11 对 5 对 10：一个不平衡小样本设计能回答什么，不能回答什么

我们重新分析 Wu 等人的乳腺癌单细胞图谱时，中途出现过一对看起来可以发表的数字：ER+ 与 TNBC 之间，内皮细胞份额差异 q = 0.0189，髓系 q = 0.0352。两个都低于 0.05，方向在留一患者下稳定，效应量都超过一个标准差。

这两个数字没有作为发现写出来。复算时它们没有变；变的是校正家族——原先那两个 q 来自一个比实际做过的检验更小的家族。改回真实的家族边界后，它们成为 0.0568 和 0.1057，15 个检验里通过 BH 的数量是 0。

这篇文章写的就是这个过程，以及一个更基本的问题：这个设计原本能回答什么。

论文和数据来源

原论文包含空间转录组、CITE-seq 与 TCGA 反卷积分析；本文只讨论公开 scRNA-seq 矩阵能直接支撑的再分析层。

队列：26 位患者，三组，谁都不平衡

26 位患者、100,064 个细胞进入分析，按临床亚型分成三组：ER+ 11 位、TNBC 10 位、HER2+ 5 位。每位患者的细胞数从 631 到 8,609 不等，跨度超过 13 倍。

统计重复是患者，不是细胞。这句话在这份数据上尤其要说清楚：如果按细胞算，任意两组之间的任何差异都会显著，因为分母是几万而不是十几。三组两两比较各测 5 个可估计部分，一共 15 个检验。

三个比较的枚举规模差别很大，因此各自的最小可达 P 也不同：

三个比较全部用穷尽枚举，没有用 Monte-Carlo 近似（monte_carlo_draws = 0），所以 P 值是精确的，不带抽样误差。三个最小可达 P 都远低于 0.05——这个设计”能达到 0.05”。后面会看到，这句话几乎不提供任何有用信息。

回收份额审计：每个类别的分母都不一样

在做任何比较之前，先看一眼每个标注类别在 26 位患者里的回收情况，因为它决定了后面哪些部分有资格进入分母。

这张表里最该注意的不是均值，是跨度。T 细胞份额从 6.6% 到 83.9%，髓系从 2.0% 到 49.7%，也就是说患者之间的回收构成本身差了一个数量级。这种组内变异直接进了后面功效计算的分母——效应量以 SD 为单位，SD 大，同样的 ΔCLR 就换不来多少功效。

均值远大于中位数（恶性上皮 25.2% vs 14.7%，髓系 10.1% vs 7.0%）说明分布是右偏的：少数患者贡献了大部分某类细胞。这类形状在小样本下让均值差异特别容易被单个患者带动，这也是我们对每个检验都做留一患者审计的原因。

一个族边界的选择，让两个结果消失

BH 校正需要先声明家族。这里有两个都能自圆其说的写法：

• 按比较分族：每个两两比较内部 5 个部分各自校正。ER+ vs TNBC 的内皮 P = 0.00378，在 5 个检验里排第一，q = 0.0189；髓系 P = 0.0141，q = 0.0352。两个都过。
• 按实际做过的检验分族：15 个检验一起校正。同样两个 P，q 变成 0.0568 和 0.1057。都不过。

我们采用后者，理由很简单：这三个比较不是三项独立研究，它们出自同一份数据、同一次分析、同一批预先声明的部分，而且我们确实全部做了、全部看过。按比较分族相当于假装只做了 5 个检验，这会让实际的假阳性率高于名义值。

产物表里两个都保留：bh_q_global_15_tests 是主结论用的值，bh_q_within_contrast 只作为对照报告，并在 bh_family_note 里写明”主 q 校正覆盖每一个执行过的检验；per-contrast 值较弱，仅供比较”。把两个都摆出来，而不是选一个，是因为读者有权知道这个选择本身就能决定结论的正负。

最终 15 个检验的全局 q 全部大于 0.05，最小的是内皮的 0.0568。0/15。

“能达到 0.05”和”有检验效能”是两件事

