1
安裝所需?R包
if?(!requireNamespace("BiocManager",?quietly?=?TRUE))install.packages("BiocManager")BiocManager::install("GEOquery")BiocManager::install("limma")BiocManager::install("affy")
2
設(shè)置路徑，并引用包
rm(list?=?ls())######需要修改setwd("/Users/xxx/GSE40732")??library(GEOquery)library(limma)library(affy)
3
GEO?數(shù)據(jù)下載
gse<-getGEO("GSE40732",destdir?=?".",?getGPL?=?T,?AnnotGPL?=?T)length(gse)?#查看?GSE有幾個(gè)平臺(tái)子集#預(yù)覽第一個(gè)平臺(tái)子集??????????gse[[1]]?exp=exprs(gse[[1]])?#使用?exprs?函數(shù)從?gse[[1]]?中獲取表達(dá)矩陣數(shù)據(jù)，并將其存儲(chǔ)在?exp?變量中cli=pData(gse[[1]])?#使用?pData?函數(shù)從?gse[[1]]?中獲取臨床數(shù)據(jù)，并將其存儲(chǔ)在?cli?變量中GPL<-fData(gse[[1]])???#獲取平臺(tái)信息
4
如果出現(xiàn)Genebank?號(hào)，可作對(duì)應(yīng)修改
library(AnnotationDbi)library(org.Hs.eg.db)accession_numbers?<-?GPL[,?2]?gene_symbols?<-?mapIds(org.Hs.eg.db,????????????????????????keys?=?accession_numbers,????????????????????????column?=?"SYMBOL",????????????????????????keytype?=?"ACCNUM",????????????????????????multiVals?=?"first")GPL$Gene_symbol?<-?gene_symbols
5
獲取需要的探針和gene數(shù)據(jù)
gpl<-GPL[,c(1,4)]???head(gpl)####Genesymbol單詞可能需要修改?#對(duì)gpl對(duì)象中的Gene?symbol元素進(jìn)行處理，并將其轉(zhuǎn)換為一個(gè)數(shù)據(jù)框。gpl$`Gene_symbol`=data.frame(sapply(gpl$`Gene_symbol`,?#使用sapply函數(shù)對(duì)gpl$Gene?symbol``中的每個(gè)元素應(yīng)用一個(gè)函數(shù)?????????????????????????????????????function(x)unlist(strsplit(x,"http:///"))[1]),??????????#使用strsplit函數(shù)將每個(gè)元素按照///進(jìn)行拆分，并提取拆分后結(jié)果的第一個(gè)部分。這樣可以獲得每個(gè)元素中的第一個(gè)基因符號(hào)。?????????????????????????????????????stringsAsFactors?=?F)[,1]?????#stringsAsFactors?=?F參數(shù)指定將字符串型數(shù)據(jù)保持為字符型而不轉(zhuǎn)換為因子。#最后，通過[,1]選擇數(shù)據(jù)框的第一列，將處理后的基因符號(hào)賦值給gpl$Gene?symbol``，完成對(duì)Gene?symbol的處理。
6
?數(shù)據(jù)轉(zhuǎn)換
exp=as.data.frame(exp)???#轉(zhuǎn)換為數(shù)據(jù)框exp$ID<-rownames(exp)????#在?exp?中加一列行名為ID，填充是?exp的?rownamesexp_symbol<-merge(exp,gpl,by="ID")?#將?exp?中的?ID轉(zhuǎn)換為?gpl?中的?genesymbol?exp_symbol<-na.omit(exp_symbol)#去掉空值head(exp_symbol[,?ncol(exp_symbol)],?6)?#查看最后一列的前?6?行tail(exp_symbol[,ncol(exp_symbol)])?#查看最后一列的最后?6?行#將?ID?列中的?genesymbol?賦值給?exp_symbol?中的?rownames####Genesymbol單詞可能需要修改table?(duplicated?(exp_symbol$`Gene_symbol`))??#合并相同?genesymbol?的數(shù)值，取平均值####Genesymbol單詞可能需要修改exp_unique<-avereps(exp_symbol[,-c(1,ncol(exp_symbol))],#avereps函數(shù)用于將具有相同標(biāo)識(shí)符的行聚合為一個(gè)平均值????????????????????#選擇exp_symbol數(shù)據(jù)框中除第一列和最后一列之外的所有列作為輸入數(shù)據(jù)進(jìn)行聚合，因?yàn)榈谝涣袨?ID，最后一列是?Gene?symbol????????????????????ID=exp_symbol$`Gene_symbol`)???????????????????#表示使用Gene?symbol列作為標(biāo)識(shí)符進(jìn)行聚合#avereps?函數(shù)會(huì)返回一個(gè)新的數(shù)據(jù)框，這個(gè)數(shù)據(jù)框的每一行是一個(gè)組的平均值，行名是該組的標(biāo)識(shí)符
7
輸出并保存文件
#保存臨床數(shù)據(jù)結(jié)果write.table(cli,file="clinicalx.xls",?sep="\t",?quote=F,?col.names?=?T,?row.names?=?F)##保存表達(dá)數(shù)據(jù)文件####保存好的基因表達(dá)矩陣列名錯(cuò)開了一位，需要修改write.table(exp_unique,file="geneEXP.xls",?sep="\t",?quote=F,?col.names?=?T,?row.names?=?T)
關(guān)注公眾號(hào)并轉(zhuǎn)發(fā)朋友圈（非強(qiáng)制）獲取原始R代碼和答疑
生信Immunity
