想看看你的差异分析跟作者结果是不是一致

生信菜鸟团

发布于 2025-07-31 14:45:03

9300

代码可运行

运行总次数：0

代码可运行

2022的文章：《Blocking ActRIIB and restoring appetite reverses cachexia and improves survival in mice with lung cancer》，其中一个火山图来源于一个转录组测序的差异分析：

Volcano plot depicting INHBA as one of the most differentially expressed genes (DEG) in the RNA-Seq of tumors from +CACS mice compared to -CACS (P < 0.01 & Log2Fold Change >1.2).
Principal Component Analysis (PCA) of RNA-Seq from tumors of -CACS (n = 9) and +CACS (n = 13)

需要自己读一下文章，然后找到了作者给出来的表达量矩阵文件：GSE165856_Marcus5291_2017_11_HTSeq.gene.counts.csv.gz，接下来就可以简单的转录组测序的差异分析一条龙啦啦！

同样的，有了差异分析结果后，文章里面的这样的火山图上面的基因列表，并不需要自己一个个手敲代码啦，截图后复制粘贴给人工智能大模型即可：

火山图上面的基因列表

可以看到，人工智能大模型很容易就分析好了我们上传的图片，给出来了图片里面的基因列表。

同样的继续让人工智能大模型帮忙绘制热图，这个时候需要载入前面的转录组测序的差异分析流程拿到的表达量矩阵和分组信息啦：

  load('deg/DEG_deseq2.Rdata')
  head(DEG_deseq2)
  DEG =DEG_deseq2 
  load(file = 'step1-output.Rdata') 
  paper_list=list(
    up_paper = c("Fst4", "F830016B08Rik", "4930565N06Rik", "Lcor", "Far1os", "Inhba", "Col10a1", "Aldh1a3", "Mmp10", "Serpinb1", "Tff3") ,
    down_paper = c("Ackr4", "Bmp6", "Ces1g", "Cyp1a1", "Gm8229", "Gm5930", "Slc5a12", "Mytt1", "Asgr1")
  )

  paper_list
library(stringr)  
  df = data.frame(
    group=rep(names(paper_list),
              times=unlist(lapply(paper_list, length))),
    gene=unlist(paper_list)
  ) 
library(stringr) 
  df=df[df$gene %in% rownames(dat),]
  cg=df$gene
  cg
library(pheatmap)
#tmp=dat[cg,]
  n=t(scale(t(dat[cg,]))) # 'scale'可以对log-ratio数值进行归一化
  n[n>2]=2
  n[n< -2]= -2
  n[1:4,1:4]
  pheatmap(n,show_colnames =F,show_rownames = F)
  ac=data.frame(group=group_list)
  ar=data.frame(gene=df$group)
  rownames(ac)=colnames(n)
  rownames(ar)=rownames(n)
  pheatmap(n,show_colnames =F,show_rownames = T,   
           annotation_row =ar,
           annotation_col=ac)

可以看到，人工智能大模型给我们的文献里面的上下调基因，确实是在我们的复现里面体现出来了同样的上下调趋势。