跑一个肝癌的单细胞转录组10x数据定量流程

前面的教程里面 能从源头解决数据分析的瑕疵吗 ,我们重现了普通单细胞转录组数据分析的从fastq文件开始的走cellranger的定量流程。接下来,继续,应粉丝要求,跑一个肝癌的单细胞转录组10x数据定量流程!

https://www.ebi.ac.uk/ena/browser/view/PRJNA793914 可以看到这个项目详情,而且前些天我们演示了如何下载这个项目的fastq格式的测序数据原始文件,详见:aspera的高速下载确实很快吗。但是从网络下载的单细胞转录组数据文件的样品名字别抹掉了,变成了顺序编号的id,而且呢,文件名字并不符合规则:

7.1G 2月 23 02:24 SRR17418283_1.fastq.gz 
 17G 2月 23 02:25 SRR17418283_2.fastq.gz
13G 2月 23 02:32 SRR17418284_1.fastq.gz
13G 2月 23 02:39 SRR17418284_2.fastq.gz
16G 2月 23 03:05 SRR17418295_1.fastq.gz
17G 2月 23 03:13 SRR17418295_2.fastq.gz
11G 2月 23 03:19 SRR17418296_1.fastq.gz
12G 2月 23 03:25 SRR17418296_2.fastq.gz
5.0G 2月 23 03:28 SRR17418297_1.fastq.gz
 12G 2月 23 03:34 SRR17418297_2.fastq.gz
12G 2月 23 03:40 SRR17418298_1.fastq.gz
12G 2月 23 03:46 SRR17418298_2.fastq.gz
9.4G 2月 23 03:51 SRR17418299_1.fastq.gz
9.7G 2月 23 03:55 SRR17418299_2.fastq.gz
18G 2月 23 04:04 SRR17418300_1.fastq.gz
18G 2月 23 04:13 SRR17418300_2.fastq.gz
5.8G 2月 23 04:16 SRR17418301_1.fastq.gz
 13G 2月 23 04:23 SRR17418301_2.fastq.gz
4.9G 2月 23 04:26 SRR17418302_1.fastq.gz
 11G 2月 23 04:32 SRR17418302_2.fastq.gz

所以是需要看网页里面的样品信息详情,然后构建一个 id.txt 的文本文件,内容如下所示:

SRR17418283,PT1
SRR17418284,PT2
SRR17418295,PT3
SRR17418296,PT4
SRR17418297,PT5
SRR17418298,TT1
SRR17418299,TT2
SRR17418300,TT3
SRR17418301,TT4
SRR17418302,TT5

然后借助于chatGPT我简单的写一个脚本批量改名,这个shell脚本会读取上面的构建好的 id.txt 的文本文件:

while IFS=',' read -r old_name new_name || [ -n "$old_name" ]; do
 if [ -z "$old_name" ] || [ -z "$new_name" ]; then
 continue
 fi
 mv "${old_name}_1.fastq.gz" "${new_name}_S1_L001_R1_001.fastq.gz"
 mv "${old_name}_2.fastq.gz" "${new_name}_S1_L001_R2_001.fastq.gz"
 echo "Renamed $old_name files to $new_name"
done < id.txt

最后拿到了如下所示的文件名字:

$ ls -lh *gz|cut -d" " -f 5-

4.9G 2月 24 09:19 PT1_S1_L001_R1_001.fastq.gz
 11G 2月 24 09:32 PT1_S1_L001_R2_001.fastq.gz
 13G 2月 21 03:54 PT2_S1_L001_R1_001.fastq.gz
 13G 2月 21 04:22 PT2_S1_L001_R2_001.fastq.gz
 16G 2月 21 05:21 PT3_S1_L001_R1_001.fastq.gz
 17G 2月 21 05:41 PT3_S1_L001_R2_001.fastq.gz
 11G 2月 21 05:56 PT4_S1_L001_R1_001.fastq.gz
 12G 2月 21 06:10 PT4_S1_L001_R2_001.fastq.gz
5.0G 2月 21 00:08 PT5_S1_L001_R1_001.fastq.gz
 12G 2月 21 00:45 PT5_S1_L001_R2_001.fastq.gz
 12G 2月 21 06:25 TT1_S1_L001_R1_001.fastq.gz
 12G 2月 21 06:40 TT1_S1_L001_R2_001.fastq.gz
9.4G 2月 21 06:53 TT2_S1_L001_R1_001.fastq.gz
9.7G 2月 21 07:02 TT2_S1_L001_R2_001.fastq.gz
 18G 2月 21 01:45 TT3_S1_L001_R1_001.fastq.gz
 18G 2月 21 02:18 TT3_S1_L001_R2_001.fastq.gz
5.8G 2月 21 02:29 TT4_S1_L001_R1_001.fastq.gz
 13G 2月 21 02:55 TT4_S1_L001_R2_001.fastq.gz
4.9G 2月 21 03:06 TT5_S1_L001_R1_001.fastq.gz
 11G 2月 21 03:28 TT5_S1_L001_R2_001.fastq.gz

也就是说,每个样品是2个文件,可以看到,有一些样品它的R1文件是远小于R2文件的,但是有一些样品的1文件和R2文件居然是大小一致的。这些文件名字就是符合要求的了,正常走cellranger的定量流程即可,代码我已经是多次分享了。参考:

差不多几个小时就可以完成全部的样品的cellranger的定量流程。基础知识非常重要,我们在单细胞天地多次分享过cellranger流程的笔记(2019年5月),大家可以自行前往学习,如下:

接下来我们就常规的跑一下它:

bin=/home/jmzeng/x10/pipeline/cellranger-6.1.2/bin/cellranger
# db=/home/bakdata/x10/pipeline/refdata-gex-mm10-2020-A 
db=/home/jmzeng/x10/pipeline/refdata-gex-GRCh38-2020-A
ls $bin; ls $db
fq_dir=/home/jmzeng/x10/jmzeng/hcc_patients_PRJNA793914/asp/raw 
$bin count --id=$1 \
--localcores=4 \
--transcriptome=$db \
--fastqs=$fq_dir \
--sample=$1 \
--expect-cells=5000

有意思的是,上面的10个样品我同一时间(2024-02-24 09:33:21)开始运行cellranger的定量流程,但是大家的耗时其实是不一样的:

log.cellranger.0.txt:2024-02-24 13:12:44 Shutting down.
log.cellranger.1.txt:2024-02-24 18:19:43 Shutting down.
log.cellranger.2.txt:2024-02-24 19:43:29 Shutting down.
log.cellranger.3.txt:2024-02-24 17:07:41 Shutting down.
log.cellranger.4.txt:2024-02-24 13:20:12 Shutting down.
log.cellranger.5.txt:2024-02-24 19:56:38 Shutting down.
log.cellranger.6.txt:2024-02-24 16:24:52 Shutting down.
log.cellranger.7.txt:2024-02-24 19:09:46 Shutting down.
log.cellranger.8.txt:2024-02-24 13:40:30 Shutting down.
log.cellranger.9.txt:2024-02-24 13:04:43 Shutting down.

也就是说有4个样品是三个半小时就跑完了cellranger的定量流程,但是另外的6个样品基本上耗时10个小时。区别在于有一些样品它的R1文件是远小于R2文件的,但是有一些样品的1文件和R2文件是大小一致的,是因为这些样品的fastq文件并没有在r1文件被trim,我们可以看2022的一个文章的示意图(The cellBC and UMI assigned Nanopore reads (median 28,120 /cell) reflect a median of 2427 genes and 6047 UMIs per cell with a good )详见:https://www.nature.com/articles/s41467-020-17800-6:

官网示意图

Comments are closed.