内容不是教程，均为个人学习笔记，以备后日查询
学习资源均来自生信技能树论坛，生信技能树B站视频，微信公众号：生信技能树
画热图我们需要给R两个数据：表达矩阵与样品分组信息
- 读入表达矩阵
rt<-read.table('Input.txt',header = T,row.names = "symbol")
#header=TURE:把读取的矩阵第一行变为列名
#row.names="x":把读取的矩阵第一列的变为矩阵的行名
#sep='':指定分割符
#comment='!':数据中前面标！的行都不读
- 这里使用统计学参数SD选取矩阵中表达具有差异的前1000个基因进行绘图
cg<-names(tail(sort(apply(rt,1,sd)),1000))
#循环函数apply(data,1 or 2,method):对数据data的每一行（1）或列（2）执行method
#函数sort():从小到大排列一个向量
#函数tail(dat,n)：从末尾开始取n个向量dat中的元素
#函数names():取名字
- 对矩阵rt中1000个差异基因（SD）画一个简单的热图
library(pheatmap)
pheatmap(rt[cg,])
- 针对得到的热图加入分组注释
annotation<-data.frame(group=c(rep('Normal',18),rep('OA',20)))
#在data frame中指定每一个样本的编号（行名）
a<-c("GSM3130531","GSM3130532","GSM3130533","GSM3130534","GSM3130535","GSM3130536","GSM3130537","GSM3130538","GSM3130549","GSM3130550","GSM3130551","GSM3130552","GSM3130553","GSM3130554","GSM3130555","GSM3130556","GSM3130557","GSM3130558","GSM3130539",
"GSM3130540","GSM3130541","GSM3130542","GSM3130543","GSM3130544","GSM3130545","GSM3130546","GSM3130547","GSM3130548","GSM3130559",
"GSM3130560","GSM3130561","GSM3130562","GSM3130563","GSM3130564","GSM3130565","GSM3130566","GSM3130567","GSM3130568")
rownames(annotation)<-a
tmp<-annotation
class(tmp)
#tmp即为热图需要的分组信息，格式必须对应，如下图
object:tmp
- 根据分组注释信息画热图
pheatmap(rt[cg,],annotation_col = annotation)
pheatmap(rt[cg,],annotation_col = annotation,show_colnames = F,show_rownames = F)
以上数据中的表达值是一个绝对表达值(0到15)，这样的情况下表达值的高低在热图中存在基因与基因之间的比较。但很多时候我们只需要看某一个基因在不同样本中的表达差异即可，不需要基因与基因之间进行对比，因此我们需要进行normalization。
- 归一化
#在矩阵rt中取出差异表达基因（cg）转置，scale,再转置回去
n<-t(scale(t(rt[cg,])))
#使用rt["TSTD1",]检查某一个基因（这里取的是TSTD1）在scale前后是否发生变化
#将大于2和小于-2的数据拉平，避免因为矩阵中某一个值过大出现其他基因的颜色没有差异
n[n>2]=2
n[n<-2]=-2
pheatmap(n,annotation_col = annotation)
Ps.现在很多数据分析都使用R包走固定的标准流程，因此在画热图时是否需要对表达矩阵进行归一化是基于所使用的R包是否已经对数据进行了处理。通过观察数据的特征来确认是否需要手动归一化后再画热图。
