Hello禾酱，大家好！此前得封，小編給大家介紹了：如何從GEO數(shù)據(jù)庫下載芯片數(shù)據(jù)及相關(guān)樣本處理信息等等（不知道的可以戳這里哦）指郁。這些芯片數(shù)據(jù)下載下來干嘛呢？下載必然是為了挖數(shù)據(jù)啦闲坎！指不定什么有意思的東西我就發(fā)現(xiàn)了呢？指不定老天爺眼睛一閉讓我發(fā)了文章了呢梗逮？
言歸正傳，拿到數(shù)據(jù)娄蔼，分析的第一步往往是進行基因差異表達分析底哗，所以，針對芯片數(shù)據(jù)涕癣，我們就來給大家介紹一款基因差異表達分析的常用方法——R包limma前标。
數(shù)據(jù)簡介與設(shè)置
為了方便演示，這里選擇了人的早幼粒細胞白血病細胞系NB4細胞的六個樣本數(shù)據(jù)（GSE2600）炼列，分析的輸入文件是下載的表達矩陣文件，而分析之前需要確保正確安裝和加載limma须蜗，同時需要對工作路徑進行設(shè)置目溉。
library('limma')
workdir="F:/GEO/20180520"
setwd(workdir)
數(shù)據(jù)處理
1、表達矩陣
數(shù)據(jù)為六個樣本柿估，讀取數(shù)據(jù)之后陷猫，大家可以利用head()簡單查看數(shù)據(jù)的情況等。
>?expreSet=read.csv2("GSE2600expressionMatrix.csv",?header?=T,?row.names?=1,check.names?=F)
>head(exprSet,3)
GSM49939GSM49940GSM49941GSM49942GSM49943GSM49944
1007_s_at23.013.826.575.994.984.6
1053_at1449.91826.72242.81508.81523.02355.5
117_at109.271.5106.7128.884.179.6
針對表達矩陣足陨，需要查看其整體分布情況娇未，可以利用boxplot()繪制box分布圖，GEO下載的表達矩陣數(shù)據(jù)基本上都是標(biāo)準(zhǔn)化的數(shù)據(jù)镊讼，可以由箱線圖的分布特點看出這些樣本的數(shù)據(jù)基本分布一致（中位數(shù)、上四分位數(shù)蝶棋、下四分位數(shù)等等），如下圖結(jié)果：
#獲取樣品數(shù)量兼贸，并設(shè)置圖片顏色
n.sample?=?ncol(exprSet)
cols?=?rainbow(n.sample)
#利用boxplot()繪圖
pdf(file=paste(workdir,"/","Probe_expressionDistribution.pdf",sep=""),?width=24,?height=18)
par(cex?=0.7)
if(n.sample>40)?par(cex?=0.5)
boxplot(exprSet,col?=?cols,?main?="expression",?las?=2)
dev.off()
2寝受、分組矩陣
確認表達矩陣之后罕偎，可以由下載保存的樣本處理信息進行分組京闰，例如此處的樣本處理分組：CONTROL/INFECTED，經(jīng)過整理蹂楣，分組信息大致如下,并基于分組信息構(gòu)建分組矩陣（design):
>group
Treatment
GSM49939???CONTROL
GSM49940???CONTROL
GSM49941???CONTROL
GSM49942??INFECTED
GSM49943??INFECTED
GSM49944??INFECTED
>?design?=?model.matrix(~?Treatment?+0,group)
>?colnames(design)?=?levels(as.factor(c("CONTROL","INFECTED")))
>?design
CONTROL?INFECTED
GSM4993910
GSM4994010
GSM4994110
GSM4994201
GSM4994301
GSM4994401
attr(,"assign")
[1]11
attr(,"contrasts")
attr(,"contrasts")$Treatment
[1]"contr.treatment"
3痊土、差異比較矩陣
基于分組矩陣的信息構(gòu)建差異比較矩陣（cont.matrix)，由差異比較矩陣顯示結(jié)果可知犯祠，是進行INFECTED 與CONTROL之間的差異分析酌呆。
>cont.matrix = makeContrasts(INFECTED-CONTROL, levels=design)
>cont.matrix
Contrasts
LevelsINFECTED-CONTROL
CONTROL-1
INFECTED1
差異表達分析
差異表達分析主要是基于lmFit()、eBayes()隙袁、topTable()完成分析過程菩收，并提取了主要的結(jié)果(tT)。
>?fit?=?lmFit(exprSet,?design)
>?fit2?=?contrasts.fit(fit,?cont.matrix)
>fit2?=?eBayes(fit2,0.01)
>tT?=?topTable(fit2,?adjust="fdr",?sort.by="logFC",?resort.by="P",n=Inf)
>
tT?=?subset(tT,?select=c("adj.P.Val","P.Value","logFC"))
>head(tT,15)
adj.P.ValP.ValuelogFC
223020_at0.999642.196175e-05746.10000
1555758_a_at0.999646.467722e-05-540.53333
218676_s_at0.999641.352768e-04-280.86667
237249_at0.999642.669173e-04-93.53333
225100_at0.999642.836527e-04-124.96667
217825_s_at0.999642.903446e-04-143.73333
222099_s_at0.999643.425427e-04493.13333
212634_at0.999644.221452e-04-166.06667
211499_s_at0.999644.391776e-04-129.56667
221098_x_at0.999644.805746e-0495.16667
208974_x_at0.999645.060448e-04947.76667
209670_at0.999645.113338e-04374.20000
202088_at0.999645.262646e-04-594.40000
219394_at0.999645.307063e-04-117.56667
212221_x_at0.999645.393084e-04347.43333
以上坡贺，就完成了分析過程，之后可以根據(jù)結(jié)果設(shè)定合適的閾值篩選出差異表達基因钧萍，并基于結(jié)果進行圖形化展示政鼠，或者進行富集分析、蛋白質(zhì)互作網(wǎng)絡(luò)分析等等公般，如此可以豐富分析結(jié)果，感興趣的同學(xué)可以去探索一下哦瞬雹！
更多生物信息課程：
1. 文章越來越難發(fā)刽虹？是你沒發(fā)現(xiàn)新思路，基因家族分析發(fā)2-4分文章簡單快速胖缤，學(xué)習(xí)鏈接：基因家族分析實操課程阀圾、基因家族文獻思路解讀
2. 轉(zhuǎn)錄組數(shù)據(jù)理解不深入？圖表看不懂初烘？點擊鏈接學(xué)習(xí)深入解讀數(shù)據(jù)結(jié)果文件肾筐，學(xué)習(xí)鏈接：轉(zhuǎn)錄組（有參）結(jié)果解讀；轉(zhuǎn)錄組（無參）結(jié)果解讀
3. 轉(zhuǎn)錄組數(shù)據(jù)深入挖掘技能-WGCNA局齿，提升你的文章檔次，學(xué)習(xí)鏈接：WGCNA-加權(quán)基因共表達網(wǎng)絡(luò)分析
4. 轉(zhuǎn)錄組數(shù)據(jù)怎么挖掘讥此？學(xué)習(xí)鏈接：轉(zhuǎn)錄組標(biāo)準(zhǔn)分析后的數(shù)據(jù)挖掘谣妻、轉(zhuǎn)錄組文獻解讀
5.?微生物16S/ITS/18S分析原理及結(jié)果解讀、OTU網(wǎng)絡(luò)圖繪制他巨、cytoscape與網(wǎng)絡(luò)圖繪制課程
6. 生物信息入門到精通必修基礎(chǔ)課，學(xué)習(xí)鏈接：linux系統(tǒng)使用捻爷、perl入門到精通份企、perl語言高級、R語言畫圖
7. 醫(yī)學(xué)相關(guān)數(shù)據(jù)挖掘課程司志，不用做實驗也能發(fā)文章骂远，學(xué)習(xí)鏈接：TCGA-差異基因分析、GEO芯片數(shù)據(jù)挖掘激才、GSEA富集分析課程、TCGA臨床數(shù)據(jù)生存分析吨述、TCGA-轉(zhuǎn)錄因子分析钞脂、TCGA-ceRNA調(diào)控網(wǎng)絡(luò)分析
8.其他課程鏈接：二代測序轉(zhuǎn)錄組數(shù)據(jù)自主分析捕儒、NCBI數(shù)據(jù)上傳、二代測序數(shù)據(jù)解讀阎毅。
