使用aspera加速从中国的GSA数据库下载测序文件

学术   2024-11-24 09:33   广东  

中国的GSA(国家基因组科学数据中心)数据库,即国家基因组科学数据中心(China National GeneBank Database),是一个存储和共享基因组数据的国家级平台。该数据库收录了大量的基因组序列数据,包括但不限于基因组组装、转录组数据、表观遗传学数据等。

值得注意的是人类的数据跟其它物种在的GSA(国家基因组科学数据中心)数据库的存储有不一样的规则:

首先如果是小鼠测序数据

比如《食管癌病人的PDX小鼠模型6个样品之CRA010501》:

  • https://ngdc.cncb.ac.cn/bioproject/browse/PRJCA016013
  • https://ngdc.cncb.ac.cn/gsa/browse/CRA010501

再比如《小鼠耳蜗单细胞转录组数据》:

  • https://ngdc.cncb.ac.cn/bioproject/browse/PRJCA006213
  • https://ngdc.cncb.ac.cn/gsa/browse/CRA004814

可以直接安装和使用这个edgeturbo软件即可:

mkdir -p ~/biosoft/
cd ~/biosoft/
wget https://ngdc.cncb.ac.cn/ettrans/download/edgeturbo-client.linux.latest.cncb.tar.gz
tar -zxvf edgeturbo-client.linux.latest.cncb.tar.gz
 ~/biosoft/edgeturbo-client/edgeturbo --help 

 ~/biosoft/edgeturbo-client/edgeturbo  dl /gsa2/CRA010501/
 ~/biosoft/edgeturbo-client/edgeturbo  dl /gsa2/CRA004814/

使用edgeturbo软件,这个速度本来就很快,无需加速。

如果是人类测序数据会存放在gsa-human

比如:

  • https://ngdc.cncb.ac.cn/bioproject/browse/PRJCA012607
  • https://ngdc.cncb.ac.cn/gsa-human/browse/HRA003340

每个项目都有一个 Excel文件可以下载,里面蕴含全部的样品的测序信息,包括fq文件地址:

全部的样品的测序信息,包括fq文件地址

打开Excel可以看到类似的ftp路径地址:

head fq.txt 
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798195/HRR798195_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798194/HRR798194_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798193/HRR798193_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798192/HRR798192_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798191/HRR798191_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798190/HRR798190_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798189/HRR798189_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798188/HRR798188_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798203/HRR798203_f1.fastq.gz
ftp://download.big.ac.cn/gsa-human/HRA003340/HRR798202/HRR798202_f1.fastq.gz

这个时候已经是可以简单的wget或者curl命令行下载,或者axel加速,但是都不够。我们可以使用aspera加速从中国的GSA数据库下载测序文件。因为可以页面的看到Aspera命令行:  帮助   信息,如下所示 :

Aspera命令行:  帮助

可以看到这个aspera01.openssh 里面的内容如下所示(未来有可能会更新,所以还是自己下载这个aspera01.openssh 文件吧 ):

