关于单细胞转录组转录因子的分析我们之前在单细胞系列讲过R语言版本的,参考:跟着Cell学单细胞转录组分析(十二):转录组因子分析,但是R语言分析起来速度非常慢,如果你动辄上万的单细胞可能要运行好几周,这显然不现实。pySCENIC则很好的解决了这个问题,分析速度很快。官方教程参考:
https://pyscenic.readthedocs.io/en/latest/
一、软件安装
老样子,还是先说一下安装和分析文件的准备,前面环境的配置和之前cellphonedb一样,如果已经操作过的,可以跳过:
#安装下载及环境设置
# 安装一个conda,为什么安装他可以理解为Rstuido之于R,后期在环境设置、软件安装上很方便。
wget https://repo.anaconda.com/miniconda/Miniconda3-latest-Linux-x86_64.sh
# bash安装,按照指引,都选yes,这样一些依赖的python包都安装了。
bash Miniconda3-latest-Linux-x86_64.sh
# 激活conda环境
source ~/.bashrc
#设置镜像
conda config --add channels r
conda config --add channels conda-forge
conda config --add channels bioconda
conda config --add channels https://mirrors.bfsu.edu.cn/anaconda/cloud/bioconda/
conda config --add channels https://mirrors.bfsu.edu.cn/anaconda/cloud/conda-forge/
conda config --add channels https://mirrors.bfsu.edu.cn/anaconda/pkgs/free/
conda config --add channels https://mirrors.bfsu.edu.cn/anaconda/pkgs/main/
conda config --set show_channel_urls yes
然后就是安装pySCENIC及环境设置:
#安装pyscenic,并创建分析环境
conda create -n pyscenic python=3.9#创建一个pyscenic 的python环境,pyscenic要求python版本3.6及以上,目前python出到3.9了,我用3.9
conda activate pyscenic #激活pyscenic 环境
#安装依赖包
conda install -y numpy
conda install -y -c anaconda cytoolz
conda install -y scanpy
#安装pyscenic
pip install pyscenic
pip安装即可,如果安装失败或者觉得下载速度太慢,可以使用下面代码,使用与其他软件的pip安装。
pip install cellphonedb -i http://pypi.douban.com/simple/ --trusted-host pypi.douban.com
二、数据文件准备
一些分析用的依赖数据库及文件:
#TF注释
#鼠的
wget https://resources.aertslab.org/cistarget/motif2tf/motifs-v9-nr.mgi-m0.001-o0.0.tbl
#人的
wget https://resources.aertslab.org/cistarget/motif2tf/motifs-v9-nr.hgnc-m0.001-o0.0.tbl
#转录组因子列表
下载地址
https://github.com/aertslab/pySCENIC/tree/master/resources
#人的文件名:hs_hgnc_tfs.txt,复制为txt文件即可
#鼠的文件名:mm_mgi_tfs.txt,复制为txt文件即可
#reference数据库,之前一些网上教程的链接文件已经不行了,因为做了更新,跑的时候会出错,我是根据报错选择了下面的文件
#鼠的
wget https://resources.aertslab.org/cistarget/databases/mus_musculus/mm10/refseq_r80/mc_v10_clust/gene_based/mm10_10kbp_up_10kbp_down_full_tx_clustered.genes_vs_motifs.rankings.feather
#人的
wget https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg38/refseq_r80/mc_v10_clust/gene_based/hg38_10kbp_up_10kbp_down_full_tx_v10_clust.genes_vs_motifs.rankings.feather
三、分析文件准备
这一步是从我们的R种seurat对象提取表达矩阵,并转化为loom文件。这里使用一个python脚本,转化一下即可。
#文件准备
#pyscenic的输入文件是行为基因名,列为细胞ID的矩阵,所以在seurat对象中导出矩阵的时候需要转置一下,可以用标准化矩阵,也可以用counts矩阵,影响不大!
#表达矩阵、meta----R中进行
write.csv(t(as.matrix(sce@assays$RNA@counts)),file = "sce_exp.csv")
#cellInfo <- sce@meta.data[,c("celltype","nCount_RNA","nFeature_RNA")]
#colnames(cellInfo) <- c('CellType', 'nGene' ,'nUMI')
#head(cellInfo)
#write.csv(cellInfo, file = "cellInfo.csv")
#转化为loom文件,Linux下的python脚本
#编辑脚本
vim trans.py
#输入以下内容
import os, sys
os.getcwd()
os.listdir(os.getcwd())
import loompy as lp;
import numpy as np;
import scanpy as sc;
x=sc.read_csv("sce_exp.csv");#R中导出的表达矩阵
row_attrs = {"Gene": np.array(x.var_names),};
col_attrs = {"CellID": np.array(x.obs_names)};
lp.create("sce.loom",x.X.transpose(),row_attrs,col_attrs)
#保存并退出
#运行trans.py
python trans.py
ls
#这样在文件夹中会出现sce.loom文件,就是接下来输入pyscenic的文件。
这样整个准备工作就完成了,后续的分析很简洁,只需要三步。然后就是可视化了,这个下回分解