Hi-C衍生指标测序深度需求的系统评估：基于生物学重复变异的指标特异性阈值

《BMC Genomics》：Systematic assessment of sequencing depth requirements for Hi-C-derived metrics

背景：Hi-C（染色体构象捕获）实验产生全基因组染色质接触图谱，可据此在多个基因组尺度量化结构特征，范围从兆碱基（megabase, Mb）尺度区室（compartment）到千碱基（kilobase, kb）尺度边界（boundary）和染色质环（chromatin loop）。尽管高分辨率分析通常依赖数亿个测序读段对（read pairs），但不同类别Hi-C衍生特征所需的最低测序深度尚未被系统建立。研究人员因此试图系统确定这些深度需求。结果：研究人员对代表多种细胞类型和实验方案的十二个Hi-C文库进行渐进随机下采样（progressive random subsampling）。来自P(s)（接触概率）对数导数（log-derivative）的环密度（loop-density）和环尺寸（loop-size）推断在所有文库中评估，而区室（compartment）、绝缘（insulation）和边界分析在具有可比全深度覆盖的九个文库子集中进行。研究人员没有将“足够深度”（sufficient depth）定义为全深度值的任意比例，而是引入基于独立生物学重复（biological replicates）间指标变异性的生物学动机标准。每个指标所需的测序深度定义为下采样离散度（subsampling scatter）首次达到独立生物学重复间观察到的变异性的点，代表额外测序无法改善的准确度。使用该标准，基因组平均环密度估计在约7百万read pairs达到其生物学下限，基因组平均环尺寸估计在约20百万，绝缘分数在约80百万，100-kb分辨率区室特征向量在约110百万read pairs。结论：不同Hi-C衍生指标达到生物学有意义准确度需要显著不同的测序深度。对于所有指标，深度诱导变异性在全测序深度达到之前早已低于生物学重复变异性。这些阈值适用于基因组平均特征，而非特定位点的局部指标。因此，它们为全基因组实验设计和测序预算分配提供实用指导，提示增加生物学重复数量可能比将单个文库测序至更深更有效地改善全局生物物理特征的推断。

论文解读：Hi-C衍生指标测序深度需求的系统评估

研究背景与问题
Hi-C（染色体构象捕获）实验可将染色质接触转化为全基因组接触矩阵，并在多个尺度揭示基因组三维组织，包括兆碱基（Mb）尺度的A/B区室（A/B compartments）、亚兆碱基尺度的结构域及边界、精细尺度的染色质环，以及接触概率P(s)（contact probability，指相距基因组距离s的两位点接触频率的基因组平均）。这些特征并非以相同方式从Hi-C图谱中提取：区室注释通常来自较低分辨率观测/期望接触矩阵的特征向量；绝缘谱是沿对角线局部接触密度计算的轨道；单个环或点对应富集像素或小像素邻域，是较难重复检测的特征；P(s)则是全基因组平均，平均程度较高。因此，不同Hi-C衍生指标对测序深度的敏感性不应相同。

测序深度是Hi-C实验设计中的核心参数。现代Hi-C文库常测序至数亿read pairs，深参考图谱可含数十亿接触。虽然单碱基成本下降，但测序深度仍是随目标read pairs数线性增加的主要可变成本。文库制备有显著固定单样本成本，而更深测序使每个文库成本线性上升。因此，降低测序深度可降低测序部分成本，使更多条件、扰动或生物学重复的设计更可行。已有研究多聚焦全接触图谱可重复性、特定环调用工具或差异环检测。HiCRep、HiC-spector和GenomeDISCO等方法用于比较Hi-C图谱并校正距离依赖、平滑、谱性质或图扩散；后续系统比较显示可重复性评分本身依赖覆盖度。Oksuz et al.发现测序深度强烈影响检测环的数量和强度，甲醛（FA）单交联方案比FA+DSG双交联方案更深度敏感，区室可用数千万有效对解析，而单个环检测可能需要数亿至十亿以上reads。Parker et al.将问题视为差异环检测的统计功效，显示检测单个环变化可能需要极深图谱。另一限制是“足够深度”的定义常依赖全深度值的任意比例。该研究指出，不同指标有不同天然生物学可重复性，因为每个Hi-C文库采样的是异质染色质构象集合中的有限亚群，而不是单一确定性基因组折叠状态。因此，研究人员采用独立生物学重复在全深度下的变异作为经验性生物学地板：当随机下采样引入的散点与该生物学重复变异相当时，即认为达到该指标的准确度极限，因为额外测序无法消除独立制备细胞群体间的真实差异。该研究同时报告技术重复地板，即同一文库独立重测序之间的变异，代表最小化生物学和建库变异后的更紧可重复性。

