ixxmu / mp_duty

抓取网络文章到github issues保存
https://archives.duty-machine.now.sh/
110 stars 30 forks source link

单细胞数据在ENA数据库和NCBI的SRA是有区别的 #2581

Closed ixxmu closed 2 years ago

ixxmu commented 2 years ago

https://mp.weixin.qq.com/s/PXJ4hyzqoth1nDUWKc028g

github-actions[bot] commented 2 years ago

单细胞数据在ENA数据库和NCBI的SRA是有区别的 by 生信菜鸟团

前些天组建了《 单细胞多组学上下游全打通 》的微信交流群,是关于  clear cell renal cell carcinoma (ccRCC)  的肿瘤微环境的单细胞多组学,包括:

  • single-cell RNA sequencing (scRNA-seq)
  • single-cell assay for transposase-accessible chromatin sequencing (scATAC-seq)

数据在PRJNA768891,我优先选择了去ENA数据库下载它,因为sratoolkit的prefetch命令下载sra文件速度太慢,可以参考:使用ebi数据库直接下载fastq测序数据  , 需要自行配置好conda环境,然后去EBI里面搜索到的 fq.txt 路径文件。

首先使用conda安装download小环境,并且配置好aspera

conda create -n download 
conda activate download 
conda install -y -c hcc aspera-cli
conda install -y -c bioconda sra-tools
which ascp 
## 一定要搞清楚你的软件被conda安装在哪
ls -lh ~/miniconda3/etc/download/etc/asperaweb_id_dsa.openssh

我们已经多次介绍过conda细节了,这里就不再赘述。脚本如下:

conda activate download 
# 需要自己搭建好 download 这个 conda 的小环境哦。
# 然后构建 step1-aspera.sh  脚本文件,内容如下所示:
cat fq.txt |while read id
do
ascp -QT -l 300m -P33001  \
-i ~/miniconda3/envs/download/etc/asperaweb_id_dsa.openssh   \
era-fasp@$id  .
done
# nohup bash step1-aspera.sh 1>step1-aspera.log 2>&1 &

这个脚本会根据你在EBI里面搜索到的 fq.txt 路径文件,来批量下载fastq测序数据文件。需要熟悉数据库网页结构哦:

  • https://www.ebi.ac.uk/ena/browser/view/PRJNA768891?show=reads
  • https://www.ncbi.nlm.nih.gov/sra/?term=PRJNA768891

fq.txt 路径文件就是纯粹的文本文件,里面的需要下载的样品的aspera链接。很快就可以下载完毕,如下所示:

$ ls -lh *gz|cut -d" " -f5-
28G 8月   2 15:14 SRR16213608_S1_L001_R1_001.fastq.gz
28G 8月   2 15:18 SRR16213608_S1_L001_R2_001.fastq.gz
29G 8月   2 15:23 SRR16213609_S1_L001_R1_001.fastq.gz
29G 8月   2 15:27 SRR16213609_S1_L001_R2_001.fastq.gz
42G 8月   2 15:33 SRR16213610_S1_L001_R1_001.fastq.gz
42G 8月   2 15:40 SRR16213610_S1_L001_R2_001.fastq.gz
41G 8月   2 15:46 SRR16213611_S1_L001_R1_001.fastq.gz
34G 8月   2 15:51 SRR16213611_S1_L001_R2_001.fastq.gz
39G 8月   2 15:57 SRR16213612_S1_L001_R1_001.fastq.gz
34G 8月   2 16:02 SRR16213612_S1_L001_R2_001.fastq.gz
44G 8月   2 16:09 SRR16213613_S1_L001_R1_001.fastq.gz
37G 8月   2 16:15 SRR16213613_S1_L001_R2_001.fastq.gz
37G 8月   2 16:21 SRR16213614_S1_L001_R1_001.fastq.gz
33G 8月   2 16:26 SRR16213614_S1_L001_R2_001.fastq.gz

前面的3个样品的单细胞atac,后面的4个样品是单细胞转录组。如下所示:

前面的3个样品的单细胞atac,后面的4个样品是单细胞转录组

去跑cellranger的时候后面的4个样品成功,但是前面的3个样品是失败的。所以我去看了看为什么前面的3个样品的单细胞atac会失败:

单细胞atac会失败

在NCBI的sra数据库可以看到任意一个单细胞atac其实是4个fastq文件,但是我们前面在ENA只能看到2个fq文件的路径。

在ENA只能看到2个fq文件的路径

所以,这个项目的atac的10x的单细胞,在ena下载是没有用的,必须去ncbi的sra下载,然后自己转fq文件。官网说的很清楚:https://support.10xgenomics.com/single-cell-atac/software/pipelines/latest/using/fastq-input

每个样品的fastq文件必须是:[Sample Name]S1_L00[Lane Number][Read Type]_001.fastq.gz的格式:

Where Read Type is one of:

  • I1: Dual index i7 read (optional)
  • R1: Read 1
  • R2: Dual index i5 read
  • R3: Read 2

Alternatively, Cell Ranger ATAC will also accept ATAC FASTQs in this format:

  • I1: Dual index i7 read (optional)
  • R1: Read 1
  • I2: Dual index i5 read
  • R2: Read 2

也就是说虽然是同样的 I1文件是可以省略的,但是这个单细胞atac必须是3个fq文件,前面的ENA只能看到2个fq文件的路径,那就是不够了。

写在文末

我在《生信技能树》,《生信菜鸟团》,《单细胞天地》的大量推文教程里面共享的代码都是复制粘贴即可使用的, 有任何疑问欢迎留言讨论,也可以发邮件给我,详细描述你遇到的困难的前因后果给我,我的邮箱地址是 jmzeng1314@163.com

如果你确实觉得我的教程对你的科研课题有帮助,让你茅塞顿开,或者说你的课题大量使用我的技能,烦请日后在发表自己的成果的时候,加上一个简短的致谢,如下所示:

We thank Dr.Jianming Zeng(University of Macau), and all the members of his bioinformatics team, biotrainee, for generously sharing their experience and codes.

十年后我环游世界各地的高校以及科研院所(当然包括中国大陆)的时候,如果有这样的情谊,我会优先见你。