ixxmu / mp_duty

抓取网络文章到github issues保存
https://archives.duty-machine.now.sh/
122 stars 30 forks source link

试试 bulkRNA 做拟时序分析? #5203

Closed ixxmu closed 4 months ago

ixxmu commented 4 months ago

https://mp.weixin.qq.com/s/iYpS96rNgzs4B-DpCo8wgA

ixxmu commented 4 months ago

试试 bulkRNA 做拟时序分析? by 老俊俊的生信笔记

引言

最近看到一篇推文 普通的 Bulk 转录组也能做“拟时序”分析吗? 。觉得挺有意思,已经有三篇文献作为 bulk 拟时序的分析参考来源,去看看最初的文献思路方法如下:

什么是 bulk 拟时序?(来自上面推文)

但是在 bulk 转录组中,拟时序的对象不再是细胞,而是基因(严格意义上来说,不清楚这种分析是否也能被称为拟时序分析,在本文中暂且将其称为拟时序分析)。在 bulk 转录组中,一般只有少数样本,最多可能 20-30 个样本点,其分辨率远远达不到单细胞转录组中的几千甚至上万个细胞的分辨率(在单细胞转录组中,1 个细胞即可被认为是一个样本点)。因此用 bulk 转录组进行拟时序分析时,首先需要将少数几个样本点变成 500 或者更多个样本点,这样才能更好地根据基因表达的峰值的时间点来对基因表达的先后顺序进行排序,进而得到“拟时序”的基因。

  • 恰好上面推文对过程进行了代码探索,这里对他的工作进行优化和整理,做成一个包来进行分享也是不错的。 感谢分享。

分析思路:

  • Step1:对样本点进行 PCA 分析,并通过 PC1 和 PC2 计算样本点之间的欧式距离获得时间线,最终将时间线缩放到 0-10 之间。
  • Step2:利用 modelr 中的 loess 函数,根据新时间线,产生时间点的表达值。
  • Step3:利用 atan2 函数对基因表达进行排序,定义基因表达的先后顺序。

安装

# install.packages("devtools")
devtools::install_github("junjunlab/bulkPseudotime")

# or
remotes::install_github("junjunlab/bulkPseudotime")

介绍

输入数据还是标准化的 tpm/fpkm 类似的表达值:

library(bulkPseudotime)
library(ClusterGVis)

# load test data
data(exps)

# check
head(exps,3)
#           zygote  t2.cell  t4.cell  t8.cell   tmorula blastocyst
# Oog4   1.3132282 1.237078 1.325978 1.262073 0.6549312  0.2067114
# Psmd9  1.0917337 1.315989 1.174417 1.064756 0.8685598  0.4845448
# Sephs2 0.9859232 1.201026 1.123076 1.084673 0.8878931  0.7174088

一个 bulkPseudotime 即可完成所有分析,输出一个 bulkPseudotimeClass 对象:

psetime_res <- bulkPseudotime(expMat = exps)

数据组成:

查看样本间的距离热图:

psetime_res@sample_dist_plot

查看原始 PCA 分布拟合产生的新的 PCA 分布

psetime_res@sample_pca_plot

可以看到趋势还是很相近的,符合细胞分化的时间点。


查看排序组合的拟时序热图(4 种排列组合):

psetime_res@heatmap_list_group

pseudotime_heatmap 可以针对选择的特定排序拟时序热图进行重新绘制,order 指定顺序,这里我们选择排序结果 1:

pseudotime_heatmap(object = psetime_res,order = 1)

标记一些基因名:

mark <- sample(rownames(exps),10,replace = F)

pseudotime_heatmap(object = psetime_res,
                   order = 1,
                   markGenes = mark)

heatmap_params 用来提供 comlexheatmap::Heatmap 的参数:

pseudotime_heatmap(object = psetime_res,
                   order = 1,
                   markGenes = mark,
                   heatmap_params = list(km = 3,cluster_rows = T))

pseudotime_line 绘制基因的拟时序曲线:

genes <- sample(rownames(exps),6,replace = F)

pseudotime_line(object = psetime_res,genes = genes)

结尾

路漫漫其修远兮,吾将上下而求索。


欢迎加入生信交流群。加我微信我也拉你进 微信群聊 老俊俊生信交流群 (微信交流群需收取 20 元入群费用,一旦交费,拒不退还!(防止骗子和便于管理)) 。QQ 群可免费加入, 记得进群按格式修改备注哦。

老俊俊微信:

知识星球:


往期回顾目录

ClusterGVis 绘制 lineplot 的优化
RNAseqQC 给你的数据来个全面的 QC 检查
如何 Pull Request 到 github 贡献你的代码
ggplot 添加分类型数据双坐标轴
富集分析流星图?
关于 ggplot2 分面的几点建议及补充
关于 footprint 可视化总结
facet_sub 为 ggplot2 添加子图
单细胞小提琴图添加注释标明细胞类型
蛋白序列无序区预测