研究人员开展的研究及总体结论
研究人员系统分析十二个公开人类Hi-C文库，来自两项研究：八个来自Oksuz et al.，涵盖两种细胞类型和两种交联方案；四个来自Rao et al.的原位Hi-C文库。研究人员将原始FASTQ read pairs从全深度下采样至200M、100M、50M、25M、10M和1M read pairs，并将每个子样本通过相同流程处理。随后量化四类指标稳定性：由P(s)对数导数短距离凹陷推断的环密度、由环尺寸相关峰推断的环尺寸、区室特征向量和绝缘谱。结论是：不同Hi-C衍生指标需要显著不同的测序深度；在全深度远未达到前，深度诱导变异性已低于生物学重复变异性；这些阈值适用于基因组平均特征，不适用于特定位点的局部指标；该结果为全基因组实验设计和测序预算分配提供指导，并提示增加生物学重复数量可能比将单个文库测序至更深更有效地改善全局生物物理特征推断。

关键技术方法（不超过250字）
十二个人Hi-C文库来自SRA：八个来自Oksuz et al.，含H1-hESC和HFFc6，采用FA或FA+DSG交联，DpnII或DdeI+DpnII消化；四个来自Rao et al.，含HMEC、IMR90、K562和KBM7。原始read pairs经distiller-nf v0.3.4流程处理至hg38，mapq30过滤，生成1 kb至10 Mb多分辨率mcool，并用cooler平衡。FASTQ用rasusa下采样至200M、100M、50M、25M、10M和1M read pairs，各深度生成独立随机子样本。P(s)用cooltools.expected_cis计算并提取凹陷和峰；区室特征向量E1在100 kb用cooltools.eigs_cis计算，采用GC含量定相，saddle图用38个分位组；绝缘分数在10 kb分辨率、100 kb菱形窗口用cooltools.insulation计算，采用Li阈值法。深度阈值通过下采样变异与全深度生物学重复和技术重复地板比较，并以对数线性插值估计。

研究结果
数据集与下采样策略：研究人员分析来自两项研究的十二个Hi-C文库。Oksuz et al.的八个文库来自H1-hESC和HFFc6，采用FA或DSG+FA交联，以及DpnII或DdeI+DpnII消化；其命名中R表示生物学重复，T表示技术重复，因此十二个文库代表十一个独立制备。Rao et al.的四个文库来自HMEC、IMR90、K562和KBM7，均为原位Hi-C。环密度和环尺寸分析使用全部十二个文库，因为这些指标来自全基因组平均P(s)，其全深度参考为每个文库内部参考。区室、绝缘和边界分析仅保留九个文库，即八个Oksuz et al.文库加KBM7，因为这些指标需标准化到每个文库全深度值，且跨文库聚合需可比全深度，约350–390M read pairs。三个被排除的Rao细胞系全深度read counts分别为IMR90 240M、HMEC 456M、K562 591M，均在Oksuz文库范围之外。每个文库下采样至六个目标深度，各深度生成独立随机子样本：1M为25个，10M和25M为20个，50M为10个，100M为5个，200M为4个，单文库全深度作为参考。由此建立可比较的深度基准。

测序深度下黏连蛋白介导的环尺寸和环密度稳定性：P(s)是相距s位点间平均接触频率，在全基因组范围近似幂律衰减，反映间期染色质皱缩组织；在中间尺度，染色质统计主要由皱缩聚合物和黏连蛋白介导环挤出描述。P(s)对数导数中的局部偏离提供环相关染色质结构读数：约10–40 kb的局部最小值（凹陷）与环密度和有限Hi-C捕获半径相关；约100 kb的局部最大值（峰）反映特征环尺寸。研究人员先建立两个噪声地板：同一条件独立生物学重复间变异和独立技术重复间变异。全深度下，环密度生物学重复变异系数（coefficient of variation, CV）为11.2%，技术重复为0.8%；环尺寸生物学CV为6.6%，技术CV为1.0%。下采样分析显示，所有十二个数据集中，平均环密度和平均环尺寸在低至10M reads时仍处于全深度值6%以内；在1M时，环尺寸保持在9%以内，而两个文库环密度偏离达15%。测序深度主要影响子样本间散点而非均值。环密度下采样CV从200M约1.6%升至50M约4.5%，再升至1M约20.6%，在7.2M read pairs达到生物学重复地板，范围1.5–20.3M；环尺寸下采样CV在20.2M read pairs跨过生物学重复地板，范围6.5–38.0M。两者均未在测试深度达到技术重复地板。协议也影响稳健性：有效片段长度v₀^eff（effective fragment length，以基因组单位表示的Hi-C捕获半径）越小，低深度越稳定；HFFc6 DdeI+DpnII文库v₀^eff约0.7 kb，在10M read pairs散点最小，而FA交联hESC文库和KBM7的v₀^eff约8.5–13 kb，不确定性较大。该部分通过下采样、P(s)对数导数凹陷和峰提取及模型拟合得出结论：环密度较稳健，环尺寸稍更需深度，且协议捕获半径与预算耦合。

