GEO数据挖掘全流程分析

声明：以下学习资料根据“生信技能树”网络系列免费教学材料整理而成，代码来自“生信技能树”校长jimmy的github。GEO数据库挖掘系列知识分享课程，于2016年首发于生信菜鸟团博客。配套教学视频在B站，特此声明。

前言：关于GEO数据

我们的目标是要从读懂文献到复刻文献实验，再到掌握GEO数据挖掘的能力。首先便是要广泛阅读，在读文献时，提炼脉络，读懂文献使用了哪个或哪些GSE数据集，对数据做了哪些处理。了解清楚后，便可下载相应的数据集，得到表达矩阵，作差异分析，注释等一系列下游分析。
一篇文章可以有一个或多个GSE数据集，一个GSE里可以有一个或多个GSM样本。多个研究的GSM样本可以根据研究目的整合为一个GDS，每个数据集有着自己对应的芯片平台（GPL），一个GSE里可能有多个平台测出的数据。
本分析是基于R语言的平台，所以需要一些R语言的基础知识。
了解GEO，表达芯片与R

第一部分：GEO芯片数据下载及整理

GEO官网
本例以GSE42872数据集为例，学习GEO数据挖掘分析，分析文献

通过阅读文献找到相应的GSE数据集，并在官网可以下相应的数据集信息及背景知识
https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE42872

这个数据集使用的是GPL6244这个芯片平台，由6个样本组成，前三个为对照组，后三个为处理组。了解了这个数据集的相关背景后，接下来就需要进行数据下载了。数据下载的方式有很多种，这里我们用R包GEOquery来下载，具体下载代码如下：

rm(list = ls())  
## 魔幻操作，一键清空~当前环境中对象全部删除
options(stringsAsFactors = F)
#在调用as.data.frame的时，将stringsAsFactors设置为FALSE可以避免character类型自动转化为factor类型
f='GSE42872_eSet.Rdata'
#把GSE42872_eSet.Rdata赋值给f，方便后面流程化处理
##根据数据集不同修改相应的GSE号

library(GEOquery)
# 这个包需要注意两个配置，一般来说自动化的配置是足够的。
#Setting options('download.file.method.GEOquery'='auto')
#Setting options('GEOquery.inmemory.gpl'=FALSE)

if(!file.exists(f)){
  gset <- getGEO('GSE42872', destdir=".",
                 AnnotGPL = F,     ## 注释文件
                 getGPL = F)       ## 平台文件
  save(gset,file=f)   ## 保存到本地
}
##这是一个函数，利用包将数据集的表达信息下载下来，赋值给了gset，而不下载注释信息和平台信息，病保存到本地，文件名为f。

load('GSE42872_eSet.Rdata') 
 ## 载入数据
class(gset) 
 #查看数据类型
length(gset) 
 ##看一下有几个元素
gset[[1]]
#取第一个元素
class(gset[[1]])
 #查看改元素的数据类型
# 因为这个GEO数据集只有一个GPL平台，所以下载到的是一个含有一个元素的list

a=gset[[1]] 
##取出第一个元素赋值给一个对象a
dat=exprs(a) 
#a现在是一个对象，取a这个对象通过看说明书知道要用exprs这个函数，该函数得到表达矩阵
#现在 得到的dat就是一个表达矩阵，只不过基因的ID是探针名
dim(dat)
#看一下dat这个矩阵的维度
dat[1:5,1:5] 
#查看dat这个矩阵的1至4行和1至4列，逗号前为行，逗号后为列
#这个表达矩阵是已经log之后的，表达量一般是0-10左右，如果是原始芯片表达的信号值一般是几千到一万，则需要log处理。

boxplot(dat,las=2) 
#画个图看一下各样本之间有没有批次效应，一般中位数都差不多，las是将横坐标样本信息竖着排列

pd=pData(a) 
#通过查看说明书知道取对象a里的临床信息用pData
View(pd)
## 查看一下，挑选一些感兴趣的临床表型，这里我们欲得到其分组title信息。
library(stringr)
#运行一个字符分割包
group_list=str_split(pd$title,' ',simplify = T)[,4]
#抽取title一列，按照空格分割，取第四个元素即Control和Vemurafenib
table(group_list)
#看一下两个分组各有几个

也可以使用中国镜像下载：

library(devtools)
install_github("jmzeng1314/A