单细胞测序这两年最不缺的就是公开数据GEO 上随便检索一个关键词动辄就是几十个 GSE 数据集。但真正让刚入门的人卡住的往往不是分析本身而是第一步把 NCBI 上的 10X 原始数据拿下来整理成 Cell Ranger 能认的格式。我自己第一次做这件事的时候光是这个 SRR 号下下来怎么有三个 fastq 文件为什么 --sample 写的名字和文件名对不上就报 No fastq files found就折腾了整整一个下午。这篇文章就把从 NCBI 找数据、判断数据形态、用 SRA Toolkit 转 FASTQ、必要时用 bcl2fastq 拆 BCL一直到 cellranger count 跑通、看 web_summary 判断数据质量这一整条链路讲清楚。适合完全没碰过 SRA 数据库的新手也适合跑过几次但总在格式问题上翻车的同学。1. 先搞清楚数据在 NCBI 的哪个角落1.1 GEO 页面上三个入口的区别很多人拿到一篇单细胞文章第一反应是打开 GEO 页面然后看到一大堆文件就懵了。其实 GEO 页面上跟能不能跑 Cell Ranger相关的入口就那么几个先分清楚它们的性质后面能省掉大量无用下载。打开一个 GSE 页面往下拉你会看到Supplementary file这一栏。这里面的东西分两种一种是GSE123456_RAW.tar这种打包的原始文件里面可能是每个样本的 BCL 压缩包也可能是已经切好的 FASTQ另一种是形如GSMxxxxxx_sample_matrix.mtx.gz、barcodes.tsv.gz、features.tsv.gz的表达矩阵文件。前者是原始数据后者是已经处理过的结果。再往下或者往页面顶部的链接区会有一个SRA Run Selector或者BioProject的入口。点进 SRA Run Selector你会看到一个表格每一行是一个 SRR 号对应一次测序运行。这才是最干净的原始数据来源。区分方法很简单矩阵文件只能做下游分析跑不了 Cell Ranger。Cell Ranger 需要的是原始 reads也就是 FASTQ 或 BCL。你要是拿 mtx 去喂 Cell Ranger它会直接告诉你找不到 fastq 文件。所以第一步判断这个数据集有没有提供 raw data。很多 2020 年以后的文章会同时提供两者也有不少老数据集只给了矩阵。1.2 从 GEO 跳到 SRA Run Selector 的正确姿势GEO 和 SRA 的关系是样本元数据和测序原始记录的关系。同一个项目在 GEO 叫 GSE 号在 SRA 叫 SRP 号样本层面 GEO 是 GSMSRA 是 SRS一次测序运行是 SRR。从 GSE 页面跳转的路径通常是这样的在页面右侧的Relations区域找到SRA链接点进去就是该项目的 SRA 页面。或者直接改 URL把https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?accGSE123456换成https://www.ncbi.nlm.nih.gov/Traces/study/?accGSE123456SRA Run Selector 会自动把该项目下的所有 run 列出来。这里有个经验优先用 SRP 号而不是 GSE 号去 Run Selector 检索。因为一个 GSE 有时候会关联多个 SRP用 GSE 检索可能漏掉一部分 run。在 Run Selector 页面顶部有个Study下拉框能看到该项目关联的所有 SRA study。Run Selector 页面上最有用的是那个可以筛选和导出的表格列包括 Run、BioSample、Experiment、LibraryLayout、Bytes、Bases、spots 等。你可以把这些行勾选后点Metadata导出成 CSV拿到本地做批量下载脚本。这一步比一个个点网页高效太多尤其是面对几十个 run 的项目。1.3 先说结论哪些项目根本跑不了 Cell Ranger踩过几次坑之后我养成了一个习惯在下载之前先花五分钟判断这个数据集到底值不值得下因为一个 10X 样本的 FASTQ 动辄 20 到 50 GB下错了纯属浪费时间和磁盘。判断标准我总结成这么几条情况能否跑 Cell Ranger说明SRA 里有 3 或 5 转录组 reads可以最常见直接 fasterq-dumpGEO 提供 RAW.tar 且内含 BCL可以需要 bcl2fastq 或 mkfastq只提供 filtered_feature_bc_matrix不可以已有细胞过滤无法回溯只提供 raw_feature_bc_matrix不可以只有计数矩阵没有 reads提供了 FASTQ 但被 trim 过视情况若 barcode 端被截短则无法用空间转录组 Visium可以但流程不同需要 spaceranger 不是 cellranger单细胞 ATAC可以但流程不同需要 cellranger-atac还有一个细节有些项目做的是单细胞核测序snRNA-seq这个用 Cell Ranger 跑完全没问题参考基因组和参数都不变只是在解读 web_summary 的时候要注意细胞数会偏低、线粒体基因比例会明显更低。2. 拿到 SRR 号之后先别急着下判断数据形态2.1 RAW.tar、SRA Run、现成 FASTQ 三种形态NCBI 上的 10X 数据落到你手里通常有三种形态每种对应的后续操作完全不同。第一种是GEO 的 RAW.tar。解压之后一般是每个样本一个GSMxxxxxx.tar再解压里面是bcl/目录或者fastq/目录。如果是 BCL你需要 bcl2fastq 或者 cellranger mkfastq如果是 FASTQ直接跳到命名整理那一步。第二种是SRA Run。这是最普遍的形态一个 SRR 号对应一次 sequencing run用 SRA Toolkit 转成 FASTQ。要注意的是SRA 里的 run 粒度不一定等于样本粒度——有的项目一个样本拆成多个 run多 lane有的项目多个样本混在一个 run 里用 index 区分。第三种是作者直接提供的 FASTQ 附件。这种情况在近两年的数据集里越来越常见因为测序成本的下降让作者更愿意直接上传大文件。这类文件通常已经在压缩包里按样本分好目录但文件名往往被改得乱七八糟需要重新整理。三种形态的判断方法下载前先读 GEO 页面上Supplementary file的文件名和大小。如果看到*_RAW.tar且几 GB 以上多半是 BCL 或 FASTQ如果看到*_fastq.tar.gz之类的那就是现成的。2.2 读 Run Browser 里的读长和 Layout在 SRA Run Selector 里点任意一个 SRR 号会进到 SRA Run Browser 页面。这个页面里有两个字段非常关键Layout和Read specification。Layout一般显示PAIRED说明是双端测序。10X 的数据基本都是双端Read1 是 cell barcode 加 UMIRead2 是 cDNA 片段。Read specification会给出每一端的长度比如28,91,8,91或者151,151。这个数字序列的解读方式是如果只有两段比如28,91那 Read1 是 28bp 的 barcodeUMIRead2 是 91bp 的 cDNA。如果是四段前两段是 Read1 和 Read2后两段是 Index Read 的两端。这个信息决定了你后面拿到 FASTQ 之后哪个文件对应 R1、哪个对应 R2。我遇到过不止一次fasterq-dump 出来的_1.fastq长度是 91_2.fastq是 28跟 Cell Ranger 的预期正好反了。判断方法很简单用zcat file.fastq.gz | head -2 | tail -1 | wc -c看一眼第一条 read 的长度就知道了。顺便说一句如果是 10X v3 化学R1 通常是 28bp16bp barcode 12bp UMIv2 是 26bp16 10。Read2 长度跟测序配置有关常见 91 或 98。这些数字记住之后看到 read 长度就能大致判断化学版本。2.3 一个容易被忽略的坑SRA Lite 与质量值这是近两年新出现的坑值得单独拎出来说。NCBI 从 2023 年前后开始把一部分旧的 SRA 数据转换成了 SRA Lite 格式。这种格式的特点是不存储原始的碱基质量值下载出来转成 FASTQ 之后质量值那一列是固定的占位字符。对 Cell Ranger 来说如果你的数据本身质量很好这个影响可能不大但严格来说它会丢失 Q30 之类的统计信息web_summary 里的质量相关指标会失真。判断方法是在 SRA Run Browser 页面看有没有SRA Lite的标记或者在 Run Selector 表格的列里找SRA-Lite字段。如果中招了可以尝试用prefetch --type all拉取完整数据或者去其他镜像数据库找同一个 run 的原始版本。我个人的做法是如果这是个关键样本宁可换个数据源也不要拿一份没有质量值的 FASTQ 硬跑。3. SRA Toolkit 把 SRA 变成 Cell Ranger 认的 FASTQ3.1 prefetch 与 fasterq-dump 的分工SRA Toolkit 装好之后你会看到两个名字像双胞胎的命令prefetch和fastq-dump以及新版的fasterq-dump。很多人搞不清该用哪个其实它们的分工非常清楚。prefetch负责下载。它把 SRA 文件从远端拉到本地缓存目录默认在~/ncbi/public/sra/下面。这个命令最大的好处是支持断点续传网络抖一下断了再跑一次它会从断掉的地方继续不会从头再来。fasterq-dump负责转换。它把本地的.sra文件解压、拆分成 FASTQ。之所以叫 faster是因为它比老的fastq-dump快好几倍支持多线程是现在的主流选择。最小可用的命令组合长这样# 下载 prefetch --max-size 100G SRR1234567 # 转换-e 指定线程数 fasterq-dump --split-files --threads 8 --outdir ./fastq_raw SRR1234567--split-files是必须加的。不加的话双端数据会被塞进一个文件里后面完全没法用。3.2 fasterq-dump --split-files 出来的文件到底哪个是 barcode跑完上面的命令./fastq_raw目录里会出现SRR1234567_1.fastq和SRR1234567_2.fastq有时候还有_3.fastq。这里就是最容易翻车的地方。默认情况下_1对应 Read1_2对应 Read2_3一般是 index read。但并非所有提交者都按照这个顺序提交数据尤其是早期项目或者作者自己用 bcl2fastq 转换后上传的情况。判断方法就是前面提到的看长度# 看第一条序列的长度 zcat SRR1234567_1.fastq.gz | sed -n 2p | wc -c zcat SRR1234567_2.fastq.gz | sed -n 2p | wc -c如果_1是 28 左右_2是 90 以上那就是标准情况不用动。如果反过来了你就需要在后续命名的时候手动交换把短的那个命名为 R1长的命名为 R2。这个交换操作一定要在命名阶段做掉不要去改文件内容。另一个情况是出现_3.fastq。这个是 index readI1Cell Ranger 不需要它因为它只看 FASTQ 内容不读 index 信息做拆分。你可以直接删掉以节省空间但如果这个 run 里混了多个样本那_3就是拆样本的唯一线索得留着配合--lanes或者其他工具处理。3.3 重命名成 _S1_L001_R1_001 规范Cell Ranger 对输入 FASTQ 的文件名有严格的格式要求这一点官网上写得很清楚但很多人第一次看会忽略。规范的命名格式是[Sample Name]_S1_L00[Lane Number]_[Read Type]_001.fastq.gz举个例子SampleA_S1_L001_R1_001.fastq.gz和SampleA_S1_L001_R2_001.fastq.gz。其中Sample Name必须和你命令行里--sample参数的值完全一致L001是 lane 号单 lane 数据固定写 L001 就行R1/R2是读端。所以转换完之后的关键一步就是重命名cd ./fastq_raw mv SRR1234567_1.fastq SampleA_S1_L001_R1_001.fastq mv SRR1234567_2.fastq SampleA_S1_L001_R2_001.fastq pigz -p 8 SampleA_S1_L001_R1_001.fastq pigz -p 8 SampleA_S1_L001_R2_001.fastq注意pigz这个工具它是 gzip 的多线程版本压缩 20G 的文件用它能快好几倍。原始 FASTQ 不压缩的话一个样本可能占 60 到 80 GB压完通常能到 15 到 25 GB磁盘压力的差别相当大。提示如果你的样本有多个 lane每个 lane 都要重命名成_L002_、_L003_这样的形式--sample只写样本名Cell Ranger 会自动把所有 lane 合并。命名不一致是后面 No fastq files found 报错的第一大原因。3.4 磁盘、内存与并发参数这条链路里最容易在硬件上翻车。几个实测数字供参考--threads参数控制 fasterq-dump 的并发。这个值不是越大越好因为每个线程都要占内存和磁盘 IO。在 16 核 64G 内存的机器上我给--threads 8比较稳妥32 核的机器可以给到 16。给太大反而会因为磁盘瓶颈导致整体变慢。磁盘空间上fasterq-dump在转换过程中会同时存在.sra文件、中间的临时文件和输出的 FASTQ峰值占用可能是最终 FASTQ 体积的两倍多。所以转换一个样本前最好预留 100 GB 以上的空间。我有个朋友就是在磁盘只剩 40G 的时候跑转换跑到一半把系统盘写满了最后连 log 文件都没写出来。prefetch还有一个--max-size参数默认是 20G超过这个大小的 run 会被跳过。10X 的 run 很容易超过这个值所以一定要显式给大一点比如--max-size 200G否则你会看到它下载了几秒钟就完成了其实什么都没下。3.5 断点续传与失败重试大规模下载最怕的就是跑到 90% 断掉。prefetch在这点上做得很好它的缓存机制让你可以直接重跑命令# 第一次 prefetch --max-size 200G SRR1234567 # 断了之后直接再来一次它会续传 prefetch --max-size 200G SRR1234567下载完成后建议做一次完整性校验vdb-validate ~/ncbi/public/sra/SRR1234567.sra这个命令会检查本地 SRA 文件的完整性输出ok才算真正下好了。别跳过这一步否则你可能在 fasterq-dump 跑到一半才发现数据是坏的。如果项目里有几十个 run写个循环是必然的while read srr; do prefetch --max-size 200G $srr fasterq-dump --split-files --threads 8 --outdir ./fastq_raw $srr done srr_list.txt这个循环可以放后台跑用nohup或者screen挂着第二天来看结果。4. 如果 GEO 给的是 BCL 原始数据4.1 判断手里的 RAW.tar 是不是 BCL解压GSE123456_RAW.tar之后你会看到一堆GSMxxxxxx.tar。随便挑一个解开如果目录结构里有Data/Intensities/BaseCalls/这样的路径并且有.bcl或者.bcl.gz文件还有RunInfo.xml、runParameters.xml那这就是标准的 Illumina BCL 输出目录。如果是fastq/目录下面直接是.fastq.gz那就跳回上一节的重命名流程。BCL 格式的好处是保留了全部原始信息理论上可以重新做 base calling坏处是必须用 bcl2fastq 转换多一道工序而且 bcl2fastq 的安装在某些系统上有点折腾人。4.2 bcl2fastq 与 mkfastq 的版本对应关系Cell Ranger 里的mkfastq实际上是 bcl2fastq 的一个封装它自己不带 bcl2fastq 的二进制需要你在系统里单独装好。版本对应关系大致是Cell Ranger 版本需要的 bcl2fastq 版本3.xbcl2fastq2 v2.204.xbcl2fastq2 v2.205.xbcl2fastq2 v2.206.xbcl2fastq2 v2.207.xbcl2fastq2 v2.20可以看到其实都是 v2.20这个版本从 2017 年发布之后就没怎么变过。安装方式有两种一是用官方提供的 rpm/deb 包二是从源码编译。第二种比较麻烦依赖一堆 boost 和 zlib个人建议能用包管理器就用包管理器。装完之后用bcl2fastq --version确认一下然后在 Cell Ranger 里用--bcl2fastq或者--bcl2fastq2参数指定路径。4.3 SampleSheet.csv 怎么写BCL 转换的核心是 SampleSheet.csv它告诉 bcl2fastq 怎么根据 index 把数据拆成不同样本。10X 的 SampleSheet 通常长这样[Header] IEMFileVersion,4 Investigator Name,xxx Experiment Name,xxx [Data] Lane,Sample_ID,Sample_Name,Index,Sample_Project 1,SampleA,SampleA,SI-GA-A1,Project1 1,SampleB,SampleB,SI-GA-A2,Project1关键点在Index这一列。10X 的 v2 和 v3 用的是成套的 index比如 SI-GA-A1 到 SI-GA-H12 这一组如果你的 SampleSheet 里 index 写错了转换出来的样本会全部为空。最保险的做法是直接用 GEO 附件里自带的 SampleSheet.csv不要自己重写。然后跑cellranger mkfastq --idrun1 \ --run/path/to/bcl_dir \ --csv/path/to/SampleSheet.csv \ --localcores16 --localmem64出来的结果在run1/fastq/下面目录结构自动就是样本名分子目录文件名也已经是规范格式可以直接喂给 count。4.4 参考基因组的版本陷阱这一步跟 BCL 无关但同样致命。Cell Ranger 需要一个参考基因组10X 官方提供了预构建的包比如refdata-gex-GRCh38-2020-A。这个包有版本而且和 Cell Ranger 的版本有对应关系2020-A 系列可以用在 Cell Ranger 3.0 到 7.02024-A 需要用 Cell Ranger 8.0 及以上如果你用的是老版本参考基因组配合新版本 Cell Ranger多数情况能跑但会提示警告反过来新参考配老软件则可能直接报错。还有一个常见误解GRCh38-2020-A 和 GRCh38-3.0.0 不是同一个东西。前者更新了基因注释包含了更多的 lncRNA 和更新过的基因名。用不同版本跑出来的基因数会有差异做横向比较的时候必须保证全项目用同一个版本。我的习惯是在项目开始就把参考基因组的 md5 记在 README 里避免半年后回来看不知道自己用的哪个版本。5. cellranger count 跑通与参数取舍5.1 目录结构与 --sample 的匹配逻辑假设你的 FASTQ 已经整理好了目录结构长这样fastqs/ ├── SampleA_S1_L001_R1_001.fastq.gz ├── SampleA_S1_L001_R2_001.fastq.gz ├── SampleB_S1_L001_R1_001.fastq.gz └── SampleB_S1_L001_R2_001.fastq.gz那么跑 SampleA 的命令是cellranger count --idSampleA_run \ --transcriptome/ref/refdata-gex-GRCh38-2020-A \ --fastqs./fastqs \ --sampleSampleA \ --expect-cells5000 \ --localcores16 \ --localmem64--sample的值是SampleACell Ranger 会自动去--fastqs目录里匹配所有以SampleA_开头的文件收集成 R1/R2 配对。这就是为什么命名必须严格——它匹配靠的是前缀字符串不是智能识别。一个非常容易犯的错--sample写成了文件名的完整前缀但带了_S1比如写--sampleSampleA_S1。这样它匹配不到任何文件报 No fastq files found。记住--sample只需要样本名前缀_S1_L001_R1_001这部分是格式约定不用写进去。5.2 expect-cells 和资源参数怎么给--expect-cells这个参数的作用是给 Cell Ranger 一个初始的细胞数量预期它会根据这个值来调整 barcode 的过滤策略。给得太少真实的细胞会被当成背景过滤掉给得太多背景 barcode 会被误判成细胞。经验值是这样先大致估一下样本上机时目标细胞数然后按这个数来给。如果完全不知道可以从 3000 到 5000 起步跑完看 web_summary 里的Estimated Number of Cells如果这个值跟你设的 expect-cells 差得很远就调整参数重跑一次。常见的判断标准是估计值落在 expect-cells 的 0.5 到 2 倍区间内比较合理。--localcores和--localmem控制资源。这两个参数填错会导致两个典型问题cores 给太多而系统实际核数不够任务会因为抢不到资源而卡住mem 给超过物理内存会在比对阶段被系统 OOM killer 干掉日志里只会留一行被杀的记录非常难排查。我的一般做法是localcores给到系统核数的 80%localmem给到物理内存的 80%。比如 32 核 128G 的机器就给 24 核 100G。5.3 跑完之后先看 web_summary 的哪几项指标cellranger count跑完大约需要 2 到 8 小时取决于数据量和机器性能。出来之后第一件事是打开outs/web_summary.html重点看这几项指标参考范围说明Estimated Number of Cells与预期相符和 expect-cells 差太多要警惕Mean Reads per Cell 20000低于 10000 说明测序深度不足Median Genes per Cell 1000太低可能是细胞活力差或 RNA 降解Fraction Reads in Cells 70%低于 50% 说明背景噪音大Sequencing Saturation60% - 90%过高说明测序过饱和加测无意义Q30 Bases in Barcode 85%低于 80% 要检查数据质量Valid Barcodes 75%太低可能是 barcode 端有污染这几项里我最看重的是Fraction Reads in Cells和Median Genes per Cell。前者反映样本的干净程度后者反映 RNA 的完整性。如果这两个指标都正常基本可以放心往下做。5.4 输出目录里哪些文件后面还会用到outs/目录里文件不少但真正高频使用的就那几个filtered_feature_bc_matrix/是过滤后的矩阵只包含判定为真实细胞的 barcode是下游分析的主力输入。raw_feature_bc_matrix/包含全部 barcode做背景估计或者空液滴分析时用。possorted_genome_bam.bam是排好序的 BAM 文件做 CNV 分析、可变剪接、或者想自己重新定量的时候需要它但体积很大一个样本可能 30 到 50 GB。molecule_info.hdf5是做aggr合并多样本时必须的文件很多人跑完 count 就把它删了等到要合并样本时又得重新跑非常亏。cloupe.cloupe是给 Loupe Browser 用的可视化文件做汇报的时候挺方便。6. 报错排查从日志第一行开始6.1 No fastq files found这是出现频率最高的报错原因就那么几种--fastqs路径写错了相对路径和绝对路径混用最容易出问题文件名不符合_S1_L001_R1_001规范--sample和文件名前缀不一致R1和R2标签写反。排查方法是从内向外一层层确认。先ls一下--fastqs指的那个目录确认文件在再用ls fastqs/SampleA_*确认前缀匹配得上最后确认每个样本都有成对的 R1 和 R2。注意Cell Ranger 对大小写敏感。SampleA和samplea是两个不同的东西命名的时候统一风格别一会儿大写一会儿小写。6.2 读长与chemistry不匹配报错信息里可能出现Read 1 is too short或者The read length does not match the expected chemistry。这是因为 Cell Ranger 会自动检测化学版本如果 read 长度偏离了它预期的范围就会报错。v2 化学的 R1 是 26bpv3 是 28bp两者差 2bp。大部分情况下 Cell Ranger 能自动识别但如果你之前在 fastq 处理中做过 trim把 R1 截短了它就会认不出来。解决办法是用--chemistry参数显式指定可选值包括SC3Pv2、SC3Pv3、SC3Pv3LT、SC5P-PE等。指定对了之后报错就消失了。但更好的做法是根本不要对 R1 做 trim因为 barcode 和 UMI 就在这 28bp 里截短任何一个碱基都会导致比对失败。6.3 一个 SRR 里混了多个样本有时候你会遇到一个 SRR 跑完Estimated Number of Cells 高得离谱比如预期 5000 结果出来 20000。这很可能是因为这个 run 里用 index 混了多个样本。判断方法是看 SRA Run Browser 里的Library Strategy和样本描述。如果多个 BioSample 共享一个 Run那就是混样了。处理这种数据比较麻烦需要根据 index 序列来拆分。可选方案是把_3.fastqindex read拿出来用demultiplex类工具按 index 拆分或者干脆用cellranger multi配合 feature barcode 的配置。这种情况我一般会优先去找作者是不是另