区室特征向量和强度：研究人员在100 kb分辨率、九个全深度可比文库中计算A/B区室。区室特征向量E1与全深度相关性在约50M read pairs仍保留为全文库的93%，25M降至83%，1M接近随机，符号一致率约51%。生物学重复间平均E1 Pearson r为97.7%，技术重复为98.3%，说明独立深文库几乎可如重测序般重现E1轨道。下采样E1相关性在111.5M read pairs跨过生物学重复地板，范围40.2–159.5M；技术重复地板约在140M被跨过。区室强度表现异常：当用同一子样本E1对bin排序时，区室强度随下采样先升高至25M约全深度的111%，再在1M降至约14%；当用全深度E1排序时，区室强度在每个深度均保持在全深度100–103%至1M。膨胀和崩溃源于用同一噪声数据计算E1并排序，随机出现的极端bin进入saddle角；用全深度E1排序可消除该效应。区室强度生物学CV为5.1%，技术CV为0.4%。10M聚合值掩盖文库间分歧：四个hESC文库降至约全深度84%，五个非hESC文库升至约117%，但用全深度E1排序后九个文库均在100–102%。该分歧与子样本E1降解程度相关，hESC在10M仅保留15–40% E1相关，其他为68–88%。该部分通过E1相关、saddle图和区室强度计算得出结论：区室可低深度定性检测，但稳定定量注释需更深测序，且标量强度比较须谨慎。

绝缘分数和边界调用：研究人员在10 kb分辨率和100 kb窗口计算全基因组绝缘分数和结构域边界，边界未按特定大小过滤，因此包括拓扑关联结构域（topologically associated domain, TAD）以及区室/亚区室结构域。绝缘轨道与全深度相关系数在200M为98%，50M为89%，10M为62%，1M几乎为零。绝缘相关生物学重复地板为平均Pearson r=0.92，技术重复为0.97，低于E1的0.98，反映亚兆碱基结构域内在随机性较高。下采样曲线在78.3M read pairs跨过生物学重复地板，范围31.7–124.1M；50 kb窗口阈值几乎相同，为78.4M，范围36.1–127.5M。边界数量和平均边界强度呈非单调膨胀：边界数量在25M升至全深度141%；平均边界强度从200M约101%升至25M约128%，10M约169%。1M时仅在九个文库中的四个仍调用边界，数量低于全深度0.1%，强度不可解释。该膨胀是选择效应：边界在同一噪声轨道上调用和评分，通过阈值的极小值更可能是偶然加深。在±2 bin容差下，边界召回从200M的89%降至10M约56%，在44.9M read pairs达到生物学地板75.8%，范围19.2–101.6M。容差增大使阈值变浅：±1为58.7M，±3为34.3M，±5为21.2M。该部分通过绝缘轨道相关、边界调用、召回和强度计算得出结论：绝缘和边界指标需中等深度，标量边界统计在低深度易膨胀。

测序成本影响：研究人员将深度阈值转换为预算。2×150 bp测序时，深度d百万read pairs产生约0.3d Gb原始序列。所分析文库测序为2×50 bp和2×101 bp；更长读段映射更可靠，因此以原始read pairs表示的阈值对2×150 bp可能偏保守。研究人员推荐单文库至少10M read pairs用于环密度，20M用于环尺寸，约相当于全文库深度和测序成本的3–6%。对捕获半径较宽方案，如仅FA交联或具有较v₀^eff的原位Hi-C，环尺寸建议约50M read pairs。成本估算假设NovaSeq X上2×150 bp约每Gb 3美元，则环密度约9美元，环尺寸约18美元，绝缘约72美元，区室E1约99美元；标准350M read pairs文库在2×150 bp需约105 Gb原始序列，测序试剂约315美元。这些值排除建库、索引、质控、人力和设施开销，应视为测序成本下限。该部分通过成本换算得出结论：指标特异性深度选择可大幅降低测序部分成本，支持更大规模重复或扰动设计。

