截屏2021-08-11 下午10.02.47的副本.png

PCA样本聚类图

rm(list = ls())  
load(file = "step1output.Rdata")
load(file = "step2output.Rdata")
#输入数据：exp和Group
#Principal Component Analysis
#http://www.sthda.com/english/articles/31-principal-component-methods-in-r-practical-guide/112-pca-principal-component-analysis-essentials

1.PCA 图----

dat=as.data.frame(t(exp))
library(FactoMineR)
library(factoextra) 
dat.pca <- PCA(dat, graph = FALSE)      #进行pca分析
pca_plot <- fviz_pca_ind(dat.pca,       #经pca分析后得到的数据
                         geom.ind = "point", # show points only (nbut not "text")
                         col.ind = Group, # color by groups跟据分组区分颜色
                         palette = c("#00AFBB", "#E7B800"),
                         addEllipses = TRUE, # Concentration ellipses
                         legend.title = "Groups"
)
pca_plot
#处理数据时需要提供的：数据本身，及分组数据
ggsave(plot = pca_plot,filename = paste0(gse_number,"_PCA.png"))
save(pca_plot,file = "pca_plot.Rdata")

2.top 1000 sd 热图----

cg=names(tail(sort(apply(exp,1,sd)),1000))  #从exp数据集中挑出来变化大的1000个基因
n=exp[cg,]   #根据1000个挑选基因的名字选子集，获得矩阵

直接画热图，对比不鲜明（同一样本不同基因进行比较）
library(pheatmap)
annotation_col=data.frame(group=Group)
rownames(annotation_col)=colnames(n) 
pheatmap(n,
         show_colnames =F,
         show_rownames = F,
         annotation_col=annotation_col
)

用标准化的数据画热图，两种方法的比较：https://mp.weixin.qq.com/s/jW59ujbmsKcZ2_CM5qRuAg
1.使用热图参数
pheatmap(n,
         show_colnames =F,
         show_rownames = F,
         annotation_col=annotation_col,
         scale = "row",            #把矩阵按“row”进行标准化，行与行之间不对比（基因与基因之间），对比同一基因不同分组表达差异
         breaks = seq(-3,3,length.out = 100)  #设置颜色分配范围，小与-3与-3颜色相同，反之成立，使颜色更加鲜明
         ) #breaks 参数解读在上面链接
dev.off()

2.自行标准化再画热图
n2 = t(scale(t(n)))  #scale只能按列进行标准化，需要对表达矩阵进行转化，基因名转置到列上，然后再转换为表达矩阵
pheatmap(n2,
         show_colnames =F,
         show_rownames = F,
         annotation_col=annotation_col,
         breaks = seq(-3,3,length.out = 100)
)
dev.off()

#对于颜色设置范围，可在进行标准化之后boxplot(矩阵)看数值分布范围，根据需要进行取值。避免离群值占比颜色鲜艳导致中间数据区别不大

关于scale的进一步探索：zz.scale.R

3.相关性热图----

pheatmap::pheatmap(cor(exp),#原始数据
                   annotation_col = annotation_col)

pheatmap::pheatmap(cor(n),#1000个表达差异高的基因
                   annotation_col = annotation_col)

pheatmap::pheatmap(cor(n2),#标准化后的表达矩阵
                   annotation_col = annotation_col
                   )

dev.off()

关于相关性背后的故事：https://mp.weixin.qq.com/s/IqMW6Qjf64dn30F4RQg5kQ