复现Paul et al.髓系祖细胞数据：连续谱系程序与小鼠级审计

科学审计更正（2026-07-19）：早期页面曾把erythroid、megakaryocyte、granulocyte、monocyte、dendritic-like和progenitor的最高分写成离散谱系身份，并据此计算组成。该结果已经撤回。本文现将8,751个细胞全部保持为未赋值状态。top score与top-two margin只描述程序分离，不用于命名。另一组marker只检查同一矩阵内的一致性。统计与分组解释采用8只小鼠；sort family、Cebpa/Cebpe扰动和Mouse_ID可能混杂，本文只做小鼠级描述。

Paul et al.的GSE72857记录了髓系祖细胞内部的转录异质性。它适合讨论谱系相关表达怎样逐渐展开，也最容易被误写成一套确定的成熟细胞类型。祖细胞处于连续分化过程，红系、巨核、粒细胞、单核/树突样和早期祖细胞程序可能重叠。一个最大分数无法证明命运已经不可逆确定。

本文从公开UMI表和experimental design开始，保留8,751个通过过滤的单细胞。证据链由数据结构、分选/扰动背景、连续程序、分开marker一致性和小鼠级汇总组成。离散身份、细胞比例和细胞级总体显著性不进入公开结论。

论文和公开文件

• 论文：Paul et al., Cell, 2015
• 题名：Transcriptional Heterogeneity and Lineage Commitment in Myeloid Progenitors
• DOI：10.1016/j.cell.2015.11.013
• PubMed：PMID 26627738
• 公开数据：NCBI GEO GSE72857
• 主要输入：UMI table与experimental design table

公开设计表包含Well_ID、Mouse_ID、plate、分选背景、测序批次、Number_of_cells以及CD34/FcgR3 index sorting读数。表达表与设计表必须按well严格对齐。该数据不是一个可以脱离实验设计直接聚类的普通细胞图谱。

分母和过滤

图1显示分选/扰动来源、小鼠、批次和index sorting字段。CMP、Unsorted myeloid、Cebpa perturbation、Cebpe perturbation等来源没有随机分配到所有小鼠。若某个程序在一只小鼠中偏高，它可能同时反映小鼠来源、分选策略和扰动背景。当前设计不允许把这些因素全部拆开。

UMI层质量检查

图2给出每细胞UMI counts、检测基因数以及两者关系。过滤保留Number_of_cells为单细胞的well，并去除表达覆盖不足的记录。这里的QC能支持归一化、降维和连续程序计算；它不能消除分选和扰动设计造成的结构。

老数据常见的技术风险在表头。UMI表的gene列需要显式确认，表达列和design中的Well_ID也要逐项对齐。错一列仍可能得到平滑UMAP，生物解释却会全部偏移。

分选与扰动来源进入表达空间

图3显示不同sort family和扰动来源在UMAP中的分布并不均一。这种结构是实验设计的一部分。它提醒我们，后续观察到的红系、粒细胞或树突样程序不能脱离样本来源解释，也不能把所有来源合并后写成自然骨髓中的总体谱系比例。

UMAP提供表达邻近关系，不提供发育时间。相邻细胞可能共享程序，远离也可能包含批次或分选影响。谱系方向需要marker、时间、扰动模型或功能实验共同支持。

图1：top score和margin只检查程序分离

图1汇总8,751个细胞的最高程序分数，并将其与top-two差值对照。所有分数都是当前数据集内部的标准化量，没有跨实验校准。最高分较高表示某组marker在当前细胞中突出，margin较小表示两个方向接近。两者都不能证明细胞已获得确定谱系身份。

祖细胞中程序共存有明确生物学可能。红系与巨核方向可能共享早期调控背景，单核与树突样程序也存在重叠。最新审计不再寻找一个任意阈值把细胞切成“可信”和“不可信”两组，所有记录均保持探索性和未分配。

图2：嵌入显示信号强弱，不显示身份边界

图2将每个细胞的最高连续分数映射回UMAP。高分区域说明当前marker面板在局部细胞上信号更强；低分区域提示程序表达弱、面板覆盖不足或候选程序难以区分。图中没有给细胞着上erythroid、monocyte等离散颜色。

将连续梯度画成硬标签，会制造清晰边界并诱导后续组成统计。图中只回答“哪里更容易读到某个预设程序”，不会把表达强度改写成命运承诺。

图3：第二marker面板揭示程序证据的强弱

