1催什、下載GSE115850數據
2、讀取數據且處理數據宰睡，數據格式是這樣的：
1c02eb8a0eb360136b8c88dd9a3ab45.png
接下來是要把行名改成基因名蒲凶，又因為行名只能是唯一的不能有重復气筋，所以要把多個Ensembl對應的同一個基因Symbol，取平均值(一般都是取平均)旋圆，這一步代碼如下：
data<-read.delim('./GSE115850.txt') ##讀取數據
head(data) ##查看數據頭6行
rownames(data)<-data$Ensembl ##把Ensembl變?yōu)樾忻?data<-subset(data,select = -1) ##把Ensembl一行可以刪掉了
data<-aggregate(.~Symbol,data,mean) ##把重復的Symbol取平均值
rownames(data)<-data$Symbol ##把Symbol變?yōu)樾忻?data<-subset(data,select = -1) ##刪除Symbol行
data<-log2(data+1) ##對表達矩陣log2處理
處理完數據后格式是這樣的：
2c2fa06f0eaffef11f179d8c121de02.png
接下來是用limma包對這個表達矩陣做差異分析宠默，列名就是它的分組信息
#加載R包
library(limma)
library(dplyr)
library(pheatmap)
library(ggplot2)
#創(chuàng)建分組信息
phe<-data.frame(sample=colnames(p),group=c(rep('ctr',3),rep('gcm',3)))
rownames(phe)<-phe$sample #改行名
phe<-subset(phe,select = -1)
design<-model.matrix(~0+factor(phe$group))
rownames(design)<-rownames(phe)
colnames(design)<-levels(factor(phe$group))
contrast.matrix<-makeContrasts('ctr-gcm',levels=design)
常規(guī)的limma包差異分析
fit<-lmFit(p,design)
fit2<-contrasts.fit(fit,contrast.matrix)
fit2<-eBayes(fit2)
DEG<-topTable(fit2,coef=1,n=Inf) %>% na.omit()
write.table(DEG,file='./GSE115850allgene.txt',sep='\t')
a<-DEG[which(DEG$P.Value<0.05&abs(DEG$logFC)>0.58),]##篩選FC>1.5和p<0.05
write.table(a,file='./GSE115850deg.txt',sep='\t')
3、輸出的差異基因列表
8f70a3d175e537c11427ff08458d7da.png
4灵巧、開始畫熱圖
exp<-p[intersect(rownames(a),rownames(p)),]
pheatmap<- pheatmap(exp,cluster_cols=F,annotation_col=phe,show_colnames=F,scale='row',show_rownames=F)
ggsave(file='DEG_heatmap.pdf',pheatmap,width=30,height=30,units='cm')
5搀矫、結果展示，效果蠻不錯
82ed70554ac9aea4295be270cdb09cb.png
自己也是學習階段刻肄，大家多多指教瓤球，謝謝！
