为避免造成环境污染,建议使用conda创建隔离环境进行分析
#创建隔离环境
conda create -n 10X
#激活隔离环境
conda activate 10X
添加必要的conda镜像
conda config --add channels conda-forge
conda config --add channels defaults
conda config --add channels r
conda config --add channels bioconda
conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/pkgs/free/
conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/conda-forge/
conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/msys2/
conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/bioconda/
conda config --add channels https://mirrors.ustc.edu.cn/anaconda/pkgs/free/
conda config --add channels https://mirrors.ustc.edu.cn/anaconda/cloud/bioconda/
conda config --set show_channel_urls yes
1.原始数据下载及转换
从GEO下载原始数据需要使用官方工具sra-tools
安装sra-tools
conda install -y -c bioconda sra-tools
数据下载
在SRA数据库中检索后勾选需要下载的数据并导出为文件
将获取的文件上传至服务器,并启用sra-tools
中的prefetch
下载
#顺序下载
prefetch --option-file SraAccList.txt
下载后会在目录下得到包含sra文件的文件夹
1.多线程下载。
cat SRR_Acc_List.txt | parallel -j 10 "prefetch {} -X 5000000000000"
2.这里还可以使用ascp进行下载,速度快,但个人使用体会不如prefetch稳定,反正是丢到服务器,不如图个稳定安心睡觉。
将所有下载的SRA文件转移到一个目录下方便后续批量分析
mkdir raw
mv SRR*/*.sra raw/
cd raw
这样我们就将所有的sra文件汇总到了raw文件夹下。接下来需要对sra文件转换为fastq。
SRA批量转换为fastq
首先在raw目录下写一个批量脚本,我这里命名为sra2fastq.sh
脚本内容:
for i in *sra
do
echo $i
pfastq-dump --gzip --split-files -t 32 $i
done
由脚本可见,这里我使用了pfastq-dump以实现多线程运算。具体安装方法:shell下运行
git clone https://github.com/inutano/pfastq-dump
cd pfastq-dump
chmod a+x bin/pfastq-dump
##将pfastq-dump添加到环境
sudo vim ~/.bashrc
#最后面添加一行,注意修改你保存的绝对路径
export PATH=$PATH:/home/xuran/Desktop/10X/pfastq/pfastq-dump/bin:$PATH
#保存后刷新bash环境
source ~/.bashrc
为脚本添加权限并运行
chmod 777 sra2fastq.sh
bash sra2fastq.sh
运行完毕后我们就讲sra文件转换为了gz压缩后的fastq文件。
这里存在三种情况:
- 从sra拆分的fastq文件只有一个:单端测序
- 从sra拆分的fastq文件有两个:双端测序
- 从sra拆分的fastq文件有三个:双端测序read+index
单细胞测序都为第2,3种情况
详见下图说明
index文件不是必须的,因此不必纠结分割出来的是两个还是三个,只需要下面分析注意一下就可以。
i7 sample index是加到Illumina测序接头上的,保证多个测序文库可以在同一个flow-cell上或者同一个lane上进行混合测序(multiplexed)。当然可以自己指定index,但更多情况下会使用10X公司提供的index序列(bundled index sets),针对不同项目使用的index也是不同的。不过共性就是:96孔板的每个孔中都加入了4种不同的index oligos混合(详见:https://kb.10xgenomics.com/hc/en-us/articles/218168503-What-oligos-are-in-my-sample-index-)。
还可以对fastq文件进行质控,对于单样本的多个run,可以舍弃其中质量较差的run,这里省略,个人认为意义不大。
2.cellranger流程
为了在下游分析中让cellranger指定识别我们的fastq文件进行下游分析,使用官网推荐的命名格式进行命名
写一个批量重命名为cellranger推荐格式的脚本,这里我保存为rename.sh
。脚本内容如下:
cat SraAccList.txt | while read i ;do (mv ${i}_1*.gz ${i}_S1_L001_R1_001.fastq.gz;mv ${i}_2*.gz ${i}_S1_L001_R2_001.fastq.gz);done
由于我这里作者上传的只有两边Read的fastq,不包含Index,因此我只批量命名了两个Read文件名称,如果你得项目包含了index可根据你的需求结合cellranger官方建议进行改写。
chmod 777 rename.sh
bash rename.sh
cellranger安装及参考基因组下载
- 下载
#参考基因组下载(human)
curl -O https://cf.10xgenomics.com/supp/cell-exp/refdata-gex-GRCh38-2020-A.tar.gz
#cellranger下载
curl -o cellranger-6.1.2.tar.gz "https://cf.10xgenomics.com/releases/cell-exp/cellranger-6.1.2.tar.gz?Expires=1649559567&Policy=eyJTdGF0ZW1lbnQiOlt7IlJlc291cmNlIjoiaHR0cHM6Ly9jZi4xMHhnZW5vbWljcy5jb20vcmVsZWFzZXMvY2VsbC1leHAvY2VsbHJhbmdlci02LjEuMi50YXIuZ3oiLCJDb25kaXRpb24iOnsiRGF0ZUxlc3NUaGFuIjp7IkFXUzpFcG9jaFRpbWUiOjE2NDk1NTk1Njd9fX1dfQ__&Signature=Zwk8l6ickwwGCQCAcsHziaoVtnS3MI0yjZMiaqV8UiyL6rW4vTjpyMHFLl04KDWpKpuMi~6D6RfHrJlZKTli---KBfC05b5u8mVcA28uDEZuSyHiOnEpI-cDdwtJ5SuGyjuNbvsUBxTCyJ~mMkrNxTivGC5XCpx6dj312qL4d4RImwEWmwMl0Nm7L8OcRGVLlujgODH51bwB03LCq1VoYIE-ECu7IhVCHlkXRzG9jLpnP98b4xf5x50WitToM4BcsW3kuQBk6w8AXjNLk3zLJwLbqKPxjAhWGZofYIKV4p-3C6G0wuJm-6XRYE0q25KrTYuDVm6Hhr~Yk46XM7Z1MA__&Key-Pair-Id=APKAI7S6A5RYOXBWRPDA"
- 解压两个文件
tar -xzvf filename
- 把cellranger添加到系统环境。
sudo vim ~/.bashrc
将下方这句添加到.bashrc文件中,注意修改目录
export PATH=~/download/cellranger-6.1.2:$PATH
- 更新一下.bashrc文件。
source ~/.bashrc
- 测试安装是否正确
cellranger testrun --id=tiny
如果运行成功会显示Pipestance completed successfully!
接下来 就进入到了cellranger 定量的正式环节,为批量运行先在目录下写一个运行脚本
在目录下编写一个cellranger运行脚本,我保存为run-cellranger.sh
。内容如下:
db=/home/xuran/Desktop/10X/cellranger/refdata-gex-GRCh38-2020-A
ls $db
fq_dir=/home/xuran/Desktop/10X/raw
cellranger count --id=$1 \
--localcores=32 \
--transcriptome=$db \
--fastqs=$fq_dir \
--sample=$1 \
--nosecondary \
--expect-cells=5000
注意修改
db后面为你所下载的参考基因组目录。
fq_dir后面为你得原始fastq地址
--localcores为指定最大使用线程数
--nosecondary为不进行聚类分群分析,因为后续使用seurat分析
--expect-cells为指定最大细胞数,根据项目决定
批量并行run-cellranger.sh进行比对定量
cat SraAccList.txt |while read id;do (nohup bash run-cellranger.sh $id 1>log-$id.txt 2>&1 & );done
批量顺序运行run-cellranger.sh进行比对定量
写一个脚本命名为muti.sh
cat SRR_Acc_List.txt |while read id;do ( bash run-cellranger.sh $id );done
后台执行脚本
nohup bash mutil.sh
SraAccList.txt
为最一开始从sra网站上保存的访问号文件。
至此,单细胞转录组从上游原始数据下载到cellranger定量的流程就结束了。
结果查看
分样本的文件夹和每个样本的运行日志,我们进入其中一个去查看。
其中最重要的就是outs文件夹
其中web_summary.html
为本次分析的质量情况
filtered_feature_bc_matrix
目录下则保存了用于下游Seurat分析的三个文件
参考来源:
https://cloud.tencent.com/developer/article/1949641
https://support.10xgenomics.com/single-cell-gene-expression/software/downloads/latest?
鸣谢:
I thank Dr.Jianming Zeng(University of Macau), and all the members of his bioinformatics team, biotrainee, for generously sharing their experience and codes.
问题交流:
Email: xuran@hrbmu.edu.cn