翻车了,承诺上传数据后2小时内,不管多少样本,都能给出表达量矩阵,然而却没有实现。
老师有一批转录组测序数据需要放到我们的云平台上进行定量分析,就是跑 Hisat2 + Stringtie 经典流程:
首先碰到的问题是数据量比较大。约1个T(即1000G),并且老师用的是校园网,他担心直接通过网站的 Upload Data 界面上传比较慢,想通过网盘或者快递给我们数据代为上传。
其实我们平台是支持数据断点续传的,并且为了应对高通量测序数据,我们对上传工具和网络带宽都做了优化。因此告诉老师不必担心,直接上传即可。
结果,从昨天晚上8,9点钟开始,到今天下午,不到24小时,已经上传了 400 多G,可见速度还是非常可以的。
解决了上传问题,老师想先分析一批,不料流程报错了。中断在了 BAM 文件的排序步骤。我们第一反应是内存暴了,因为这个项目的基因组比较大(有10多个G),因此加大内存重试数次,结果问题依旧。
进一步排查,才发现是因为基因组的部分染色太长,超过了 BAI 索引格式所支持的最大长度。网上搜到的信息:
The original BAM index format (BAI) has a fixed bin size system, and suffers from a maximum reference sequence length of 512Mbp making it unsuitable for large chromosomes (e.g. many plants, some marsupials).
解决办法是用 CSI 代替 BAI 索引格式。
流程中原来用到的排序命令是:
sambamba sort -t $Cpu --tmpdir `pwd` -o $Sample.final.bam $Sample_$Lib.markdup.bam
该命令排序的时候默认生成 BAI 索引,并且似乎不能生成 CSI 索引,因此只能换排序命令为:
samtools sort -t $Cpu -m 4G -o $Sample.final.bam $Sample_$Lib.markdup.bam
samtools index -@ $Cpu -c $Sample.final.bam
其中 samtools index 的 -c 参数就是指定创建 CSI 索引的。可以看其帮助文档的描述:
Usage: samtools index -M [-bc] [-m INT] <in1.bam> <in2.bam>...
or: samtools index [-bc] [-m INT] <in.bam> [out.index]
Options:
-b, --bai Generate BAI-format index for BAM files [default]
-c, --csi Generate CSI-format index for BAM files
-m, --min-shift INT Set minimum interval size for CSI indices to 2^INT [14]
-M Interpret all filename arguments as files to be indexed
-o, --output FILE Write index to FILE [alternative to <out.index> in args]
-@, --threads INT Sets the number of threads [none]
经过上述修改后,我们的经典流程变成了:
fastp,进行质量控制,去除测序数据中不合格的序列,保留高质量序列用于后续分析;
hisat2,将经过质控得到的 Clean data,比对到参考基因组上;
samblaster,去除PCR或光学重复;
samtools,去除不合格的比对序列,将 SAM 文件转换成 BAM 文件,并对其进行排序;
stringtie,结合 hisat2 得到的比对文件(BAM)和基因组注释文件(GTF),进行定量分析,得到每个样本的表达量文件(FPKM / TPM);
prepDE.py,将 FPTM / TPM 转化为 Counts 计数的表达量矩阵,用于后续差异基因等分析。
大型基因组的转录组数据分析,请原谅我,大意了。
关于简说基因
生信平台
Galaxy中国(UseGalaxy.cn)致力于打造中国人的云上生物信息基础设施。大量在线工具免费使用。无需安装,用完即走。活跃的用户社区,随时交流使用心得。
生信分析
我们能够承接所有 NGS 组学数据分析业务,包括但不限于 WGS / WES / RNA-seq 等。基因组组装、注释,以及各种重测序业务都可以与简说基因合作。
生信培训
简说基因的生信培训班,荣获学员的一致好评。如果你也对生物信息学感兴趣,欢迎来跟简说基因,学真生信。
联系方式
QQ交流群(免费):925694514
微信交流群(免费):加微信好友,邀请入群
客服微信:usegalaxy