GEO芯片数据下载及清洗

1. 获得平台探针注释信息

#source("https://bioconductor.org/biocLite.R")
#biocLite("GEOquery")
setwd("/Users/baiyunfan/desktop/GEO")
library(GEOquery)
GPL <- getGEO("GPL17843")
probe<-GPL@dataTable@table

这个probe矩阵里就有详细的探针信息啦～其中第15行是gene_symbol，接下来会用到

2. 获取表达谱数据和表型矩阵

expr.df <- read.table(file = "GSE84957_series_matrix.txt", header =TRUE, 
                      comment.char = "!", row.names=1)

Data <- getGEO(filename="GSE84957_series_matrix.txt")
pData <- pData(phenoData(Data))

这里的expr.df是读取我们预先下载好的表达谱矩阵，而getGEO是直接在线读取数据，二者选取一种方法即可～

3. 将探针和symbol对应上

symbol<-sapply(1:dim(expr.df)[1],function(i){as.character(probe[which(probe[,1]==rownames(expr.df)[i]),15])})
expr.df<-cbind(symbol,expr.df)

4. 由于多个探针对应着一个基因，所以把相同的基因的表达谱数据取平均数

更正：感谢@生信技能树@土豆学生信 的友情指导。经深入学习后发现，多个探针对应同一基因的情况可以取最大值，平均数和中位数，各有优缺点，用什么方法目前尚无定论。平均数是比较保守，本实验暂且用平均值来计算。

uni<-as.character(unique(expr.df[,1]))
uni<-uni[-2]
test<-matrix(data=NA,ncol=19,nrow=1)
test<-as.data.frame(test)
colnames(test)<-colnames(expr.df)
expr.df[,1]<-as.character(expr.df[,1])
expr.df<-expr.df[-which(expr.df[,1]==""),]
for(i in 1:length(uni)){if(length(which(expr.df[,1]==uni[i]))!=1){test<-rbind(test,c(uni[i],colMeans(expr.df[which(expr.df[,1]==uni[i]),-1])))}else{test<-rbind(test,expr.df[which(expr.df[,1]==uni[i]),])}}
test<-test[-1,]

就得到最终矩阵啦～

最后安利一波生信技能树的B站良心课程，免费的，免费的，免费的，质量好，质量好，质量好。