上面那张最小可达 P 的表容易给人一个错觉：既然 0.000229 远小于 0.05，这个设计显然够用。这个推理是错的，因为它问的是”P 最小能到多少”，而不是”在实际存在的效应量下，我有多大概率检出它”。

我们对每个检验做了功效模拟，关键是在观测到的效应量处读值，而不是在网格上取最大值：

三个比较的最大功效都在 0.94 以上，看着很体面。但那个最大值对应的是一个远大于实际观测到的效应量。落到观测处，15 个检验没有一个达到 0.8，最高的是内皮的 0.753。表里因此对每一行标了 test_is_underpowered_for_this_effect = True，15/15 全部为真。

这里有两个限定必须跟着数字一起给出，否则 0.753 会被读得比它应有的分量重。

第一，模拟只在 0.5、1.0、1.5、2.0 这四个效应量格点上跑（每格 200 次试验、每次 2,000 个排列），观测处的值是在格点之间内插出来的。以 ER+ vs TNBC 为例，四个格点的功效是 0.245、0.675、0.865、0.985，观测效应量 1.206 落在 1.0 与 1.5 之间，因此 0.753 是这两点之间的插值，不是直接测出来的。它足以支持”低于 0.8”这个判断，但不适合当成一个精确的功效估计去引用。

第二，这些功效是在未校正的 0.05 阈值下算的，表里的 note 字段写明了这一点。而主结论用的是覆盖 15 个检验的 BH 校正。校正之后每个检验实际需要跨过的门槛更严，对应的功效只会更低。也就是说 0.753 是这个检验在最宽松判据下的上界，真实的、与主分析一致的功效比它还小。两个限定叠加起来，结论的方向反而更稳：15/15 欠功效这个判断不依赖于插值的精度。

这带来一个必须一起说清的结论：这一层的零发现是功效受限（zero_discoveries_is_power_limited = True），而不是分辨率受限（zero_discoveries_is_resolution_limited = False）。最小可达 P 完全够小，问题在于效应量相对组内变异不够大、样本量不够多。所以”0/15 通过 BH”不能读成”三个亚型的组成没有差别”，只能读成”这个规模的队列没有把任何一个差别确立下来”。

这个区分在实践中决定了下一步该做什么。分辨率受限的设计，加样本没用——如果 5 对 5 的最小可达 P 在 8 个检验的 BH 下天花板就是 0.0635，那么把测序做得再深、把细胞抓得再多，都不会让任何一个结果通过，唯一的出路是改设计（增加供者、改成配对、换检验族）。功效受限的设计相反：最小可达 P 已经足够小，缺的是把效应从组内变异里分出来的样本数，所以加患者是有效的，而改设计不一定。这份队列属于后者，因此”再多招 10 位患者”是一个有意义的建议；而对一个 5 对 5 的分选设计说同样的话就是浪费钱。

顺着这条线还能算出一个对规划有用的数字。ER+ vs TNBC 的内皮效应量是 1.206 个标准差，在 11 对 10 的规模下功效 0.753。要把功效推到 0.8 以上，在效应量不变的假设下需要的是每组再多几位患者，而不是每位患者多几千个细胞——因为细胞数增加只会让每位患者的份额估计更精确，不会改变患者之间的变异，而后者才是这个检验的分母。这是”n 是患者不是细胞”这句话在预算层面的直接后果。

逐部分实况

15 个检验的效应量与全局 q，按比较列出：

ER+ vs TNBC：内皮 +0.921（P = 0.0038，q = 0.0568）、髓系 −1.014（P = 0.0141，q = 0.1057）、PVL +0.579（P = 0.0747，q = 0.2089）、CAFs −0.300（P = 0.483，q = 0.6035）、T 细胞 −0.186（P = 0.701，q = 0.8087）。

ER+ vs HER2+：PVL +0.916（P = 0.0341）、内皮 +0.827（P = 0.0217）、T 细胞 −1.024（P = 0.0836）、髓系 −0.442（P = 0.272）、CAFs −0.276（P = 0.427）。全局 q 分别是 0.1279、0.1087、0.2089、0.5039、0.6035。