讨论部分总结与结论翻译
讨论指出，自Hi-C出现以来测序深度增加数个数量级，Hi-C从揭示广泛折叠模式的全基因组实验转变为可跨尺度分析染色质组织的定量高维读数。同时，Hi-C解释日益压缩为少量生物物理读数，如黏连蛋白介导随机环参数、区室特征向量、区室强度、绝缘谱和边界位置。这些指标反映不同尺度并依赖不同平均程度。该研究核心问题是每个Hi-C衍生指标需要多少测序深度。由于文库制备固定成本显著，测序深度是主要线性可变成本，降低深度可不比例降低总成本，但能显著降低测序部分，使更多条件、扰动或生物学重复设计更可行。该研究强调没有单一通用深度需求，深度需求取决于基因组尺度、平均程度、Hi-C协议，以及该指标在独立生物学重复间的可重复性。该研究避免将足够深度定义为全深度任意比例，而采用经验性生物地板。结果显示深度层级：P(s)衍生环密度凹陷最稳健，约1.5–20M read pairs达到生物地板；环尺寸峰约6.5–38M；绝缘谱约32–124M；区室特征向量尽管代表大尺度组织，但在该标准下最苛刻，约40–160M跨过生物学重复地板。排序反映指标在下采样下降解程度与其独立全深度生物学重复可重复性的组合。均值本身不是足够深度定义：平均许多低深度子样本可重建更深实验，但实际研究者只获得一个随机reads样本，关键问题是单个低深度文库能偏离全深度估计多远，因此子样本间散点比均值更合适。环密度是全局稳态参数，由每基因组长度染色质结合黏连蛋白总数决定；环尺寸是涌现长度尺度，由黏连蛋白持续性和环挤出屏障分布及动力学塑造。协议与预算耦合：较小v₀^eff的紧密捕获使低深度P(s)推断更稳定，较宽捕获需更深测序。区室E1看似需最深，因为生物学重复E1相关约r=0.98且接近技术重复，深度诱导噪声必须很小才能与生物地板相当；区室可低深度定性检测，但稳定定量注释需更深。绝缘和边界因生物地板较低，尽管更局部且对覆盖度敏感，仍在较浅深度达到其生物地板。标量汇总统计需特别谨慎：区室强度、边界数量和平均边界强度在接触图谱恶化时仍可看似改善，这是在同一噪声数据上选择特征所致，应在相同测序深度比较或使用共同参考图谱定义区室和边界。基因组当量直觉显示，10M 2×150 bp read pairs约3 Gb原始序列，约一个单倍体基因组当量；KBM7近单倍体，K562非整倍体，相同read pairs下每基因组拷贝覆盖不同。限制包括：仅限人类细胞系bulk Hi-C并映射至人类基因组；生物地板来自至少100M read pairs深文库，重复组为2–7个，组间变异存在，IMR90最可变，HMEC因缺深独立重复未纳入；随机下采样模拟read深度损失，但不模拟文库复杂度损失，真正浅实验可能变异更大；子样本与全深度参考共享reads，使观察散点偏小、阈值略浅，区室和绝缘阈值约分别占文库30%和20%，环指标为2–6%；未直接测试Micro-C、capture Hi-C、单细胞Hi-C等；读长为2×50和2×101 bp，阈值对2×150 bp可能保守但不完全转移；阈值跨协议平均；P(s)估计器平滑和模型约束可能贡献稳健性；未测试其他指标，绝缘和边界阈值依赖调用设置。定性层级比精确数值阈值更稳健。结论翻译：基因组平均Hi-C指标的测序深度需求具有指标特异性，并取决于基因组平均程度和生物学可重复性。使用全深度重复间变异性作为经验性生物地板，研究人员发现清晰层级：全基因组P(s)衍生环密度最稳健，环尺寸推断需要稍深测序，绝缘和边界指标需要中等深度，100 kb区室特征向量因较高生物学重复可重复性需要最深测序。这些结果反对单一通用Hi-C深度标准，并提供根据下游指标、实验方案和内在生物变异选择测序深度的实用框架。它们不扩展到特定位点单个环或边界的检测，后者需要显著更深测序。