cat aspera01.openssh 
-----BEGIN RSA PRIVATE KEY-----
MIIEogIBAAKCAQEA2ZwvCa5s/iDOZdt47Z+81WiNFwY+FvMDP0zRixuiTbVeudyI
6KtHITsVxSl2gA0RDAujwbswUm3m5vt+xsZGMPsIZdaEDeq0PsgkZQngSjjKnIbw
04J0r9DDtvsgZTEK9cWQ9074mSuEo5VVUBZWWltrqEE0Mb2z4nM4G3KJw7DaAB5x
azmYQwhVq3hj4jhdFhGqWQdtjk65Ib7gGbnPD5P9LBz9xeSystGpxDDSG+2TCBrz
wCqgFtOaFEvu99iEI6rUYFTwJi7iFlmSTo5DyWqTEZ8n5xoLNNLMlJaK9vwLI0v8
GpYWZjSTDhCLNcQ+Ox2h4vAFJOqdnaAHHarGxwIDAQABAoIBAG0LuAKI9rMuK1+d
aG7jMhRbDQCxryF+62yoCQVFdsKsFWjb23uEgONlIVaonWozogxANPJa5C8aRbAP
7QqcxwW6dg3aPlhe+a3QS3SaS1vGM7nWYQ4MfH6yBbFHDIfI2K4qL6fOWgMfbJsw
tQyp9OHYmA07h3U5k8/xXvdaITD/HGuANKHz8k6jAc8LGB/xbyLaDcvvfFkpsMM9
4MGPW/y0O3fWEzTcSzM8WI17fDhlZuSGiAOI43HDhwNHoFTf+kDaePGmzRdMkmKu
O/yf0FEZbGX7Fq9nHeqOLdcDSebHTvJcaRvn/AfRX7COh6dCsB7GVFC6vHE4uMFD
4phVwWkCgYEA+0DpR9nNbHb7krWQRIsNhj0aB/9ROa4CTCYtHOE4YbYMsiCzovww
WZfd9IaB9IEaKpIkGNvWlnUg7oBSVhOgQ17/IA4pX0Vw3jll5fndXyw9kDYq9Zfi
Im3lC/JkfwWRRAh0su1BmyCpfcbISFhrZS5ezvo8waCjYqzGtskzj6UCgYEA3biR
sob+9Y/QfHIFiDtFoWjfLv+FE1OQnD7AKsKUc2fyoXv6Eis/t4AO9CApoDAXCVBf
cinsirLI+0zk2387OtBnj5be5R5/eCoebJ1P4A/yE3krRpCu+DzEW3OEaSQhb3t5
db9HDw2oJbAQEH5+QvNsXcFlk3rCVMgjQHEWMPsCgYBXIUiEwTQqxAwz8UDJdEtJ
XREU6uZ5ES38yFRmqnfJ8r7uWsbQ45HE0BSJgse4SbkQQEPCVyJQZOf5rYhgD0hJ
dL0WmbKhNkyQ0+jRWtf49DeZNh+psXUbKW7/uJw/LGgW5rPJhNt2d1ovouq7o+YC
XZVFQ6QkJZfjVTVIF/gIfQKBgCZdZXiKu1sqQ10FLcfg2a+QQe4T8KbjcsZWZVIA
0UcW2XjvxtmV+jR6SBwd5JX/PD00Vw+eCXwGa5hwOblxToJixUinRnJG0K+uOg15
OvT+TVjpQn+3UU9K1H8ugd7fjZmIt/+T1WvZZRsAWAdCm/5huIKQkE7wkuewqcjg
yII5AoGAW+UYcTuWd1aTP6hObznA5Urmf4RGNEOieOgvuxxFx3m4Z7jx6J6hLWzG
hgFVibMsgcRPM8ok8hzmiQ+u0v7Dol3kdxLBz/PwDSX8KIg9w5gZAbXZyxd2QB5V
Z2rKkpsDFrIf/iIYJUoAqhAiuO28u8v5toC2WX1qFHUxeaZ+8d8=
-----END RSA PRIVATE KEY-----

很简单的构建命令行测试一下单独的下载一个文件看看 :

ascp -P33001 -i aspera01.openssh  -QT -l100m -k1 -d  aspera01@download.cncb.ac.cn:gsa-human/HRA003340/HRR798256/HRR798256_r2.fastq.gz  ./ 

### HRR798256_r2.fastq.gz       61.6Mb/s                            

所以写一个批量函数:

sed 's|ftp://download.big.ac.cn/|aspera01@download.cncb.ac.cn:|g' fq.txt > tmp
mv tmp fq.txt 
cat fq.txt |while read id
do
ascp -QT -l 300m -P33001 -k 1 -i  aspera01.openssh $id  .
done
# mamba activate download 
# nohup bash step1-aspera.sh 1>step1-aspera.log 2>&1 &
# which ascp 

当然了,这样的超高速下载如果是单个样品通常是没什么意外,但是如果要下载成百上千个文件,很容易因为网络波动导致其中的某几个失败,就需要重新针对失败的样品进行独立的下载。

可以看到的是每个10x技术的单细胞转录组样品其实会有很多fq文件,如下所示