HER2+ vs TNBC：T 细胞 +0.838（P = 0.173）、髓系 −0.571（P = 0.302）、PVL −0.337（P = 0.452）、内皮 +0.094（P = 0.814）、CAFs −0.024（P = 0.969）。这个比较里没有任何 P 低于 0.1，且 CAFs、内皮两项在留一患者下方向会翻转。

留一患者审计的结果值得单列：15 个检验里 12 个方向稳定，3 个不稳定（全部在 HER2+ vs TNBC 这个比较里，n = 5 那一侧只要删掉一位患者就少了 20% 的信息）。方向稳定并不等于结果成立——ER+ vs TNBC 的内皮方向完全稳定，q 仍然是 0.0568。

五个部分之外：为什么恶性上皮不在检验里

进入 CLR 推断的只有 5 个部分：CAFs、内皮、髓系、PVL、T 细胞。9 个标注类别里被排除了 4 个，硬零门禁是唯一依据——一个部分必须在每一位患者里至少出现 1 个细胞：

• 恶性上皮：26 位患者中 6 位为零，被排除。
• 正常上皮：13 位为零。
• 浆细胞样：12 位为零。
• B 细胞：6 位为零。

被排除的是恶性上皮这件事本身需要解释，因为它是这份数据里第二大的类别（占全部细胞的 24.5%，平均回收份额 25.2%）。它在 6 位患者里完全没有被回收到，中位回收份额只有 14.7%，而最高的患者是 70.0%——这是解离、捕获与标注的产物，不是”这 6 位患者没有肿瘤细胞”。乳腺肿瘤组织的上皮细胞在酶解和液滴捕获中的存活率本来就低且高度可变。

如果把它强行放进几何中心，这 6 位患者的对数比值会由伪计数（0.5）而不是观测丰度决定，且因为几何中心由所有部分共享，这个虚构的数值会同时移动其余每一个部分。排除的理由到此为止，与这个部分是否”重要”无关。

伪计数带来的位移有多大，我们逐部分算了出来

“由伪计数决定”这句话可以量化。对每个部分，我们取它最稀的那位患者，算出把观测计数换成伪计数后对数丰度会偏移多少：

进入分析的 5 个部分，最坏情况下的位移是内皮的 0.154，相对于它 +0.921 的效应量约占 17%——不小，但有界，而且是可以报告的量。四个被排除的部分位移是无穷，因为观测计数为零时对数比值完全由伪计数虚构，没有任何观测信息参与。硬零门禁与”稀疏但非零”的分界就落在这里：后者引入有界偏差，可以量化并披露；前者引入的是一个纯粹的人造数值。

还有一层需要指出：最稀患者往往就是总细胞数最少的那几位（CID3941 和 CID4461 各只有 631 个细胞）。门禁的压力因此集中在浅测序的患者身上，而这些患者在患者等权的检验里与深测序患者的权重完全相同。

每位患者一票，无论他贡献了 631 还是 8,609 个细胞

患者层检验里，每位患者提供一个 CLR 值，然后按组求均值。这意味着贡献 631 个细胞的 CID3941 与贡献 8,609 个细胞的患者在检验里权重相同。

这是有意的选择，不是疏漏。按细胞数加权会让少数深测序患者主导结果，而细胞数反映的是样本处理与上机决策，不是这位患者在生物学上更有代表性。代价是浅测序患者的 CLR 估计本身噪音更大，这部分噪音进入了组内变异，也就直接压低了前面那张表里的功效——这份队列 13 倍的细胞数跨度，有一部分是通过组内变异转化成功效损失的。

如果要用加权方案，正确的做法是在看效应量之前声明并同时报告两种结果，而不是在拿到 P 值之后挑一个。本次分析只声明并执行了等权一种。

排除那 6 位患者会怎样：一个反方向的诚实结果

既然 6 位患者的恶性上皮完全缺失，一个自然的敏感性问题是：只留下 20 位有恶性上皮的患者，结论会不会变？

