本来用小洁老师的转录组流程代码跑的非常丝滑,找了一个公司给的转录组测序数据本想大展拳脚,然后就error了,记录一下debug的过程:

出现问题的位置

首先正常地取gene symbol的表达矩阵,去重,DEG,挑选DEG画热图的时候error了:

  1. cg1 = rownames(DEG1)[DEG1$change !="NOT"]
  2. h1 = draw_heatmap(dat[cg1,],Group,n_cutoff = 2)

image.png

自己试图debug的过程

dat[head(cg1),]#no error
dat[tail(cg1),]#no error
table(cg1 %in% rownames(dat))#all TRUE
#####改成data.frame可以取子集了...
dat0 <- dat %>% as.data.frame()
dat0[cg1,]
draw_heatmap(as.matrix(dat0[cg1,]),Group,n_cutoff = 2)

image.png
我蕉绿了,,,
优秀的小洁老师来拯救世界了

小洁老师debug的过程

table(is.na(rownames(dat)))
table(cg1 %in% rownames(dat))
table(is.na(cg1))
dat[cg1,]
d2 = dat[rownames(dat) %in% cg1 ,]
dim(d2)
table(!duplicated(cg1))
table(!duplicated(rownames(dat)))
table(!duplicated(rownames(d2)))
head(sort(cg1))
#""      "AASS"  "ABCA1" "ABCA4" "ABCG1" "ABCG2"

所以是行名里的””导致matrix不能取子集,但data.frame可以取子集,,,
只是用data.frame取了子集后,再as.matrix, pheatmap不能识别””,一整个大无语啊,,,

“”当然是从表达矩阵gene symbol里来的

#cg1是从DEG来的,DEG是从表达矩阵来的
head(sort(rownames(dat)))
#""      "AASS"  "ABCA1" "ABCA4" "ABCG1" "ABCG2"

再往前翻,在原始的gene.count.matrix里,因为我直接使用了公司给到的matrix里的gene symbol,直接!duplicated,所以

head(sort(gene.count.matrix$gene_name))
# "" "" "" "" "" ""
#所以在filter gene的时候应该加上
dat <- dat[!dat$gene_name == "",]

总结

head(sort()) yyds!

致谢

对小洁老师的感谢值得一个一级标题!!!