$ ls -lh raw/*f1*gz |cut -d" " -f5-
6.9G 11月 21 01:56 raw/HRR798188_f1.fastq.gz
6.5G 11月 21 00:33 raw/HRR798189_f1.fastq.gz
5.6G 11月 20 23:36 raw/HRR798190_f1.fastq.gz
6.3G 11月 20 23:28 raw/HRR798191_f1.fastq.gz
141M 11月 20 23:10 raw/HRR798192_f1.fastq.gz
137M 11月 20 23:10 raw/HRR798193_f1.fastq.gz
115M 11月 20 23:10 raw/HRR798194_f1.fastq.gz
128M 11月 22 19:00 raw/HRR798195_f1.fastq.gz

5.8G 11月 21 05:06 raw/HRR798196_f1.fastq.gz
5.6G 11月 21 04:31 raw/HRR798197_f1.fastq.gz
6.7G 11月 23 01:03 raw/HRR798198_f1.fastq.gz
6.1G 11月 21 03:16 raw/HRR798199_f1.fastq.gz
401M 11月 21 02:05 raw/HRR798200_f1.fastq.gz
391M 11月 21 02:02 raw/HRR798201_f1.fastq.gz
455M 11月 22 19:02 raw/HRR798202_f1.fastq.gz
416M 11月 21 01:59 raw/HRR798203_f1.fastq.gz

需要自己去找到对应的关系,然后进行文件名修改 :

HRX584075       HRR798195_f1.fastq.gz
HRX584075       HRR798194_f1.fastq.gz
HRX584075       HRR798193_f1.fastq.gz
HRX584075       HRR798192_f1.fastq.gz
HRX584075       HRR798191_f1.fastq.gz
HRX584075       HRR798190_f1.fastq.gz
HRX584075       HRR798189_f1.fastq.gz
HRX584075       HRR798188_f1.fastq.gz
HRX584076       HRR798203_f1.fastq.gz
HRX584076       HRR798202_f1.fastq.gz
HRX584076       HRR798201_f1.fastq.gz
HRX584076       HRR798200_f1.fastq.gz
HRX584076       HRR798199_f1.fastq.gz
HRX584076       HRR798198_f1.fastq.gz
HRX584076       HRR798197_f1.fastq.gz
HRX584076       HRR798196_f1.fastq.gz

首先呢,上面的这些fastq文件名字是需要改名的。。。。如果你熟悉10x单细胞转录组数据,就知道:

  • 首先,1-26个cycle就是测序得到了26个碱基,先是16个Barcode碱基,然后是10个UMI碱基;通常是R1文件
  • 然后,27-34这8个cycle得到了8个碱基,就是i7的sample index;通常是I1文件
  • 最后35-132个cycle得到了98个碱基,就是转录本reads(目前很多测序仪都是150bp了),通常是R2文件

也就是说R2 文件是真正的测序reads,肯定是文件最大。而且I1文件是可以省略的。。。

眼尖的小伙伴们肯定是发现了下面的r1和r2文件其实没什么文件大小差异,是因为虽然是r1是可以碱基数量很少但是如果作者并没有处理它其实它仍然是测序仪设置的100bp或者150bp这样的初始化长度,所以就跟r2文件是一样的大小了。

128M 11月 22 19:00 raw/HRX584075_S1_L001_R1_001.fastq.gz
115M 11月 20 23:10 raw/HRX584075_S2_L001_R1_001.fastq.gz
137M 11月 20 23:10 raw/HRX584075_S3_L001_R1_001.fastq.gz
141M 11月 20 23:10 raw/HRX584075_S4_L001_R1_001.fastq.gz
6.3G 11月 20 23:28 raw/HRX584075_S5_L001_R1_001.fastq.gz
5.6G 11月 20 23:36 raw/HRX584075_S6_L001_R1_001.fastq.gz
6.5G 11月 21 00:33 raw/HRX584075_S7_L001_R1_001.fastq.gz
6.9G 11月 21 01:56 raw/HRX584075_S8_L001_R1_001.fastq.gz
416M 11月 21 01:59 raw/HRX584076_S1_L001_R1_001.fastq.gz
455M 11月 22 19:02 raw/HRX584076_S2_L001_R1_001.fastq.gz
391M 11月 21 02:02 raw/HRX584076_S3_L001_R1_001.fastq.gz
401M 11月 21 02:05 raw/HRX584076_S4_L001_R1_001.fastq.gz
6.1G 11月 21 03:16 raw/HRX584076_S5_L001_R1_001.fastq.gz
6.7G 11月 23 01:03 raw/HRX584076_S6_L001_R1_001.fastq.gz
5.6G 11月 21 04:31 raw/HRX584076_S7_L001_R1_001.fastq.gz
5.8G 11月 21 05:06 raw/HRX584076_S8_L001_R1_001.fastq.gz

这样的名字就是符合cellranger定量程序的规则啦, 接下来就完完全全参考 小鼠的5个样品的10x技术单细胞转录组上游定量(文末赠送全套代码),走cellranger流程即可。正常走cellranger的定量流程即可,代码我已经是多次分享了。参考:

差不多几个小时就可以完成全部的样品的cellranger的定量流程。

学徒作业

完成上面的HRA003340项目的单细胞转录组测序fq文件的下载以及cellranger的定量流程。


生信技能树
生物信息学学习资料分析,常见数据格式及公共数据库资料分享。常见分析软件及流程,基因检测及癌症相关动态。
 最新文章