答案是会变，而且是往”更显著”的方向变。在 20 位患者的子队列上（ER+ 9、TNBC 8、HER2+ 3），ER+ vs TNBC 的内皮变成 +1.167（P = 0.00128，q = 0.0191），髓系变成 −1.162（P = 0.00485，q = 0.0364）。两个都低于 0.05。15 个检验的符号全部保持。

这个结果我们照实报告，但不作为发现，理由有三条：

1. 子队列不是预先声明的主分析。 排除标准来自”哪些患者缺少某个部分”，这是看过数据之后才知道的信息。把它当主结论就是在事后选择一个更有利的样本。
2. 子队列更小，却给出更小的 P。 这本身就说明这两个 q 对具体是谁进了队列高度敏感，而不是说明证据更强。HER2+ 一组只剩 3 位患者。
3. 主分析已经声明了族边界。 换队列重跑不等于换掉那个声明。

放在一起看，正确的读法是：内皮与髓系在 ER+ 与 TNBC 之间可能确实存在组成差异，但这份队列的规模无法把它确立下来，而且是否确立取决于队列定义这种本该在看数据前定好的选择。这是”待验证的线索”，不是结果。

百分点和 CLR 回答的不是同一个问题

组成数据最常见的误读是把两种量当成一回事。原始回收百分比回答的是”占回收到的细胞多少”，CLR 回答的是”相对于 5 个可估计部分的几何中心如何变化”。分母不同，两者方向可以不一致，这不矛盾。

在这份数据上尤其要小心，因为恶性上皮占了平均 25.2% 的细胞却不在参照集合里。也就是说，CLR 的参照是”5 个非上皮部分的几何中心”，而不是”整个肿瘤组成”。任何一句”某个亚型里髓系更多”都必须带上这个限定，否则读者会以为分母是全部细胞。

任何时候主语都只能是”解离、捕获、QC 和标签之后回收到的份额”，不能是组织内的绝对细胞数量。没有组织学计数或空间验证，回收份额与组织组成之间的差距无法估计。

门槛阶梯与下采样

硬零门禁的阈值 1 是我们选的，所以要检查阈值位置本身有多大影响。阈值排成 1、5、20 三档，每一档重新决定可估计集合、重算效应、P 和 q。这份数据的三档给出三个不同的可估计集合（n_distinct_estimable_sets = 3），没有任何两档共享同一组结论——与另外两篇不同，这里的阶梯真的起了区分作用，说明这份队列的稀疏结构离阈值很近。

下采样把每位患者压到 500 个细胞、重复 100 次。1,500 行结果里 1,442 行方向与全量一致，58 行翻转；1,020 行带有地板事件标记，也就是说在 500 细胞的深度上，某个部分在某位患者里掉到了门槛以下的情况相当常见。这一层说明的是：这些效应量的估计对测序/回收深度有真实依赖，不是深度无关的稳定量。

为什么组成层和程序层没有合并成 30 个检验

上一节说”族边界要覆盖每一个执行过的检验”，那么一个合理的追问是：组成层 15 个、程序层 15 个，为什么不合并成 30 个一起校正？

我们的处理是两层各自成族，理由是它们回答的是两个不同的问题，而且在分析计划里是分开声明的。组成层问”某个部分在两个亚型之间的相对丰度是否不同”，程序层问”某个转录程序在母区室内部的平均打分是否不同”。一个部分完全可以丰度不变而程序打分改变，反之亦然；把它们放进同一个 BH 序列，等于宣称我们在这 30 个检验里寻找”任意一个差异”，而实际的分析计划是两个独立的问题各自寻找。

这个选择是可以质疑的，所以要说明它的方向：分成两族比合并成一族更宽松。如果按 30 个检验校正，组成层内皮的 q 会从 0.0568 进一步变大，结论（0 个通过）不变；程序层那 4 个通过的最小 q 是 0.000926，即使乘以 2 也仍然远低于 0.05，结论同样不变。也就是说这个族边界的选择在本次数据上不改变任何结论——但这句话必须由计算说出来，不能靠”分开更自然”来假设。

