前言
想要白嫖子宫内膜癌的单细胞测序数据用来做验证,搜到一篇做的很全面的单测文章,但作者上传的是SRA文件,只好拾起3年前跟生信技能树曾老师启蒙的上游数据分析技能。
网上一顿搜索,装了好多包,结果本次白嫖前的工作准备真正用上的就3个,sra-tool,parallel-fastq-dump(基于 sra-tool,所以 sra-tool 是必须的)和star solo,嗯还有迅雷(10年没用过了我)。
这篇笔记旨在记录一些我记不住的设置(主要是star使用的g++编译,在非intel芯片上安装需要进行的环境设置),也希望能给想用mac进行上游数据分析的同道一些参考。
一、搭建测序数据处理环境
1. conda环境搭建参考
https://blog.csdn.net/weixin_47614014/article/details/118070452
2. 上游数据处理包安装参考
https://mp.weixin.qq.com/s/mo_oDD_ZWrC18EKHkLNR9g
二、参考基因组下载
http://ftp.ebi.ac.uk/pub/databases/gencode/
①genome.fa.gz是基因组文件,解压为genome.fa文件
②gtf或gff3.gz是注释文件(有全基因组注释的、外显子注释的、非编码RNA注释的,按需)
③可用复制链接使用迅雷下载
三、上游测序数据下载(本例数据量大,勿作练习用)
https://www.ncbi.nlm.nih.gov/sra/?term=SRP349751
1. 网上教程教我的
①点All runs
②下载Accession list
③使用Sra toolkit下载
2. 速度实在太慢了,我用的迅雷
①15个数据,也不是很多,链接慢慢改吧,速度喜人!
②下载 normalize 的数据,跟使用 sra toolkit 下载的一致,可以用 fast-dump 直接转,lite.sra的数据需要特殊处理
③需添加后缀.sra
四、sra文件解压为fastq文件
使用parallel-fastq-dump,因为快
1. 安装 parallel-fastq-dump
https://anaconda.org/bioconda/parallel-fastq-dump/files
下载noarch版本(python环境,兼容Mac系统)或者osx版本,因为在终端使用conda下载容易出现不兼容OSX的版本,所以我都是进官网下载,再解压安装(反正用Macbook,也是因为没租用linux服务器是吧?不需要在乎服务器可视化问题)
# 需要使用python环境的包,我都放在miniforge2/pkgs文件夹里面了,我也不知道不放里面解压有没有影响
tar -vxzf parallel-fastq-dump-0.6.3-py36_1.tar.bz2
# 配置环境
vim ~/.zshrc
# export PATH="/Users/用户名/miniforge2/pkgs/parallel-fastq-dump-0.6.3-py36_1/bin:$PATH"
source ~/.zshrc
# 查看是否安装成功
conda list
2. 使用parallel-fastq-dump解压sra文件
双端的单测数据一般都生成3个fastq文件,分别是I1(index文件),R1(细胞标签read文件)和R2(单个细胞内测序文件),后续分析主要用到R1和R2具体可参考https://www.jianshu.com/p/dadd202c34be
# 即便使用parallel-fastq-dump 8线程,20几G的sra文件也用了半个多小时
parallel-fastq-dump --sra-id SRR17165228 --threads 4 --outdir ../rawfq/ --split-files --gzip
# 如果测序文件只是几个G(细胞量少吧),可以使用批量程序挂机
ls *.sra | while read id
do ( nohup parallel-fastq-dump -t 72 -O ../rawfq/ --split-files --gzip -s $id &)
done
# 没有试fasterq-dump,虽然它也支持多线程,但不支持gzip压缩,本例数据太大,应该不够空间存放。。fastq-dump、fasterq-dump和parallel-fastq-dump区别和使用可以参考https://www.jianshu.com/p/97e7a70aaf79
五、fastq文件比对并输出Counts文件
1. 安装STAR
https://github.com/alexdobin/STAR
我遇到的主要问题包括:
①git clone速度慢→可以通过GitHub Desktop进行下载
②gcc编译→修改编译环境,参考https://blog.csdn.net/qq_33957603/article/details/131757260
③-mavx2报错→source/Makefile文件里修改为CXX_SIMD_FLAGS=-march=native
ps. 更新了2.7.11a版本似乎可以使用zip下载了,同时makefile文件里面也不存在“CXX_SIMD_FLAGS”了,应该已经解决Mac安装该软件时候-mavx2报错问题。
2. 建立比对索引
STAR --runMode genomeGenerate --genomeDir ~/Downloads/Reference/hg38 --genomeFastaFiles ~/Downloads/Reference/hg38/GRCh38.p14.genome.fa --sjdbGTFfile ~/Downloads/Reference/hg38/gencode.v44.annotation.gtf
3. 下载barcodes
STARsolo与普通的转录组比对区别在于你需要在比对时加上whitelist,参考https://www.jianshu.com/p/b2076d670558
barcodes是细胞的标签序列文件,使用10XGenomics试剂盒产生的数据,可以从10XGenomics处下载barcodes
https://github.com/10XGenomics/cellranger/tree/master/lib/python/cellranger/barcodes
本例用到的是3M-february-2018.txt
4. 运行单细胞比对
先解压fastq.gz文件解压到移动硬盘(其实上面可以用faster-dump生成fastq文件,这里就不需要解压了,主要还是硬盘容量不够,fastq缓冲文件会占用硬盘空间)
尽量删除占用电脑硬盘的文件,因为bam文件比fastq文件大几倍,本例生成的bam文件200G左右(请大佬们教我怎么将缓冲文件也保存到移动硬盘,可能可以解决这个电脑硬盘不够的问题)
STAR --runThreadN 16 \ #多线程,所以比cellrange快,但是Mac本来也用不了cellrange,不知道修改编译环境后能不能安装,未尝试
--genomeDir ../../Reference/hg38 \ #比对索引目录
--readFilesIn *_3.fastq *_2.fastq \ #见上述3为R2文件,2为R1文件,starsolo需要先比对R2文件,再通过R1文件分配到单细胞里面去
--soloType CB_UMI_Simple \ #10X的单细胞测序用这个one UMI and one Cell Barcode of xed length in read2, e.g. Drop-seq and 10X Chromium
--soloCBwhitelist ../../Reference/barcodes/3M-february-2018.txt \ #barcodes文件
--soloBarcodeReadLength 0 #默认是1,但本例报错,可能与添加了接头序列导致与barcodes文件的碱基数不匹配
starsolo还可以根据不同的比对文件使用单细胞测序数据进行多种分析,例如ncRNA、exon可变剪接等,参考文章先保存在这了https://mp.weixin.qq.com/s/XoXBq6OMkghlhsg4yurHJA
5. 生成进行下游分析的文件啦
后话
这次白嫖花了半个月时间,希望后面再白嫖需要上游分析的测序数据,可以少花些时间,毕竟代码这东西用着用着就又报错了😭😭😭