图3给出两组marker分数的Spearman相关。dendritic-like约0.60，erythroid约0.47，monocyte约0.41，granulocyte约0.35；megakaryocyte约0.18，progenitor约0.12。相关较低提示当前程序定义或marker覆盖较弱，不能靠换一个图形把证据补足。

第二组marker没有参与生成最高程序摘要，因而能避免原样重复assignment marker。它仍来自同一批8,751个细胞和预设基因知识，只能称为一致性检查。外部参考、作者标签、克隆追踪或功能分选才有资格提供更独立的验证。

图4：小鼠级中位数暴露设计混杂

图4按Mouse_ID汇总程序中位数。Mouse 1的megakaryocyte程序明显偏高，Mouse 5的granulocyte程序较高，Mouse 8的dendritic-like程序较高，Mouse 3的monocyte程序相对上移。每个方向主要由少数小鼠贡献，不能把所有细胞合并后写成稳定总体效应。

这些小鼠同时分布在不同sort family和Cebpa/Cebpe扰动背景中。图4只描述数据结构，不把程序差异归因于基因扰动、自然发育或个体差异。要拆分这些因素，需要实验设计中存在可比较的独立重复。

连续谱系程序的正确读法

erythroid程序可由Gata1、Klf1、Alas2等相关基因支持；megakaryocyte方向关注Pf4、Itga2b等；granulocyte方向关注Mpo、Elane、Prtn3等；monocyte相关方向常见Csf1r、Ms4a7、Fcer1g。成组marker有助于判断方向，仍不能替代细胞命运观测。

程序分数适合表示“更偏向哪条表达轴”“该轴在不同小鼠中的中位水平如何”。它不适合回答某个细胞最终会成为哪类成熟细胞，也不能从UMAP顺序推出发育方向。trajectory、RNA velocity和谱系追踪各有额外假设，不能由当前score直接替代。

与源论文的对应范围

原论文研究髓系祖细胞转录异质性和谱系承诺，结合更完整的分析框架讨论连续分化。本文复核GEO公开UMI表中可直接重建的表达结构、分选背景和marker程序。二者的主线相容，证据深度不同。

本文没有复刻原文全部分支模型、转录因子网络和功能结论。Cebpa/Cebpe样本的程序变化只能结合实验设计描述，不能从这份小鼠级热图单独估计扰动效应。原文taxonomy也不能由最高分marker模块重新命名。

要把描述性差异推进到可检验结论，至少需要确认每种扰动对应多少只独立小鼠、分选门是否一致、批次是否与处理条件重合，并在小鼠层建立效应量和不确定区间。细胞数再多，也不能补足实验单位不足。若不同小鼠的方向不一致，正文应优先呈现异质性，而不是只报告合并后的平均趋势。

组成结论为什么撤回

为每个细胞选择最大程序，再统计各sort family中的标签比例，会同时受到程序共存、缺乏绝对阈值、不同小鼠细胞数和实验来源混杂影响。比例图看起来完整，却无法说明自然骨髓中的谱系丰度，也无法给出扰动的因果效应。

新流程不发布离散标签计数。top score与margin留作诊断，另一组marker给出程序一致性，小鼠热图承担生物重复层摘要。若某程序只有弱一致性或主要由一只小鼠贡献，正文会直接写出这一限制。

后续分析的最低证据要求

若要做离散注释，应引入作者逐细胞标签或经过验证的造血参考映射，并保留无法映射的细胞。若要讨论Cebpa/Cebpe扰动，需要核对每个条件下的小鼠数、分选门和批次，建立以小鼠为单位的比较。若要推断发育顺序，需要使用专门轨迹方法，并以时间、克隆或未参与建模的基因证据约束方向。

对于项目复现，优先级应放在元数据对齐和统计单位上。一个表达图可以很漂亮，Mouse_ID与实验条件混在一起时，结论仍需停留在描述层。

结论审计

小结

GSE72857能够重建髓系祖细胞中的连续表达程序，也保留了分选、扰动和小鼠来源造成的结构。8,751个细胞在最新审计中全部标记为探索性、未分配，不再生成谱系身份比例。分开marker面板显示各程序的一致性差别明显；8只小鼠的中位数又说明若干信号集中在少数动物和实验来源中。

这篇复现可以支持marker相关谱系方向的描述，不能证明单细胞命运、总体谱系丰度或扰动因果效应。数据的连续性和实验设计都需要留在正文里。

这也是后续复现必须保留的边界。