不能接受的是另一个方向：在同一层里把 15 个缩成 5 个。那是把已经做过的检验从家族里删掉，而且恰好是让结论从 0 变成 2 的那一种缩法。

HER2+ 只有 5 位患者，这是整个设计的约束点

不平衡设计有一个容易被忽略的性质：决定分辨力的是最小的那一组，不是总样本量。26 位患者听起来不算少，但涉及 HER2+ 的两个比较，枚举规模分别只有 4,368 和 3,003，比 ER+ vs TNBC 的 352,716 小两个数量级。

这个差距在三个地方同时显现。一是留一患者的杠杆：n = 5 的一侧删掉一位就少掉 20% 的信息，所以 3 个方向不稳定的检验全部落在 HER2+ vs TNBC。二是功效上限：这个比较的网格最大功效是 0.940，比 ER+ vs TNBC 的 0.985 低，且在观测效应量处只有 0.100–0.268。三是程序层的覆盖：上皮程序要求母区室至少 20 个细胞，HER2+ 通过这一关的只剩 3 位患者，任何涉及 HER2+ 的程序结论都建立在 3 个数据点上。

结论是：这份队列对”ER+ 与 TNBC 有什么不同”有一定分辨力（尽管仍不足以通过校正），对”HER2+ 与另外两组有什么不同”基本没有分辨力。把三个比较并列陈述会掩盖这一点，所以逐部分实况那一节把三个比较分开列，而不是排成一张按 P 值排序的总表。

程序层：4 个通过 BH，但它们是阳性对照

程序打分层也做了 15 个检验（5 个程序 × 3 个比较），其中 4 个通过全局 BH：

• ER+ vs TNBC 管腔程序 +1.084，q = 0.000926
• ER+ vs TNBC 增殖程序 −0.201，q = 0.000926
• ER+ vs TNBC 基底程序 −0.245，q = 0.00555
• HER2+ vs TNBC 管腔程序 +0.529，q = 0.0227

这四条不是发现。ER+ 的定义就是雌激素受体阳性，而管腔程序的打分基因里第一个就是 ESR1；TNBC 的定义是三阴性，基底样表型和高增殖是它的已知特征。这四条通过 BH 说明的是分析流程能把已知的强信号找出来，也就是一个阳性对照。如果它们没通过，那才是流程有问题的信号。

属于探索性的两个程序都没通过：髓系 LAM 样程序在三个比较里 q 分别是 0.816、0.412、0.753；肌成纤维 CAF 程序是 0.698、0.400、0.400。

程序层的另一个限制是患者覆盖。三个上皮程序都要求母区室至少 20 个细胞，因此 26 位患者里只有 20 位有可用均值（HER2+ 只剩 3 位）；肌成纤维 CAF 程序剩 23 位。留出基因验证是分开做的：每个程序的打分基因与留出基因重叠为 0，细胞层 Spearman ρ 分别是增殖 0.545、管腔 0.526、CAF 0.399、LAM 样 0.233、基底 0.188。这些是细胞层相关，只能说明程序打分不是纯噪音，不能当作患者层结果。

要把内皮那条线索确立下来，需要什么

既然功效模拟已经跑过，就可以反过来回答一个更有用的问题：这条线索要多大的队列才立得住。

ER+ vs TNBC 的内皮效应量是 1.206 个标准差，在 11 对 10 处观测功效 0.753。功效随样本量增长，而这个比较的功效曲线在网格上最高能到 0.985，也就是说样本量不是无法逾越的障碍——差距不大。按同一模拟框架，要在这个效应量下把观测功效推过 0.8，需要的是每组再多两三位患者，而不是翻倍。

这件事值得单独说，因为它把”零发现”从一个泄气的结论变成一个可执行的判断：情形属于”队列差几位患者”，而非”效应太小、不值得追”。同一句话对 HER2+ vs TNBC 不成立——那里内皮的效应量只有 0.123 个标准差，观测功效 0.100，把它推过 0.8 需要的样本量远超任何单中心队列的规模，那个比较才是应该停下来的。

