首页
学习
活动
专区
圈层
工具
发布
33 篇文章
1
CellChat三部曲1:使用CellChat对单个数据集进行细胞间通讯分析
2
CellChat三部曲2:使用CellChat 对多个数据集细胞通讯进行比较分析
3
CellChat 三部曲3:具有不同细胞类型成分的多个数据集的细胞通讯比较分析
4
多个单细胞亚群合并
5
如何读取单细胞数据
6
纯生信单细胞数据挖掘-全代码放送
7
单细胞测序流程(单细胞rna测序)
8
单细胞亚群比例变化和表达量差异分析
9
生信中各种ID转换
10
单细胞功能注释和富集分析(GO、KEGG、GSEA)(2021公开课配套笔记)
11
细胞亚群的生物学命名
12
scRNA包学习Monocle2
13
单细胞转录组基础分析六:伪时间分析
14
Seurat包的findmarkers函数只能根据划分好的亚群进行差异分析吗
15
​cytoscape的十大插件之二--MCODE插件
16
从零到壹:Cytoscape插件使用心得~MCODE篇
17
cytoscape的cytohubba及MCODE插件寻找子网络hub基因
18
上下调基因各自独立进行GO数据库的3分类富集(求美图代码)
19
拟时序分析的热图提取基因问题
20
单细胞亚群合并与提取(2021公开课配套笔记)
21
单细胞转录组之Seurat包全流程-数据过滤、降维分群及可视化
22
CellChat细胞通讯(二)可视化篇
23
GWAS全基因组关联分析流程(BWA+samtools+gatk+Plink+Admixture+Tassel)
24
WGCNA分析,简单全面的最新教程(在线做,但也需要懂原理)
25
统计遗传学:第九章,GWAS+群体分析+亲缘关系分析
26
干货:把知识经验整理为电子书
27
如何在箱线图添加显著性--代码分享
28
ANNOVAR 软件用法还可以更复杂
29
3DSNP 数据库 | 注释 SNP 信息
30
使用FUSION进行TWAS分析
31
R包”gwasrapidd”------快速获取GWAS Catalog数据库的信息
32
连锁不平衡小工具-----LDlink的使用教程
33
🤩 CMplot | 完美复刻Nature上的曼哈顿图(一)
清单首页生信文章详情

拟时序分析的热图提取基因问题

昨天我在单细胞天地讲解了使用monocle2进行拟时序分析的方法,基本上跟着我的代码走一波就可以学会了,当然具体参数理解需要自行发力哦,见:使用monocle做拟时序分析(单细胞谱系发育) 用法只是最基础的知识而已,更多的时候,我们需要活学活用,比如课程学员提到的问题,就是因为做不到活学活用,他想知道下面的拟时序分析的热图提取基因,学员把基因按照发育顺序绘制了热图,而这些基因被他分成了3组,想拿基因去做GO/KEGG等数据库进行功能注释,不知道如何获取基因名字。

我这里不能拿学员真实项目数据来演示,所以还是用我们的老朋友,拿scRNAseq包的表达矩阵测试,见:使用monocle做拟时序分析(单细胞谱系发育)

首先根据细胞发育谱系来绘制热图

因为前面的教程 使用monocle做拟时序分析(单细胞谱系发育) 我们已经把细胞发育情况做出来了,就是NPC细胞跟另外3种细胞从生理上就不一样,所以是单独的发育轨迹,而 “GW16” and “GW21” ,“GW21+3” 这种孕期细胞,就可以很清晰的看到时间被反映在我们的拟时序分析结果了。

大家可以重温教程 使用monocle做拟时序分析(单细胞谱系发育) 里面的4个绘图代码:

代码语言:javascript
代码运行次数:0
复制
plot_cell_trajectory(cds, color_by = "Biological_Condition")  
# 可以很明显看到细胞的发育轨迹 
plot_cell_trajectory(cds, color_by = "State")
plot_cell_trajectory(cds, color_by = "Pseudotime")
plot_cell_trajectory(cds, color_by = "State") +
  facet_wrap(~State, nrow = 1)

我们前面构建细胞发育谱系,使用的是不同Biological_Condition的细胞类型之间使用monocle找到的两千多个基因。

简单的一个函数就可以绘制热图:

代码语言:javascript
代码运行次数:0
复制
plot_pseudotime_heatmap(cds[ordering_genes,],
                        num_clusters = 3,
                        cores = 1,
                        show_rownames = T)

可以看到, 这个图就和学员求助的图一模一样啦,因为基因数量两千多个,所以画出来肯定是看不清晰的啦。

既然基因被分成了3组,想拿基因去做GO/KEGG等数据库进行功能注释,就需要获取基因名字。

这个做不出来,不怪学员,因为正常人很难想到,这个绘图函数其实是可以调整返回数据对象的,而且官网例子也没有提到。

需要看函数的帮助文档,如下:

代码语言:javascript
代码运行次数:0
复制
hclust_method    
The method used by pheatmap to perform hirearchical clustering of the rows.
hclust_method = "ward.D2"
return_heatmap    
Whether to return the pheatmap object to the user.

很明显,这个函数其实就是pheatmap的一个包装罢了,本质上也是调用 hclust 而已,使用的是ward.D2距离。

然后解析热图函数返回对象

根据帮助文档,我们修改参数,这样monocle的plot_pseudotime_heatmap函数就有返回值了,是一个对象。

代码语言:javascript
代码运行次数:0
复制
p=plot_pseudotime_heatmap(cds[ordering_genes,],
                        num_clusters = 3,
                        cores = 1,return_heatmap=T,
                        show_rownames = T)

从pheatmap的对象里面提取基因名字就很简单了,就是在p$tree_row里面

代码语言:javascript
代码运行次数:0
复制
> p$tree_row

Call:
hclust(d = d, method = method)

Cluster method   : ward.D2 
Number of objects: 2200 

就可以拿到基因名对应的cluster啦,代码如下:

代码语言:javascript
代码运行次数:0
复制
clusters <- cutree(p$tree_row, k = 3)
clustering <- data.frame(clusters)
clustering[,1] <- as.character(clustering[,1])
colnames(clustering) <- "Gene_Clusters"
table(clustering)

实际上学员提问是有问题的

因为学员之间丢出这个热图,然后咨询如何在图片里面提取基因名字,所以大家只能是问是pdf还是png的图片呢?是不是可以AI或者PS解析它,拿到基因名字呢?

如果学员是直接问:使用monocle的plot_pseudotime_heatmap函数绘制的热图里面的基因聚集成为3类,该如何提取基因名字,其实就很简单了。

下一篇
举报
领券