基于 python 的单细胞转录因子分析_TTS56的博客

相关文章推荐
读研的充电器 · mysql日期只比较年月日 ...· 1 年前 ·
不羁的烤红薯 · Charles + Android + ...· 2 年前 ·
果断的火柴 · VS ...· 2 年前 ·
爱逃课的硬盘 · Common Event Format ...· 2 年前 ·
知识渊博的冰棍 · 无法启动Zabbix ...· 2 年前 ·
Download reference datas
wget -c https://github.com/aertslab/pySCENIC/archive/refs/heads/master.zip
x master.zip
cd master
mv resources/* ../../
wget -c https://resources.aertslab.org/cistarget/motif2tf/motifs-v9-nr.hgnc-m0.001-o0.0.tbl
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/tc_v1/gene_based/encode_20190621__ChIP_seq_transcription_factor.hg19-tss-centered-5kb.max.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/tc_v1/gene_based/encode_20190621__ChIP_seq_transcription_factor.hg19-500bp-upstream.max.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/tc_v1/gene_based/encode_20190621__ChIP_seq_transcription_factor.hg19-tss-centered-10kb.max.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg38/refseq_r80/tc_v1/gene_based/encode_20190621__ChIP_seq_transcription_factor.hg38__refseq-r80__10kb_up_and_down_tss.max.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg38/refseq_r80/tc_v1/gene_based/encode_20190621__ChIP_seq_transcription_factor.hg38__refseq-r80__500bp_up_and_100bp_down_tss.max.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc8nr/gene_based/hg19-500bp-upstream-7species.mc8nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-500bp-upstream-7species.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc8nr/gene_based/hg19-500bp-upstream-10species.mc8nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-500bp-upstream-10species.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg38/refseq_r80/mc9nr/gene_based/hg38__refseq-r80__500bp_up_and_100bp_down_tss.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc8nr/gene_based/hg19-tss-centered-10kb-7species.mc8nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-tss-centered-10kb-7species.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg38/refseq_r80/mc9nr/gene_based/hg38__refseq-r80__10kb_up_and_down_tss.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc8nr/gene_based/hg19-tss-centered-10kb-10species.mc8nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-tss-centered-10kb-10species.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc8nr/gene_based/hg19-tss-centered-5kb-7species.mc8nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-tss-centered-5kb-7species.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-tss-centered-5kb-7species.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc8nr/gene_based/hg19-tss-centered-5kb-10species.mc8nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-tss-centered-5kb-10species.mc9nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc8nr/region_based/hg19-regions-9species.all_regions.mc8nr.feather
wget -c https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/region_based/hg19-regions-9species.all_regions.mc9nr.feather
The pipline of pyscenic only 3 steps 
Step.1 
pyscenic grn \
        --num_workers 6 \
        -o /data/expr_mat.adjacencies.tsv \ # input Count data
        # csv (rows=cells x columns=genes) or loom (rows=genes x columns=cells).
        /data/expr_mat.tsv \
        /data/allTFs_hg38.txt
Step.2 
pyscenic ctx \
        /data/expr_mat.adjacencies.tsv \ # First Step out put file
        /data/hg19-tss-centered-5kb-7species.mc9nr.feather \
        /data/hg19-tss-centered-10kb-7species.mc9nr.feather \
        --annotations_fname /data/motifs-v9-nr.hgnc-m0.001-o0.0.tbl \
        --expression_mtx_fname /data/expr_mat.tsv \ # the same to the first input data
        --mode "dask_multiprocessing" \
        --output /data/regulons.csv \
        --num_workers 6
Step.3 
pyscenic aucell \
        /data/expr_mat.tsv \
        /data/regulons.csv \
        -o /data/auc_mtx.csv \
        --num_workers 6
pyscenic
micromamba activate SCpip install pyscenic -i https://mirrors.aliyun.com/pypi/simple/
安装docker
需要有root权限或者在docker的用户组
#1.Update the apt package index and inst
                                    chromVAR是一个用于分析稀疏染色质开放的R包。chromVAR的输入文件包括，ATAC-seq处理后的fragments文件(过滤重复和低质量数据), DNAse-seq实验结果，以及基因组注释(例如motif位置)
chromVAR先根据所有细胞或者样本的平均情况来计算期望开放性， 然后用它来计算每个注释，每个细胞或样本的偏差，最后对开放进行纠正。
安装加...
关于单细胞转录组转录因子的分析我们之前在单细胞系列讲过R语言版本的，参考：跟着Cell学单细胞转录组分析(十二):转录组因子分析，但是R语言分析起来速度非常慢，如果你动辄上万的单细胞可能要运行好几周，这显然不现实。pySCENIC则很好的解决了这个问题，分析速度很快。
                                    转录因子分析可以了解细胞异质性背后的基因调控网络的异质性。转录因子分析也是单细胞转录组常见的分析内容，R语言分析一般采用的是SCENIC包，具体原理可参考两篇文章。1、《SCENIC : single-cell regulatory networkinference and clustering》。2、《Ascalable SCENIC workflow for single-cell gene regulatory network analysis》。但是说在前头，SCENIC的计算量超级大，非常耗费内存
                                    假设你完成了上面的步骤，那接下来的分析至少在代码上很简单，三个步骤，可能会等待一段时间，尤其是第一、二步骤，不过相比于R简直是神速。参考基因组的情况根据实际情况自行下载，当然我下载的也可以用，具体深入的原理有兴趣的可以去了解，我只是参考文献使用的！为了减轻大家的负担，文件我已经下载好了，包括人的、鼠的，以及转化文件的py脚本，已上传QQ群文件，群成员可在群里免费获取！上一节说了pySCENIC的分析环境配置及安装，除了这些，还有一些必要条件，例如相关文件的下载，一些数据转化等等。
                                    文章目录一、安装二、使用1、准备工作2、预处理过滤低质量细胞样本3、检测特异性基因4、主成分分析（Principal component analysis）5、领域图，聚类图（Neighborhood graph）6、检索标记基因7、保存数据8、番外
如果没有conda 基础，参考： Conda 安装使用图文详解（2021版）
pip install scanpy
conda install -y -c conda-forge leidenalg
1、准备工作
# 载入包
import
                                    ILC1、 ILC2、ILC3差异基因表达分析，结果表明，ILC1 cells共有79个上调表达基因，参与干扰素γ的调控，ILC2 cells共有58个上调表达的基因，在前列腺素、Notch信号通路及环境感应中发挥作用，ILC3 cells共有371个上调表达的基因，根据GO注释有85个与免疫相关，且存在未知功能的基因。GSEA基因集富集分析。基因表达趋势分析（需要三个或以上有处理梯度的样品）：利用STEM软件将基因按照表达模式分为不同的趋势、各趋势基因的GO/KEGG功能富集分析、特定基因集聚类分析；
                                    第一次接触基因评分是在一篇文章中，也不知道这样的叫法对不对，作者选定了几个炎症基因，利用seurat包的一个打分函数AddModuleScore，依据基因的平均表达水平进行分析，最后得到的score称为炎症分数，其实这样的叫法有点欠缺，但是有这样的做法。（演示数据没有意义）
选择基因进行计算。
DefaultAssay(immune) <- “RNA”
cd_features <- list(c(
‘TNF’,
‘CCL2’,
‘CCL3’,
‘CCL4’,
‘CXCL10’,
‘S100A8’
                                    在2019/08/07的Nature刊中，中科院景乃禾课题组发表了文章——Molecular architecture of lineage allocation and tissue organization in early mouse embryo ，我在这篇文章中发现了一个被汤神组 （就是Hemberg-lab单细胞转录组数据分析（二）- 实验平台中开辟了单细胞转录组领域的人）反...