所以功效模拟的用途不止是给零发现一个解释。它把 15 个检验分成了两类：值得用更大队列再做一次的，和不值得的。这个分类在设计阶段就能算出来，代价只是一次模拟。

结论的可复现层

和另外两篇一样，这篇的每个数字都能从冻结产物里取回。结果树里 8 张图和 13 张表全部登记在产物清单中，逐条带 SHA256；图表清单另外记录每张图用到的源表及其哈希，因此任何一张图都能追到具体哪几行数字。

分析脚本自身被复制进 provenance 目录，与结果同树保存。这意味着读到”内皮 q = 0.0568”时，可以打开同一棵树里的脚本确认这个 q 是按 15 个检验算的，而不必相信文章的描述。运行指标（耗时、峰值内存）单独存放并显式排除在产物清单之外，因为它们不可能逐字节复现，混进清单只会让哈希校验永远失败。

这一层不产生任何新结论，但它决定了上面所有数字是可核对的还是只能被相信。对一篇主要内容是”我们没有确立任何发现”的文章来说，这一点尤其重要：阴性结论同样需要能被复查，否则它和没做过没有区别。

这份数据能写什么、不能写什么

可以写的：这是一个 11/5/10 的不平衡三组设计；三个比较全部用穷尽枚举给出精确 P；在覆盖全部 15 个检验的族边界下 0/15 通过 BH；15 个检验在各自观测效应量处全部低于 0.8 功效；恶性上皮因 6 位患者硬零而不在参照集合里；排除这 6 位后两个结果会跨过 0.05，但那是事后队列。

不能写的：任何把亚型间组成差异当作发现的表述；任何”相对于健康乳腺组织”的比较——这个队列里没有正常乳腺对照，26 位患者全是肿瘤样本，“正常上皮”是肿瘤组织内的一个标注类别，不是对照组；任何与治疗反应、预后或分子分型指导的关联，本次分析没有任何临床结局数据。

三层本次没有做，一并披露：没有做常规 QC 过滤；没有做 ambient RNA 校正，也没有给出其影响的定量估计；细胞标注沿用原始注释，没有用外部参考数据集独立重做。

复现边界

这次分析能支持的，是一个 11/5/10 不平衡三组设计的分辨力实测：三个比较全部用穷尽枚举给出精确 P，在覆盖全部 15 个检验的族边界下 0/15 通过 BH，且 15 个检验在各自观测效应量处全部低于 0.8 功效。恶性上皮因 6 位患者硬零而不在参照集合里，因此 CLR 的参照是五个非上皮部分的几何中心，不是整个肿瘤组成。

不能支持的是任何把亚型间组成差异当作发现的表述，任何”相对于健康乳腺组织”的比较（本队列没有正常乳腺对照），以及任何与治疗反应、预后或分子分型指导的关联。程序层 4 个通过 BH 的结果是阳性对照而非发现，它们的打分 marker 与亚型定义本身重叠。三层未做的工作同属边界：没有常规 QC 过滤，没有 ambient RNA 校正及其定量估计，细胞标注沿用原始注释、没有用外部参考独立重做。

可复用的判断

如果你手上是一个小样本多组设计，四条可以直接抄走：

1. 先把族边界写下来，再看 P 值。 家族定义能让同一个 P 决定两种相反的结论，事后再定就失去了错误率的含义。
2. 功效要在观测效应量处读。 曲线的最高点只说明”如果效应足够大”，与你手上的数据无关。
3. 区分零发现的两种原因。 分辨率受限（最小可达 P 都不够小）和功效受限（P 能小但效应/样本不够）需要完全不同的补救：前者要改设计，后者要加样本。
4. 把”排除某些样本后结论变强”当成警告而不是升级。 它衡量的是结论对队列定义的敏感度。

这次分析留下的不是一条亚型组成的结论，留下的是一句可以放进任何小样本项目的话：先算清这个设计在你预期的效应量下有多大概率检出它，再决定要不要做这个比较。算出来往往会省下后面全部的争论。