我们继续完成pySCENIC的分析!本来想这一节可视化也讲了,但是不着急,我发现有些伙伴没搞明白原由,或者太会“衣来伸手饭来张口”,所以这里着重整理了需要下载的文件!!!
上一节说了pySCENIC的分析环境配置及安装,除了这些,还有一些必要条件,例如相关文件的下载,一些数据转化等等。
为了减轻大家的负担,文件我已经下载好了,包括人的、鼠的,以及转化文件的py脚本,已上传QQ群文件,群成员可在群里免费获取!
假设你完成了上面的步骤,那接下来的分析至少在代码上很简单,三个步骤,可能会等待一段时间,尤其是第一、二步骤,不过相比于R简直是神速。注意:建议还是用服务器(别开玩笑用免费的2G内存的服务器😂),除非你的数据不大可用≥64G内存的本机。
分析第一步:GRN---运行完得到sce.adj.csv文件
pyscenic grn --num_workers 10 \
--sparse \
--method grnboost2 \
--output sce.adj.csv \
sce.loom \
hs_hgnc_tfs.txt
#这一步的目的
#推断转录因子与提供的表达矩阵基因的共表达模块,基于grnboost2,R中时GENIE3
参考基因组的情况根据实际情况自行下载,当然我下载的也可以用,具体深入的原理有兴趣的可以去了解,我只是参考文献使用的!数据库更新了,用之前的文件会出错!
鼠的下载地址:
Index of /cistarget/databases/mus_musculus/mm10/refseq_r80
人的下载地址:
Index of /cistarget/databases/homo_sapiens/hg38/refseq_r80
分析第二步:RcisTarget---运行完得到sce.regulons.csv文件
pyscenic ctx --num_workers 10 \
--output sce.regulons.csv \
--expression_mtx_fname sce.loom \
--all_modules \
--mask_dropouts \
--mode "dask_multiprocessing" \
--min_genes 10 \
--annotations_fname motifs-v9-nr.hgnc-m0.001-o0.0.tbl \
sce.adj.csv \
hg38__refseq-r80__10kb_up_and_down_tss.mc9nr.genes_vs_motifs.rankings.feather
#这一步的目的
#进行TF-motif富集分析,识别直接靶标
#得到转录因子(TF)与其对应的直接作用的靶点,称为regulon(每一个regulon是1个TF和其调控的靶基因)
分析第三步:AUCell---运行完得到sce_SCENIC.loom文件,即分析结果
pyscenic aucell --num_workers 3 \
--output sce_SCENIC.loom \
sce.loom \
sce.regulons.csv
#这一步的目的
#使用AUCell对每个细胞的每个regulon活性进行评分。
以上就是pyscenic的分析内容了!
更多精彩内容请至我的公众号---KS科研分享与服务
我们也可以提取数据,用热图的方式呈现,这里我是用ggheatmap做的,也可以用pheatmap、complexheatmap或ggplot2做。RSS分析,查看细胞类型特异性转录因子,需要先加载seurat对象,提取metadata信息,并进行分析!上面我们展示的是转录因子在不同细胞的评分,按照这个道理,我们依然可以选定某种细胞,看样本间转录因子的差别!当然了,全部展示没有啥意义,还是可以提取数据,可视化需要的TF!更多精彩内容请至我的公众号---KS科研分享与服务。先加载需要的R包,都加载了,没毛病。
关于单细胞转录组转录因子的分析我们之前在单细胞系列讲过R语言版本的,参考:跟着Cell学单细胞转录组分析(十二):转录组因子分析,但是R语言分析起来速度非常慢,如果你动辄上万的单细胞可能要运行好几周,这显然不现实。pySCENIC则很好的解决了这个问题,分析速度很快。
**pySCENIC全部往期精彩系列首先说一句,我们之前也发过R语言版本的SCENIC,但是后来我们感觉容易出错,而且费时,所以就没有再探究过。可是总是有小伙伴喜欢跑R,然后说这里错了,那里找不见,其实我们的帖子写于2022年,但是数据库已经更新了,去官网下载新的数据库,不能无脑跑代码。回到pySCENIC,之前我们写过整个系列4篇帖子,分析可视化都是很完善了。可是近期跑的时候发现在第一步有点问题,要么跑不动,要么出错,怀疑是软件和数据库没有更新的缘故,故而更新一下测试。这个帖子主要有两部分内容。
今天我们复现一篇cell子刊的图表,这篇文章有一副关于转录因子的图表,观察这个图有什么特点呢?第一是热图是二项值热图,只有0,1两个值,我们知道,在R语言版本的SCENIC分析中,最后可以得到二项值热图,那么pyscenic的分析结果中也是可以进行二项值分析并做热图的。第二是热图左侧注释有auc值的注释,热图是按照分组split的,而且两组的celltype分布是对称的。第三则是左侧行名的展示。我们的复现结果如下,因为数据是随意挑选的,所以结果看上去没有原图那么明显,但是所有的元素我们都完成了。
最近公众号小伙伴好像扎堆做单细胞转录因子富集分析,这里还是建议使用服务器,因为自己电脑可能跑起来比较费劲。我们也在上一篇内容里面对分析进行了更新,我们提供的方法都在liunux终端的conda环境中运行。我们是在R语言里面提取的seurat单细胞的矩阵去python中分析,所以最后可视化的时候需要将R里面的文件转化为python可读可操作的对象。我们封装了两个函数,第一个是R里面数据提取seurat_to_adata.R,第二个是python中数据构建函数seurat_to_adata.py。
SCENIC (Single-Cell rRegulatory Network Inference and Clustering)是一种能够从单细胞 RNA-seq 数据(SCENIC)或单细胞 RNA-seq+单细胞 ATAC-seq 的组合(SCENIC+, 今年新发的Nature Methods, 2023, 20, 1355–1367)中同时进行转录因子推断、基因调控网络重建的方法。3.制作一个trans.py,主要目的是把提取的表达矩阵转成loom文件,python trans.py运行。
这个protocol解释了如何使用软件容器和Nextflow 管道对单细胞 RNA 测序数据执行快速 SCENIC 分析以及标准最佳实践步骤。SCENIC 重建调节子(即转录因子及其目标基因),评估单个细胞中这些发现的调节子的活性,并使用这些细胞活动模式来寻找有意义的细胞簇。在这里,我们展示了 SCENIC 的改进版本,具有多项进步。SCENIC 已经用进行了重构和重新实现,速度提高了十倍,并且已打包到容器中以方便使用。现在还可以使用表观基因组轨迹数据库以及基序来完善调节子。
接上一篇,pySCENIC分析完之后就可以进行可视化了,个人认为最重要的有三个图:rss点图,rank图,转录因子表达水平热图,3个图结合着看。3.可视化rank图,可以自己改一下plotRSS_oneSet参数让图更好看一点。4. regulon表达水平热图,转录因子表达水平同理,得到的图类似只是输入不一样。2.可视化点图:寻找cluster特异性转录因子。1.loom文件读入R,提取数据。
作为 TF 结合基序数据库由于 TF 的覆盖范围不同,并且预测算法以不同的方式对绑定进行建模,因此即使使用相似的建模策略,GRN 推理方法之间的结果也可能会有所不同。频率主义方法将事件的概率定义为该事件在大量相同实验中发生的次数的比例,而贝叶斯概率将其定义为基于观察到的数据和先前的数据对所述事件发生的置信度的度量信息。GRN 中的相互作用可以是有向的或无向的(分别表示基因之间的因果关系或缺乏因果关系)、有符号的(表示正向或负向的调节模式)和/或加权的(表示相互作用的强度)。