总览 · ← 返回运行 20261003-153851-search-t1-D-smoke2h

节点 n5 在终选来历上

OT 全细胞耦合的精确祖先位移：投到 30 个输入 PCA 程序模块并阻尼 λ=3 外推，再叠加 0.5×原始全基因位移，addnz 解码，均匀抽样输出，关闭生长重采样。

方法说明?节点程序自带的 METHOD.md：这个程序做了什么、为什么。

来自 programs.git 8bf3e4f97f:solution/METHOD.md

调研员的计划

代码改动?这个节点的程序和父节点程序的逐行差别：绿色是新增，红色是删除。

对比：父节点版本 3554b5ef18。改动的文件：solution/METHOD.md +79 −81、solution/README.md +4 −3、solution/run.py +305 −68

diff --git a/solution/METHOD.md b/solution/METHOD.mdindex 06c26ae..dde565a 100644--- a/solution/METHOD.md+++ b/solution/METHOD.md@@ -1,81 +1,79 @@-# ot_moscot — Waddington-OT / moscot coupling, one-step displacement extrapolation--Seed (2026-10-02) made from the G37 candidate `modeling/candidates/T1/ot_moscot/` (commit 227eeb2). Same method and-hyper-parameters. Changes for the seed contract only: the dev-only environment overrides (`G37_LAMBDA`, `G37_DECODE`,-`G37_GROWTH`, `G37_JAX_GPU`) and the unused `knn` decode branch are removed; JAX forced to CPU, torch threads fixed-at 8 (`EXECUTION.json {"gpu": false}`); `stage_pair` reads every input stage of the view the same way (no reference to-manifest `mode` / `source`). Output depends only on the view's data, the time differences between stages and `--seed`.--Contract: `python run.py --data <view> --out <pred.h5ad> --seed <int>`; `g37_common.py` must stay next to `run.py`.--## Method--Two input stages `prev` (time t0) and `last` (t1); target time t2 (final view: E8.5, E9.5 -> E10.5).--1. **Embedding.** Genes measured in both stages; top 2000 by variance (both stages pooled); z-score, clip at 10;-   PCA, 30 components (randomized, `random_state=seed`). Fitted on the input stages only.-2. **Growth prior (Waddington-OT).** Proliferation / apoptosis scores (`scanpy.tl.score_genes`, moscot's mouse-   gene lists) -> birth = generalised logistic(prolif; 1.7, 0.3, 0.25, 0.5), death = logistic(apopt; 1.7, 0.3, 0.1,-   0.2), per-day growth g = exp(birth - death) (Schiebinger 2019; same formula and defaults as moscot's-   `BirthDeathProblem.estimate_marginals`).-3. **Coupling.** `moscot.problems.time.TemporalProblem` prev -> last on the PCA (`joint_attr="X_pca"`),-   source marginal ∝ g^(t1-t0), target uniform; entropic unbalanced Sinkhorn, `epsilon=1e-3`, `tau_a=0.95`,-   `tau_b=1`, `scale_cost="mean"` (moscot tutorial settings). JAX on CPU unless `G37_JAX_GPU=1`.-4. **Output cells.** n = number of latest-stage cells clipped to `[min_cells, max_cells]` (as copy_last),-   drawn without replacement from the latest stage with probability ∝ g^(t2-t1) (the WOT birth-death model-   continued over the target interval).-5. **Displacement extrapolation.** For output cell j: ancestor mean a_j = Σ_i π_ij x_i / Σ_i π_ij (barycentric-   projection of the coupling, gene space, all panel genes); smoothed position s_j = mean expression of its 30-   nearest latest-stage cells in the PCA; step_j = λ (t2-t1)/(t1-t0) (s_j - a_j), λ = 1 (continue the last-   observed displacement at the same rate). The step is added **to the cell's non-zero entries only** and clipped-   at 0; genes not measured in both stages (external stages) get no step. Each cell keeps its own residual.--**One input stage (proxy view, E8.5 only):** steps 1, 3, 5 need two stages; only steps 2 + 4 run, i.e. a-growth-weighted copy of the latest stage (g^(t2-t1) resampling). This is the only thing the proxy can test.--Sources:-- Schiebinger G. et al. Optimal-transport analysis of single-cell gene expression identifies developmental-  trajectories in reprogramming. *Cell* 176, 928–943 (2019). doi:10.1016/j.cell.2019.01.006 (WOT: unbalanced-  entropic OT between snapshots, growth from proliferation/apoptosis signatures, birth-death logistic).-- Klein D., Palla G., Lange M. et al. Mapping cells through time and space with moscot. *Nature* 638, 1065–1075-  (2025). doi:10.1038/s41586-024-08453-2 (TemporalProblem; code moscot 0.5.2, BSD-3).-- Cuturi M. Sinkhorn distances. NeurIPS 2013; Chizat L. et al. Scaling algorithms for unbalanced optimal-  transport problems. *Math. Comp.* 87, 2563–2609 (2018).-- Extrapolating the barycentric displacement one more step is our use of the coupling (WOT/moscot interpolate,-  they do not extrapolate); listed in agent/knowledge/T1_methods_landscape.md §1.--## Data / knowledge used--Only the view's input stages. Generic knowledge: moscot's built-in mouse proliferation (97) and apoptosis (193)-gene lists (`moscot.utils.data`, from the WOT paper; stage-agnostic gene-function annotation). No held-out stage,-no information from (E9.5, E13.5], no pre-trained weights.--## Hyper-parameters--| Name | Value | Where it came from |-|---|---|---|-| `N_HVG`, `N_PCS`, `K_SMOOTH` | 2000, 30, 30 | a priori; not tuned |-| `EPSILON`, `TAU_A`, `TAU_B`, `scale_cost` | 1e-3, 0.95, 1, mean | a priori (moscot tutorial settings) |-| growth prior | on | a priori (WOT / moscot defaults); G37 saw it cost ~3 points on the old proxy, kept on |-| `LAMBDA` | 1 | a priori (continue the observed displacement at the same rate); G37 also ran 0.5 on X3 (50.4 vs 50.2), not changed |-| decode `addnz` | — | **chosen on the X3 ruler (G37)** against `add` (25.6) and `knn` (48.7): the dense step destroys the zero pattern. Also a first-principles choice (the scorer compares sparse log-expression), but the evidence that picked it was X3 |--No re-tuning for the seed.--## Resources (Spark, CPU)--Final view (16.8k x 17.1k coupling): 69 s, max RSS 6.7 GB; proxy 4 s / 1.6 GB; X3 ~ 8 s / 2.3 GB; proxy2 ~ 20 s / 4.3 GB.-The dense coupling (n_prev x n_last float32, ~1.1 GB on final) and the dense earlier stage (~2.2 GB) dominate memory.--## Findings (G37, local scorer: fast engine, truth half B, scorer seed = program seed)--- Dense gene-space step (`add`) is destructive: X3 25.6 at λ=1, still 32.9 at λ=0.25 (cell_state, covariation-  collapse) — the step makes every zero slightly positive. `addnz` fixes it (X3 50.2 at λ=1, 50.4 at λ=0.5);-  `knn` decode 48.7 (λ=1). Growth resampling has no effect on X3/proxy2 (all cells are kept there).-- proxy2 (E8.5 official -> Qiu E9.0 heart, then +0.5 d): every variant ≈ copy of the Qiu cells (~27.5), the-  cross-dataset step is batch effect.-- Full eval (seeds 0-2, half B): proxy 46.8 / 46.8 / 47.9 (copy_last 50.0 / 50.3 / 50.1) — the WOT growth-  resampling alone costs ~3 points on E8.5 -> E9.5 (direction, cell_state); X3 50.2 / 50.1 / 49.8 (copy_last 50.0,-  pseudobulk_shift 40.5-40.7); proxy2 27.7-27.9 (copy_last 27.4-27.6). `G37_GROWTH=0` turns the proxy into copy_last.-- Final-view prediction (`~/vec/scratch/g37/ot_moscot/final.h5ad`, seed 0): board-cap cells, only the last input's types,-  composition close to the last input's; no new states (the method cannot create them).+OT 全细胞耦合的精确祖先位移：投到 30 个输入 PCA 程序模块并阻尼 λ=3 外推，再叠加 0.5×原始全基因位移，addnz 解码，均匀抽样输出，关闭生长重采样。++# 节点 5（improve，父节点 2 ot_moscot）：local_ot — 模块投影阻尼位移 + 弱原始残差++## 最终方法（run.py 默认配置）++两输入阶段（X3: E8.75+E9.0 → E9.5；final: E8.5+E9.5 → E10.5）：++1. **输出细胞**：从最新阶段均匀无放回抽取 `target_n_cells` 个。父节点的 g^dt 生长重采样**默认关闭**（在本管线里它降低 X3，见下表；`--ablate growth_on` 保留为对照）。+2. **联合嵌入**：沿用父节点 `g37_common.embed`（两阶段合并，只在两阶段都实测的基因里选 HVG 2000，z-score clip10，30 PCs）。30 个 PCA 载荷 = 数据驱动的“程序模块” W（只用输入阶段拟合）。+3. **OT 耦合**：moscot TemporalProblem，全部细胞（非父节点之外的子集），ε=1e-3，τ_a=0.95，τ_b=1，scale_cost="mean"，源边际 = Waddington-OT 增殖/凋亡生长先验（同父节点）。每个输出细胞 j 用它自己的耦合列（精确祖先重心 q_j），不做近邻平滑映射。+4. **位移**：d_j = kNN30 平滑的最新阶段位置 s_j − 祖先重心 a_j（原始 log 空间）。+   - 模块投影步：c_j = W·d_j/sd /(1+α)（W 正交，岭回归解析解，α=1），step_proj = λ·(dt_out/dt_in)·sd·Wᵀc_j，λ=3；+   - 混合步：step = 0.5·(dt_out/dt_in)·d_j（**全部**两阶段实测基因）+ step_proj（HVG 上叠加）；+   - 只加到该细胞的非零元上（addnz），clip ≥0，零模式保留细胞自身残差。+5. **单输入退路**（proxy10）：耦合需要两个阶段 → 均匀照抄最新阶段（与 copy_last 同流程）。分支只看输入阶段个数，不看视图身份/绝对时间。++实现但**默认关闭**的机制（PLAN 机制 B）：按局部检测率趋势的受约束基因开启/关闭（候选 = 细胞 top-3 模块的 top-50 载荷基因、|Δp|≥0.05、每细胞上限 M、开启值从近邻非零值抽取）。默认 `M_CAP=0`，`--mcap 30` 可启用。关闭原因见“被证否的部分”。++实现但默认不用的稳定性机制（PLAN 机制 A 的收缩部分）：`--subset` 跑两个 80% 子集耦合并按 SNR_k=|两次平均|/|两次差| 收缩模块系数。X3 上精确全细胞耦合一致更好（下表 #16 vs #8），故默认 `--exact`（w_k 全部=1，收缩退化为恒等）。++## 查分记录（全部 20 次额度；X3 = A 半，seed 0，均按 vec-score score_parts）++| # | 配置（X3） | de_score | de_dir | mmd_u | variogram | 总分 |+|---|---|---|---|---|---|---|+| 1 | 父节点原码复跑 | （截断） | 12.26 | 14.96 | 10.69 | ≈50.6 |+| 2 | 父节点原码 @proxy10 | 12.40* | 13.61* | 14.25* | 9.41* | ≈49.7* |+| 3 | 子集耦合 hvg λ0.5 M30（PLAN 默认） | 11.23 | 12.57 | 15.08 | 9.95 | 48.83 |+| 4 | 同上 M0（关开关） | 11.04 | 12.58 | 15.25 | 10.00 | 48.88 |+| 5 | `--ablate displacement`（均匀照抄） | 10.86 | 12.53 | 15.16 | 10.05 | 48.60 |+| 6 | 子集 hvg λ1 M0 | 11.32 | 12.63 | 15.36 | 9.99 | 49.30 |+| 7 | 子集 hvg λ1 M30 | 11.62 | 12.63 | 14.94 | 9.88 | 49.07 |+| 8 | 子集 hvg λ2 M0 | 11.62 | 12.65 | 15.62 | 9.95 | 49.84 |+| 9 | 子集 full λ1（原始全基因步） | 11.62 | 12.39 | 15.82 | 9.30 | 49.13 |+| 10 | 子集 full λ2 | 11.83 | 12.24 | 15.41 | 8.81 | 48.30 |+| 11 | 子集 hvg λ2 M0 + growth | 11.83 | 12.66 | 15.53 | 9.19 | 49.20 |+| 12 | 子集 hvg λ0.5 M30 + growth | 11.42 | 12.58 | 14.95 | 9.17 | 48.12 |+| 13 | 子集 full λ1 + growth | 11.72 | 12.46 | 15.44 | 8.54 | 48.17 |+| 14 | 精确 full λ1 + growth（≈父节点机制） | 11.94 | 12.19 | 13.52 | 10.10 | 47.74 |+| 15 | 精确 full λ1 | 12.15 | 12.25 | 14.29 | 10.72 | 49.42 |+| 16 | 精确 hvg λ2 M0 | 11.62 | 12.71 | 15.65 | 10.14 | 50.12 |+| 17 | 精确 hvg λ3 M0 | 11.83 | 12.70 | 17.57 | 10.21 | 52.31 |+| 18 | 精确 hybrid λ2 + 0.5·raw | 11.83 | 12.67 | 17.73 | 10.59 | 52.82 |+| 19-20 | **精确 hybrid λ3 + 0.5·raw（最终默认）** | 11.73 | 12.68 | **20.03**(raw 0.0218) | 10.50 | **54.94** |++（* #2 的 score_parts 取自父节点变化量表 proxy10 行。#19 与 #20 是同一预测文件：一次输出解析失败、一次重查确认，额度因此耗尽。）++λ 方向在 hvg 投影步上单调向好（48.88 → 49.30 → 49.84/50.12 → 52.31），加 0.5 倍原始残差再 +0.5~2.6。**未验证 λ4**（额度耗尽），最终配置取已验证的最优点而不是外推点。++## 机制证据（PLAN mechanism_evidence 对应项）++1. **型内位移方差占比 0.965**（精确模式，X3 seed 0，stderr 诊断；标签仅用于诊断）：位移在类型内部高度异质，不是每型常向量。+2. 开关机制（默认关）：λ0.5 时每细胞开启 7.2 / 关闭 6.7 个基因，λ1 时 12.8 / 12.1 个，上限 M=30 未触发；模块 w_k（子集模式）0.57–0.85。+3. 零值比例：最新阶段 0.9394 → 输出 0.9422（+0.28pp < 2pp 护栏）；方差/离群未塌缩（100% 细胞被位移触及，step max 9.4 但中位小）。+4. 被改动细胞比例：100%（addnz 步作用于每个细胞的非零元）。+5. **消融对照**（`--ablate`，harness 用第一个标识符 `displacement`）：+   - `displacement` → 均匀照抄，输出显著不同（#5: 48.60 vs 默认 54.94）→ mechanism_active 应判 yes；+   - `opening` → M=0，与默认输出**逐位相同**（默认即 M=0，如实说明：开关机制不在最终配置中）；+   - `growth_on` → 恢复 g^dt 重抽样，输出不同（#11/#14：降分）。++## 被证否的部分（如实报告）++- **PLAN 机制 B（检测率开关基因）在已测幅度下有害或无效**：同幅度 M30 vs M0 两次对照（#3 vs #4：48.83/48.88；#6 vs #7：49.30/49.07），mmd_u 与 variogram 均被拉低，de_score 无稳定收益。最终配置 M=0，代码保留可启用。+- **生长重采样在本管线的 X3 上有害**（#11 vs #8、#12 vs #3、#14 vs #15：variogram/mmd 下降），默认关闭，只留对照。+- **子集耦合 + 近邻映射的祖先平滑不如精确耦合列**（#9 vs #15、#8 vs #16），SNR 收缩在精确模式下不可辨识（单次耦合，w≡1）。PLAN 的收缩机制被更强的“精确祖先 + 模块投影去噪”替代。++## 验证过 / 未验证++- 已验证：X3 与 proxy 两视图 `vec-check` ok；seed 0 重跑逐字节相同（确定性）；seed 1、2 正常运行且诊断量接近；**伪装视图**（时间统一 +1 天、manifest 键序打乱、路径更换）输出与真实视图逐位相同（view-independent）；X3 全程 13–15 s、CPU-only（EXECUTION.json gpu=false）。+- 未验证（额度耗尽）：proxy10 上默认配置的查分（退路 = 均匀照抄，预期 ≈ 节点 1 copy_last 的 52.75，抽样差异 ≤0.5）；最终配置在 seed 1/2 的分数；λ=4；开关机制在 λ3 混合步下的表现；final 视图实跑（父节点同构型全耦合在 final 约 70 s / 6.7 GB，本方法内存轮廓与父节点相当）。+- B 半风险：X3 A 半 54.94，正式分用 B 半，噪声约 ±2；但相对父节点 A 半复跑（≈50.6）的 +4.3 超出噪声，方向由 8 个配置的单调趋势支持。++## 知识来源++- 未使用任何保留阶段/保留基因型的测量值；未使用 `uns.celltype_palette`；无硬编码类型名单、阶段名、细胞数、比例。程序只用视图内表达、标签集合外的时间**差**（dt_in、dt_out）。+- 发育事件表（G42.2 附件）：**未使用**（本方法是纯数据驱动的 OT+线性模块外推，不含标记程序或事件顺序）。+- Waddington-OT 生长先验（增殖/凋亡基因评分 → 出生/死亡率）：Schiebinger et al. 2019, Cell, DOI 10.1016/j.cell.2019.01.006（经 moscot 实现；marker 列表是通用细胞状态标记，非阶段特异测量）。父节点 2 继承。+- 熵正则非平衡 OT（Sinkhorn，τ_a/τ_b 边际松弛）：Peyré & Cuturi 2019（moscot/otter 实现）。父节点继承。+- 评分规则简报 §3 事实 5（组成/位移需与真实变化同向才得分）用于解读：λ 增大时 mmd_u/de_dir 同向改善，说明外推方向与真值一致。diff --git a/solution/README.md b/solution/README.mdindex dc1cbaf..d311d94 100644--- a/solution/README.md+++ b/solution/README.md@@ -1,4 +1,5 @@-# ot_moscot+# local_ot 改进（节点 5，父节点 2 ot_moscot） -Waddington-OT / moscot：两个最新输入阶段在联合 PCA 上做非平衡熵 OT 耦合（WOT 增殖/凋亡先验生长），最新阶段按 g^Δt 重抽样，每个细胞沿“kNN 平滑位置 − 耦合祖先均值”再走一个等比例步长，只加在非零基因上，夹到 ≥0。-来自 G37 候选 modeling/candidates/T1/ot_moscot；addnz 解码是在 X3 上选的（详见 METHOD.md）。纯 CPU，final 约 70 s、6.7 GB。+moscot 全细胞 OT 耦合（WOT 生长先验边际）给出每个输出细胞的精确祖先重心；位移 = kNN30 平滑位置 − 祖先重心，投到 30 个输入 PCA 模块上、阻尼 λ=3·(dt_out/dt_in) 外推，另加 0.5·(dt_out/dt_in) 的原始全基因位移；只作用于非零元（addnz），clip ≥0；输出细胞从最新阶段均匀抽取（生长重采样默认关，`--ablate growth_on` 为对照）。单输入视图退化为均匀照抄。检测率开关基因机制已实现但默认关（M=0，实测有害，见 METHOD.md）。+纯 CPU（EXECUTION.json gpu=false），X3 约 15 s；内存轮廓与父节点相当（final 预估 <8 GB）。+详细方法、20 次查分记录、消融与验证见 METHOD.md。diff --git a/solution/run.py b/solution/run.pyindex a8a575c..536dfce 100644--- a/solution/run.py+++ b/solution/run.py@@ -1,27 +1,26 @@ #!/usr/bin/env python3-"""ot_moscot: Waddington-OT / moscot TemporalProblem coupling, extrapolated one step past the latest input stage.--Two input stages (final: E8.5, E9.5):-  1. joint PCA of both stages (HVG, z-score, 30 PCs; fitted on the inputs only);-  2. moscot TemporalProblem prev -> last on the PCA, source marginals from Waddington-OT prior growth-     (proliferation / apoptosis gene scores), entropic unbalanced Sinkhorn (epsilon 1e-3, tau_a 0.95);-  3. output cells = latest-stage cells resampled with weights g^dt_out (prior growth rate continued for the target-     interval, Waddington-OT birth-death model);-  4. each output cell j moves by LAMBDA * dt_out / dt_in * (kNN-mean(x_j) - ancestor_mean_j) in gene space, where-     ancestor_mean_j is the coupling-weighted (barycentric) mean of its ancestors in the earlier stage: the last-     observed displacement continued for the target interval; the cell keeps its own residual. The step is applied-     to the cell's non-zero entries only (DECODE "addnz": a dense step turns every zero into a small positive value-     and wrecks cell_state / covariation, see METHOD.md). Clipped at 0.---ablate (G39.7): `growth` turns step 3 off (uniform resampling); any other name turns the OT displacement (steps 1, 2,-4) off, leaving the growth-weighted copy of the latest stage.-One input stage (proxy: E8.5 only): steps 1, 2, 4 need two stages; only the growth resampling (3) runs, i.e.-a growth-weighted copy of the latest stage.--Seed version (agent/seeds/T1__val/ot_moscot, 2026-10-02) of modeling/candidates/T1/ot_moscot (G37): same method and-hyper-parameters; the dev-only environment overrides (G37_LAMBDA / G37_DECODE / G37_GROWTH / G37_JAX_GPU) and the-`knn` decode branch are removed; JAX and torch on CPU, fixed thread count (EXECUTION.json gpu false).-Parameter provenance (METHOD.md): DECODE = addnz was chosen on the X3 ruler (G37); LAMBDA = 1 a priori, also checked-on X3 (0.5 vs 1 within 0.2).+"""local_ot improved: damped module-projected OT displacement + weak raw residual, uniform output sampling.++Chosen configuration (all hyper-parameters selected on the X3 ruler, see METHOD.md):+  0. output cells = uniform sample without replacement of the latest stage. Growth resampling is OFF by default+     (`--ablate growth_on` restores the parent's g^dt_out resampling; it lowered X3 in this pipeline);+  1. joint embed of the two input stages (HVG 2000 among genes measured in both, z-score clip10, 30 PCs);+     modules W = the 30 PCA loadings; one moscot TemporalProblem coupling prev -> last on ALL cells of both+     stages (entropic unbalanced Sinkhorn, epsilon 1e-3, tau_a 0.95, tau_b 1, source marginals from the+     Waddington-OT growth prior), exact per-cell ancestor barycentre q_j (column of the coupling of cell j);+  2. displacement d_j = kNN30-smoothed latest-stage position s_j - a_j (raw log space);+     projected part: c_j = W d_j,z / (1+ridge alpha) (z-units), step_proj = LAMBDA*(dt_out/dt_in)*sd*W^T c_j;+     hybrid step = RES_LAMBDA*(dt_out/dt_in)*d_j on ALL masked genes  +  step_proj on the HVG;+  3. step added to the cell's non-zero entries only ("addnz"), clipped at 0; zeros keep the cell's own pattern;+  4. implemented but OFF by default (M_CAP=0): detection-rate-constrained gene opening/closing on module+     top-loading genes (`--mcap 30` enables; harmful at the tested amplitudes, METHOD.md);+  5. --subset (dev): two 80%-subset couplings + per-module SNR shrinkage w_k instead of the exact coupling.+--ablate displacement (any unrecognised name too): no step -> uniform copy of the latest stage.+--ablate opening: M=0 (identical to the default configuration).+--ablate growth_on: g^dt_out resampling of output cells, rest unchanged.+One input stage (proxy): the coupling needs two stages -> uniform copy of the latest stage (copy_last flow).+All quantities are computed from the view at run time; only time DIFFERENCES are used, no absolute stage times,+no hard-coded stage names, cell counts or type lists (view-independent per CONTRACT). """  from __future__ import annotations@@ -29,6 +28,7 @@ from __future__ import annotations import argparse import os import sys+import time from pathlib import Path  os.environ.setdefault("XLA_PYTHON_CLIENT_PREALLOCATE", "false")@@ -36,93 +36,330 @@ os.environ.setdefault("XLA_PYTHON_CLIENT_MEM_FRACTION", "0.1") os.environ["JAX_PLATFORMS"] = "cpu"  # CPU only: deterministic, no GPU slot  import numpy as np+from scipy import sparse  sys.path.insert(0, str(Path(__file__).resolve().parent))-from g37_common import embed, growth_rates, knn_mean, stage_pair  # noqa: E402+from g37_common import embed, growth_rates, stage_pair  # noqa: E402  from src.task1_temporal.view_io import load_manifest, panel_genes, target_n_cells, write_prediction  # noqa: E402  N_HVG = 2000-N_PCS = 30+N_PCS = 30         # K modules = joint-input PCA loadings EPSILON = 1e-3 TAU_A = 0.95 TAU_B = 1.0-K_SMOOTH = 30-LAMBDA = 1.0     # 1 = continue the observed displacement at full rate-GROWTH = True    # resample output cells by g^dt_out (WOT birth-death model)+K_SMOOTH = 30      # neighbours for smoothing + detection rates + opening values+K_MAP = 10         # nearest subset cells used to map an output cell onto the subset coupling+SUB_CAP = 4000     # coupling subset cap per stage+SUB_FRAC = 0.8     # always keep two DIFFERENT subsets so the SNR shrinkage is estimable+RIDGE_ALPHA = 1.0+LAMBDA = 3.0       # damping of the module-projected extrapolated step (chosen on X3: 0.5/1/2/3 monotone better)+RES_LAMBDA = 0.5   # weight of the raw (all-masked-gene) displacement added to the projected step+M_CAP = 0          # per-cell cap on openings/closings; 0 = mechanism B off (harmful at tested amplitudes, METHOD.md)+TOP_MODS = 3       # modules whose top-loading genes are opening/closing candidates+TOP_LOAD = 50      # genes per module+DP_MIN = 0.05      # |p_next - p_last| threshold for candidates N_THREADS = 8  +def _t0():+    return time.time()+++def _log(msg, t0):+    print(f"[{time.time() - t0:7.1f}s] {msg}", file=sys.stderr, flush=True)++ def weighted_rows(w: np.ndarray, n: int, rng: np.random.Generator) -> np.ndarray:     p = w / w.sum()     return np.sort(rng.choice(len(w), size=n, replace=n > len(w), p=p))  +def subset_size(n: int) -> int:+    m = min(n, SUB_CAP)+    return int(min(m, max(round(SUB_FRAC * n), min(n, 1000))))+++def run_coupling(Zp, Zl, gr_prev, gr_last, sub_p, sub_l, tag, eps):+    import anndata as ad+    import pandas as pd+    from moscot.problems.time import TemporalProblem++    obs = pd.DataFrame({"time": np.r_[np.zeros(len(sub_p)), np.ones(len(sub_l))],+                        "growth": np.r_[gr_prev[sub_p], gr_last[sub_l]]},+                       index=[f"p{i}" for i in range(len(sub_p))] + [f"l{i}" for i in range(len(sub_l))])+    obs["time"] = obs["time"].astype(float)+    small = ad.AnnData(X=np.zeros((len(obs), 1), dtype=np.float32), obs=obs)+    small.obsm["X_pca"] = np.vstack([Zp[sub_p], Zl[sub_l]])+    tp = TemporalProblem(small)+    tp = tp.prepare(time_key="time", joint_attr="X_pca", a="growth")+    tp = tp.solve(epsilon=eps, tau_a=TAU_A, tau_b=TAU_B, scale_cost="mean")+    sol = tp[(0.0, 1.0)].solution+    print(f"ot[{tag}]: converged={getattr(sol, 'converged', None)} cost={getattr(sol, 'cost', None)}",+          file=sys.stderr, flush=True)+    return np.asarray(sol.transport_matrix, dtype=np.float32)+++def ancestor_weights(Zl, Zl_sub, P, rows, sub_l=None):+    """(n_out, n_prev_sub): coupling-weighted ancestor distribution of each output cell.++    The coupling was solved on subsets; each output cell (a full latest-stage cell) is mapped through its+    K_MAP nearest cells of the latest-stage subset, averaging their column-normalised coupling rows.+    With --exact the subset is every latest-stage cell and the cell's own (exact) coupling column is used.+    """+    from sklearn.neighbors import NearestNeighbors++    if sub_l is not None and len(sub_l) == Zl.shape[0]:+        q = np.asarray(P[:, rows].T, dtype=np.float32)  # (n_out, n_prev): exact coupling column per output cell+        q /= np.maximum(q.sum(axis=1, keepdims=True), 1e-30)+        return q+    Q = P / np.maximum(P.sum(axis=0, keepdims=True), 1e-30)+    nn = NearestNeighbors(n_neighbors=min(K_MAP, Zl_sub.shape[0])).fit(Zl_sub)+    idx = nn.kneighbors(Zl[rows], return_distance=False)  # (n, k) into subset+    q = np.asarray(Q[:, idx].mean(axis=2), dtype=np.float32).T  # (n, n_prev_sub)+    return q+++def cap_events(ev: np.ndarray, prob: np.ndarray, M: int) -> np.ndarray:+    """Keep at most M events per row, those with the highest probability."""+    if M <= 0:+        ev[:, :] = False+        return ev+    cnt = ev.sum(axis=1)+    for j in np.where(cnt > M)[0]:+        idxs = np.where(ev[j])[0]+        order = idxs[np.argsort(-prob[j, idxs], kind="stable")]+        ev[j, order[M:]] = False+    return ev++ def main() -> None:+    t0 = _t0()     parser = argparse.ArgumentParser()     parser.add_argument("--data", required=True)     parser.add_argument("--out", required=True)     parser.add_argument("--seed", type=int, default=0)-    # G39.7 mechanism-off control: `growth` -> no growth resampling (uniform); any other name (the OT displacement,-    # e.g. transport / displacement / mechanism) -> no displacement step: the growth-weighted copy of the latest stage+    # G39.7 mechanism-off controls: `displacement` (any unrecognised name too) -> uniform copy of the latest+    # stage; `opening` -> M=0; `growth_on` -> restore g^dt resampling.     parser.add_argument("--ablate", default=None)+    parser.add_argument("--lam", type=float, default=None, help="dev only: override LAMBDA")+    parser.add_argument("--mcap", type=int, default=None, help="dev only: override M_CAP")+    parser.add_argument("--epsilon", type=float, default=None, help="dev only: override EPSILON")+    parser.add_argument("--proj", default="hybrid", choices=["hvg", "full", "hybrid"],+                        help="dev only: hvg = module-projected step; full = raw displacement on all masked genes;"+                             " hybrid = lam*proj + res_lam*raw (default)")+    parser.add_argument("--reslam", type=float, default=None, help="dev only: override RES_LAMBDA")+    parser.add_argument("--exact", action="store_true", default=True,+                        help="single full-cell coupling with exact per-cell ancestor barycentres (default)")+    parser.add_argument("--subset", action="store_true",+                        help="dev only: two 80%%-subset couplings with SNR shrinkage instead of --exact")     args = parser.parse_args() +    lam = LAMBDA if args.lam is None else float(args.lam)+    m_cap = M_CAP if args.mcap is None else int(args.mcap)+    eps = EPSILON if args.epsilon is None else float(args.epsilon)+    res_lam = RES_LAMBDA if args.reslam is None else float(args.reslam)+    args.exact = args.exact and not args.subset++    ablate = args.ablate+    displacement_off = ablate is not None and ablate not in ("opening", "growth_on")+    growth_on = ablate == "growth_on"+     manifest = load_manifest(args.data)     genes = panel_genes(args.data, manifest)     prev, last, mask, dt_in, last_e = stage_pair(args.data, manifest, genes)     dt_out = float(manifest["target"]["time"] - last_e["time"])     rng = np.random.default_rng(args.seed)     n = target_n_cells(manifest, last.n_obs)-    growth_on = GROWTH and args.ablate != "growth" -    if prev is None or (args.ablate and args.ablate != "growth"):+    # ---- output cell selection (growth resampling only under --ablate growth_on) ----+    if growth_on:         (g_last,) = growth_rates([last], genes, mask, args.seed)-        rows = weighted_rows(g_last ** dt_out if growth_on else np.ones(last.n_obs), n, rng)+        sel_w = g_last ** dt_out+        _log(f"diag: growth last min/median/max {g_last.min():.3f}/{np.median(g_last):.3f}/{g_last.max():.3f}", t0)+    else:+        sel_w = np.ones(last.n_obs)+    rows = weighted_rows(sel_w, n, rng)++    if prev is None or displacement_off:+        _log(f"fallback path (prev={prev is None}, displacement_off={displacement_off}); uniform copy", t0)         write_prediction(last.X[rows], genes, args.out, seed=args.seed)         return -    import anndata as ad-    import pandas as pd     import torch-    from moscot.problems.time import TemporalProblem+    torch.set_num_threads(N_THREADS)++    ratio = dt_out / dt_in+    _log(f"dt_in={dt_in} dt_out={dt_out} ratio={ratio} n_out={n} lam={lam} M={m_cap} eps={eps}", t0)++    # ---- joint embedding (modules W = PCA loadings on HVG, z-scored) ----+    (Zp, Zl), info = embed([prev.X, last.X], mask, N_HVG, N_PCS, args.seed)+    hvg, sd, W = info["hvg"], info["sd"], info["components"]  # W: (K, H) orthonormal+    H = len(hvg)+    _log(f"embed done: {Zp.shape} {Zl.shape}, hvg {H}", t0) -    (Zp, Zl), _ = embed([prev.X, last.X], mask, N_HVG, N_PCS, args.seed)     g_prev, g_last = growth_rates([prev, last], genes, mask, args.seed)+    gr_prev, gr_last = g_prev ** dt_in, g_last ** dt_in -    obs = pd.DataFrame({"time": np.r_[np.zeros(prev.n_obs), np.ones(last.n_obs)],-                        "growth": np.r_[g_prev ** dt_in, g_last ** dt_in]},-                       index=[f"p{i}" for i in range(prev.n_obs)] + [f"l{i}" for i in range(last.n_obs)])-    obs["time"] = obs["time"].astype(float)-    small = ad.AnnData(X=np.zeros((len(obs), 1), dtype=np.float32), obs=obs)-    small.obsm["X_pca"] = np.vstack([Zp, Zl])-    tp = TemporalProblem(small)-    # source marginals = WOT prior growth over the input interval (normalised by moscot); target uniform-    tp = tp.prepare(time_key="time", joint_attr="X_pca", a="growth")-    tp = tp.solve(epsilon=EPSILON, tau_a=TAU_A, tau_b=TAU_B, scale_cost="mean")-    sol = tp[(0.0, 1.0)].solution-    print(f"ot: converged={getattr(sol, 'converged', None)} cost={getattr(sol, 'cost', None)}", file=sys.stderr)-    P = np.asarray(sol.transport_matrix, dtype=np.float32)  # (n_prev, n_last)--    rows = weighted_rows(g_last ** dt_out if growth_on else np.ones(last.n_obs), n, rng)-    print(f"diag: growth last min/median/max {g_last.min():.3f}/{np.median(g_last):.3f}/{g_last.max():.3f}", file=sys.stderr)-    Pc = P[:, rows]-    del P-    Pc /= np.maximum(Pc.sum(axis=0, keepdims=True), 1e-30)-    factor = LAMBDA * dt_out / dt_in-    torch.set_num_threads(N_THREADS)-    Xp = torch.from_numpy(prev.X.toarray())-    anc = (torch.from_numpy(Pc).T @ Xp).numpy()  # barycentric ancestor mean, (n, n_genes)-    del Xp, Pc-    smooth = knn_mean(Zl, last.X, rows, K_SMOOTH)-    step = (smooth - anc) * factor-    del smooth, anc-    step[:, ~mask] = 0.0-    X = last.X[rows].toarray()-    step *= X > 0  # "addnz": move only measured (non-zero) entries, keeps each cell's zero pattern-    X += step+    # ---- latest-stage neighbourhoods of the output cells ----+    from sklearn.neighbors import NearestNeighbors+    nn_last = NearestNeighbors(n_neighbors=min(K_SMOOTH, last.n_obs)).fit(Zl)+    nb_idx = nn_last.kneighbors(Zl[rows], return_distance=False)  # (n, K) into last+    Wk = sparse.csr_matrix((np.full(nb_idx.size, 1.0 / nb_idx.shape[1], dtype=np.float32),+                            nb_idx.ravel(), np.arange(0, nb_idx.size + 1, nb_idx.shape[1])),+                           shape=(n, last.n_obs))+    Xl_hvg = np.asarray((last.X[:, hvg]).todense(), dtype=np.float32)  # (n_last, H)+    s_hvg = np.asarray(Wk @ Xl_hvg, dtype=np.float32)                  # smoothed position (n, H)++    # ---- two subset couplings -> two ancestor means / coefficient sets ----+    nsub_p, nsub_l = subset_size(prev.n_obs), subset_size(last.n_obs)+    _log(f"coupling subsets: prev {nsub_p}/{prev.n_obs}, last {nsub_l}/{last.n_obs}", t0)+    if args.exact:+        nsub_p, nsub_l = prev.n_obs, last.n_obs+    qs, Cs, prev_subs, A_full = [], [], [], []+    for r in ((0,) if args.exact else (0, 1)):+        rgr = np.random.default_rng(args.seed + r)+        sub_p = np.arange(prev.n_obs) if args.exact else np.sort(rgr.choice(prev.n_obs, size=nsub_p, replace=False))+        sub_l = np.arange(last.n_obs) if args.exact else np.sort(rgr.choice(last.n_obs, size=nsub_l, replace=False))+        P = run_coupling(Zp, Zl, gr_prev, gr_last, sub_p, sub_l, f"{args.seed}+{r}", eps)+        q = ancestor_weights(Zl, Zl[sub_l], P, rows, sub_l=sub_l)+        del P+        Xps = np.asarray((prev.X[sub_p][:, hvg]).todense(), dtype=np.float32)+        a_hvg = q @ Xps  # (n, H) raw log space+        d_z = (s_hvg - a_hvg) / sd+        Cs.append((d_z @ W.T) / (1.0 + RIDGE_ALPHA))  # (n, K) ridge on orthonormal W+        qs.append(q)+        prev_subs.append(sub_p)+        if args.proj in ("full", "hybrid"):+            A_full.append(q @ np.asarray(prev.X[sub_p].todense(), dtype=np.float32))+        _log(f"coupling {r} done, |c| mean {np.abs(Cs[-1]).mean():.4f}", t0)++    C0, C1 = Cs[0], Cs[-1]+    Cbar = 0.5 * (C0 + C1)+    num = 0.5 * (np.abs(C0) + np.abs(C1)).mean(axis=0)+    den = np.abs(C0 - C1).mean(axis=0)+    snr = np.where(den > 1e-12, num / np.maximum(den, 1e-12), np.inf)+    w_shrink = (snr ** 2) / (1.0 + snr ** 2)+    w_shrink = np.where(np.isfinite(w_shrink), w_shrink, 1.0)+    _log(f"module shrinkage w: min {w_shrink.min():.3f} med {np.median(w_shrink):.3f} max {w_shrink.max():.3f} "+         f"n>0.5 {(w_shrink > 0.5).sum()}/{len(w_shrink)}", t0)++    # ---- step: module-projected, shrunk, damped displacement (raw log space, HVG only) ----+    factor = lam * ratio+    step = factor * ((Cbar * w_shrink) @ W) * sd  # (n, H)++    # ---- opening / closing pool: top-50 |loading| genes of every kept module ----+    opening_off = (args.ablate == "opening") or m_cap == 0+    Aw = np.abs(Cbar) * w_shrink  # (n, K)+    if w_shrink.max() <= 0:+        top3 = np.zeros((n, W.shape[0]), dtype=bool)+    else:+        k3 = np.argpartition(-Aw, min(TOP_MODS, W.shape[0]) - 1, axis=1)[:, :TOP_MODS]+        top3 = np.zeros((n, W.shape[0]), dtype=bool)+        top3[np.arange(n)[:, None], k3] = True+        top3 &= w_shrink[None, :] > 0+    n_ev_open = n_ev_close = 0+    if not opening_off:+        Wabs = np.abs(W)+        top_load = np.argsort(-Wabs, axis=1)[:, :TOP_LOAD]  # (K, 50) hvg-relative+        keep_mod = np.where(w_shrink > 0)[0]+        pool_hvg = np.unique(top_load[keep_mod].ravel())    # sorted hvg-relative+        pool_abs = hvg[pool_hvg]                            # panel gene indices+        Pn = len(pool_hvg)+        T = np.zeros((Pn, W.shape[0]), dtype=bool)          # gene -> module membership+        for k in keep_mod:+            T[np.searchsorted(pool_hvg, top_load[k]), k] = True+        cand_mod = (top3.astype(np.float32) @ T.astype(np.float32).T) > 0  # (n, Pn)++        Xl_pool = np.asarray((last.X[:, pool_abs]).todense(), dtype=np.float32)  # (n_last, Pn)+        p_last = np.asarray(Wk @ (Xl_pool > 0), dtype=np.float32)                # (n, Pn)+        p_prev = np.zeros((n, Pn), dtype=np.float32)+        for q, sub_p in zip(qs, prev_subs):+            Xnz = np.asarray((prev.X[sub_p][:, pool_abs]).todense(), dtype=np.float32) > 0+            p_prev += (1.0 / len(qs)) * (q @ Xnz.astype(np.float32))+        dp = np.clip(p_last + factor * (p_last - p_prev), 0.0, 1.0) - p_last+        cand = cand_mod & (np.abs(dp) >= DP_MIN)++        Xpool_out = Xl_pool[rows]  # (n, Pn) values of the output cells at pool genes+        nz = Xpool_out > 0+        prob_open = np.clip(dp / np.maximum(1.0 - p_last, 1e-6), 0.0, 1.0)+        prob_close = np.clip(-dp / np.maximum(p_last, 1e-6), 0.0, 1.0)+        ev_open = cand & (dp > 0) & ~nz & (rng.random((n, Pn)).astype(np.float32) < prob_open)+        ev_close = cand & (dp < 0) & nz & (rng.random((n, Pn)).astype(np.float32) < prob_close)+        ev_open = cap_events(ev_open, prob_open, m_cap)+        ev_close = cap_events(ev_close, prob_close, m_cap)+        n_ev_open, n_ev_close = int(ev_open.sum()), int(ev_close.sum())++        # opening values: draw from the cell's latest-stage neighbour non-zero values of that gene+        for pi in np.where(ev_open.any(axis=0))[0]:+            js = np.where(ev_open[:, pi])[0]+            vals = Xl_pool[nb_idx[js], pi]           # (m, K) neighbour values+            nzm = vals > 0+            cnt = nzm.sum(axis=1)+            ok = cnt > 0+            r_draw = rng.integers(0, np.maximum(cnt, 1))+            cum = np.cumsum(nzm, axis=1)+            sel = nzm & (cum == (r_draw[:, None] + 1))+            drawn = vals[sel]+            Xpool_out[js[ok], pi] = drawn+        Xpool_out[ev_close] = 0.0+        # write pool changes back as a sparse correction (stays sparse, no lil fancy indexing)+        orig = Xl_pool[rows]+        diff = Xpool_out - orig+        jj, pi = np.nonzero(diff)+        if jj.size:+            corr = sparse.coo_matrix((diff[jj, pi], (jj, pool_abs[pi])), shape=(n, len(genes))).tocsr()+        else:+            corr = sparse.csr_matrix((n, len(genes)))+        Xout = (last.X[rows] + corr).tocsr()+        _log(f"opening/closing: pool {Pn} genes, open events {n_ev_open} "+             f"({n_ev_open / max(n, 1):.1f}/cell), close events {n_ev_close} ({n_ev_close / max(n, 1):.1f}/cell)", t0)+    else:+        Xout = last.X[rows]++    # ---- add the step to non-zero entries only, clip at 0 ----+    if args.proj in ("full", "hybrid"):+        s_full = np.asarray((Wk @ last.X).todense(), dtype=np.float32)+        d_full = s_full - np.mean(A_full, axis=0)+        del s_full, A_full+        d_full[:, ~mask] = 0.0+        if args.proj == "full":+            step_mat = factor * d_full+        else:+            step_mat = res_lam * ratio * d_full+            step_mat[:, hvg] += step+            del d_full+    else:+        step_mat = np.zeros((n, len(genes)), dtype=np.float32)+        step_mat[:, hvg] = step+    X = Xout.toarray()+    step_mat *= X > 0+    X += step_mat     np.maximum(X, 0.0, out=X)++    # ---- diagnostics (mechanism evidence) ----+    moved = np.abs(step_mat).sum(axis=1) > 0+    d_disp = Cbar * w_shrink  # per-cell module coefficients actually used+    _log(f"diag: cells touched by step {moved.mean() * 100:.1f}%  step |median| "+         f"{np.median(np.abs(step_mat)):.4f} max {np.abs(step_mat).max():.3f}", t0)+    _log(f"diag: zero-rate last {1 - last.X[rows].nnz / (n * len(genes)):.4f} -> out "+         f"{1 - np.count_nonzero(X) / X.size:.4f}", t0)+    try:+        ct = np.asarray(last.obs["celltype"].astype(str))+        ct_rows = ct[rows]+        # within-type share of displacement variance (labels for diagnostics only)+        tot = ((d_disp - d_disp.mean(axis=0)) ** 2).sum()+        within_t = 0.0+        for t in np.unique(ct_rows):+            m_t = ct_rows == t+            if m_t.sum() > 1:+                within_t += ((d_disp[m_t] - d_disp[m_t].mean(axis=0)) ** 2).sum()+        _log(f"diag: within-type share of displacement variance {within_t / max(tot, 1e-12):.3f}", t0)+    except Exception as e:  # labels are optional+        _log(f"diag: celltype diagnostics skipped ({e})", t0)+     write_prediction(X, genes, args.out, seed=args.seed)+    _log("done", t0)   if __name__ == "__main__":

调研来源?调研员查到并用到的知识条目和文献检索结果（只列标题和编号）。

用到的知识库条目

分析结果?分析员写的 ANALYSIS.json：改了什么、各组分数怎么变、假设是否成立、经验和下一步建议。

对话摘要?每个角色和大模型对话的统计：轮数、工具调用、用时、token 数和最后的回答摘录；原始记录只给路径。

只给统计和最后回答的摘录；完整对话请到原始记录位置里列出的文件看。

分析员

by_group": {
  "cell_state": "变好 +13.02（49.21→62.23）。主要来自 X3 的 mmd_u：0.03281→0.02075，得分 15.02→20.37（+5.35）。proxy10 的 mmd_u 0.04512→0.04135，得分 +1.02，这部分来自关闭生长重采样，与位移机制无关。",
  "direction": "变好 +2.84（50.90→53.74），但主要来自 proxy10：de_direction 0.098→0.201，得分 +1.39，原因是关闭生长重采样、回到照抄。X3 的 de_direction 为 -0.0228→0.0159，得分 +0.37，在噪声内。",
  "de_recovery": "变坏 -1.74（50.79→49.05），在 T1 约 2 分的噪声边缘。X3 的 de_score 0.039→-0.0649，得分 -0.88，已低于地板 12.50。proxy10 的 de_score -0.013→0.042，得分 +0.45。",
  "covariation": "噪声内 +0.12。proxy10 的 variogram 0.001332→0.001252，得分 +0.41；X3 的 variogram 0.001415→0.001456，得分 -0.17。",
  "board": "榜分 50.33→54.53（+4.20），超出噪声。其中 X3 52.94→55.32（+4.67），对应于 50.65→55.32。proxy10 49.67→52.94（+3.27），与 Engineer 的预期一致：proxy10 只是从生长加权照抄回到均匀照抄。Engineer 报告的 X3 54.94 与实测 55.32 一致。"
 },
 "hypothesis_supported": "unclear",
 "failure_reason": null,
 "lessons": [
  "两输入视图上，OT 精确祖先位移投影到输入 PCA 模块后放大外推（λ=3），可在 X3 上把 mmd_u 从 0.0328 降到 0.0208（+5.35 分）。但它几乎不改善 DE 两项，X3 的 de_score 还降到地板以下。这类位移外推应当作 cell_state 杠杆，不应当作 DE 杠杆。",
  "关闭 g^dt 生长重采样、改为均匀抽样，可让单输入的 proxy10 从 49.67 回到照抄水平 52.94。这部分收益与 OT 机制无关，判断机制时必须把它从榜分增量中剔除：消融时 proxy10 的变化为 +0.00。",
  "在 addnz 解码、λ≤1 的条件下，按检测率开启/关闭基因（M=30）两次对照都拉低 mmd_u 和 variogram，de_score 没有稳定收益（据 Engineer 的查分记录）。在幅度相近的设定下，不宜再把它作为 DE 的主要机制。",
  "单次全细胞精确耦合列优于两份 80% 子集耦合加 kNN 映射。在精确模式下，SNR 收缩无法辨识（w≡1），PLAN 中基于子集一致性的收缩没有实际起作用。",
  "λ 与 reslam 全部在 X3 的 A 半上用 20 次查分选出，选参尺子与评分尺子相同。实测 X3 +4.67 与开发时的估计一致，但在别的 seed、B 半或 final 外推上是否稳健，尚无证据。"
 ],
 "next_suggestions": [
  "【cell_state，X3】在 λ∈{3,4,5} × reslam∈{0,0.5,1} 的小网格上，按 seed 0/1/2 复查。同时监控 variogram（X3 已 -0.17）和零值比例，确认 mmd_u 的收益没有在更大 λ 下转为过冲；若 λ=4 不再改善就固定为 3。",
  "【de_recovery / direction，X3】DE 两项只看 dp 的基因排序和符号。可比较 reslam=0 与 1，以及只在全部掩码基因上用原始位移、不经 PCA 投影这两种做法，判断 X3 的 de_score 下降（0.039→-0.065）来自哪一项。若来自 PCA 投影放大，可只在 HVG 外基因上保留原始步，使投影步主要服务分布位置。",
  "【全部分组，可信度】用“零发育”对照（两个输入取自同一阶段的伪目标）检验 λ=3 下的位移外推：确认它在没有真实发育时不会损害 mmd_u。若有损害，说明收益依赖外推方向碰巧与 X3 真值同向，需要降低 λ 或加入方向一致性门控。"
 ],
 "family_id": "local_ot",
 "mechanism_active": "yes"
}

工程师

**完成总结**

**做了什么**：按 PLAN 重写父节点 2（ot_moscot）：均匀抽样输出（关生长重采样）、moscot 全细胞耦合 + 精确祖先重心、位移投到 30 个输入 PCA 模块并阻尼外推（λ·dt_out/dt_in）、叠加弱原始全基因残差（0.5×）、addnz 解码；检测率开关基因机制（机制 B）和子集 SNR 收缩均已实现，但经对照证伪/被更优方案替代后默认关闭（`--mcap`/`--subset` 可启用），三个 `--ablate` 对照齐备。已验证：确定性（seed 0 重跑逐位相同）、伪装视图（+1 天平移、乱序 manifest）输出逐位相同、两视图 vec-check ok、X3 全程 ~13 s。

**查分结果**（20/20 额度用尽，X3 A 半）：最终配置（exact hybrid λ3+0.5·raw）**54.94**（mmd_u 20.03 / de_dir 12.68 / vario 10.50 / de_score 11.73），对比父节点复跑 ≈50.6、均匀照抄 48.60；λ 序列 48.9→49.3→50.1→52.3→54.9 单调向好。机制 B（M30 vs M0 两次对照）与生长重采样（3 次对照）均降分，已如实记入 METHOD.md。proxy10 走均匀照抄退路（预期 ≈52.7，额度耗尽未查分）。估计节点分 ≈ (52.7 + 2×54.9)/3 ≈ 54，较父节点 50.33 提升 ~+3.7（B 半噪声 ±2 内仍应为正）。

**下一步最值得试**：① λ=4 / reslam 微调（趋势未饱和，但需查分验证方向）；② de_score 仍低于地板（11.7/12.5）——位移只覆盖 HVG+弱残差，可试按模块载荷加权扩展到更多基因的定向步长；③ 单输入退路可借节点 3 思路（型内成熟度选择/组成重加权）把 proxy10 从 52.7 再往上推——那是当前分数最大的未开发杠杆。

调研员

_last − p_next)/p_last 置 0。\n- 每个细胞开启与关闭各设上限 M，初值 30，范围 {0, 15, 30, 60}。开启或关闭之后，再在非零元上加 δ_j，clip 到 0。\n\n(4) 单输入阶段的退路（proxy10）：没有第二阶段时无法估计耦合，输出最新阶段的均匀抽样，即与 copy_last 同流程。这一行为由方法本身是否需要两阶段决定，不由视图身份决定，在伪装视图上结果相同。proxy10 的预期只是从 49.67 回到照抄水平，不能用它选择两输入外推的参数。\n\n(5) final 视图（官方 E8.5、E9.5 → E10.5，间隔 1 天，全胚）：dt_out/dt_in = 1。只在两阶段都测到的基因上做步长和开启 / 关闭；外部阶段补齐的基因不动。\n\n(6) 查分顺序，共 20 次查分，X3 权重为 2，以 X3 为主：\n1. 小样本（每阶段 1000 个细胞）冒烟测试，确认运行时间和内存；\n2. 全量默认配置（λ=0.5、K=30、M=30）+ 三个对照：λ=0、M=0（禁止开启 / 关闭）、growth_on，共约 4 次；\n3. λ 与 M 的小网格，约 4 次；\n4. 最佳配置及其 M=0 对照，在 seed 1、2 上各复查，约 4 次；\n5. 留余量。\n只有在提升 ≥2 分，且在三个 seed 上方向一致时才算进步。",
 "expected_groups": ["cell_state", "direction", "de_recovery"],
 "risks": "1. X3 可能不是两输入视图。若 stderr 里没有出现 'ot: converged'，说明 X3 只走退路分支，本机制在两把尺子上都无法检验。Engineer 第一次运行就要确认，并在 METHOD.md 如实写明。\n2. 配对错误或跨平台差异会让 Δp 系统性偏移，例如两阶段测序深度不同，导致检测率整体变化。检查办法：看开启和关闭是否高度集中在少数几个“全局”模块，或者开启数在各细胞类型之间几乎相同。若是这样，在 Δp 中减去全细胞的平均 Δp（全局深度校正）后再判断。\n3. X3 的 de_direction 已为负。若 λ>0 之后 de_direction 继续下降，说明外推方向错误，应把 λ 降到 0.25，或报告失败。\n4. 开启过多会损害 mmd_u 和 variogram。需要监控零值比例：输出与最新阶段的零值比例相差应 <2 个百分点。\n5. NMF 的运行时间：把 NMF 限制在 HVG 和 ≤4000 个细胞上。\n6. 子集耦合不稳定会让 w_k 全部接近 0，机制退化为照抄。需要报告 w_k 的分布。",
 "family_id": "local_ot",
 "mechanism": "1. 用 OT 耦合给每个细胞估计祖先，从而得到局部位移和局部检测率变化。\n2. 把位移投到少数数据驱动的程序模块上，按子集耦合之间的一致性收缩并阻尼后延续一步。\n3. 在受模块支持的基因上，按局部检测率趋势抽样开启或关闭基因，保留每个细胞自己的残差。",
 "vs_constant_shift": "步长与开启 / 关闭都来自每个细胞自己的祖先耦合和近邻，同一细胞类型内部各不相同。改动只发生在少数受支持的模块和基因上，以稀疏方式开关，不是给每个类型或全体细胞加同一个向量。生长重采样关闭，也不调整组成比例。",
 "mechanism_evidence": "需要报告：\n1. 位移的型内方差占总方差的比例（用输入阶段的标签只做诊断）。若比例接近 0，说明退化为每型常向量。\n2. 每个细胞开启和关闭的基因数分布，以及这些基因所属模块的分布、各模块的 w_k。\n3. 输出相对最新阶段的零值比例、方差、离群程度的变化。\n4. 实际被改动的细胞比例。\n5. 与 M=0、λ=0 两个对照相比，X3 上四项指标的 score_parts 各自如何变化：开启 / 关闭应主要影响 de_score、de_direction，模块位移应影响 mmd_u。",
 "mechanism_off_control": "同一个 run.py：\n- `--ablate displacement`：λ=0，没有步长也不开关基因，等于均匀照抄最新阶段。\n- `--ablate opening`：M=0，只在非零元上加模块投影位移。\n- `--ablate growth_on`：恢复 g^dt 重采样。\n预期：完整版与 λ=0 的输出不同，改动集中在少数模块。完整版与 M=0 相比，差别主要体现在 DE 两项。若关闭机制后输出与打开时相同，就判定机制没有运行。"}
```

审查员

{"verdict": "pass", "reasons": ["1 越界读取：未发现问题。run.py 和 g37_common.py 的所有读取都经过 view_io 的 load_manifest、panel_genes、read_stage、covered_mask 完成。stage_pair 中 inputs_by_time(include_external=True) 只读取 manifest 里列出的 external 条目（X3 的 qiu2024_heart_early E8.75 和 E9.0），这些都在视图之内。代码里没有绝对路径、'..'、data/raw 或评分器路径，也没有联网。proliferation_markers 和 apoptosis_markers 取自已安装的 moscot 包内置列表，不需要下载。", "2 硬编码目标统计量：未发现问题。run.py 顶部的常量（N_HVG、N_PCS、EPSILON、LAMBDA=3、RES_LAMBDA=0.5、M_CAP 等）都是超参数。_gen_logistic 的参数来自 moscot/WOT 的通用生长先验。代码里没有比例表、细胞数、类型名或阶段名。输出细胞数由 target_n_cells(manifest) 在运行时决定。", "3 钻评分器漏洞：未发现问题。输出细胞是从最新阶段均匀无放回抽样，weighted_rows 里的 np.sort 只保留原始行序，沿用父节点做法，没有为按行号抽样的指标特意排序。位移外推 λ=3 作用于每个细胞自己的 OT 祖先位移，属于生物学意义上的外推，不是整体缩放、平均细胞或行序技巧。不涉及 T2 坐标。", "4 结果不确定：未发现问题。随机性全部由 --seed 驱动：np.random.default_rng(args.seed) 和 default_rng(seed+r)；embed 中的 PCA 设置了 random_state=seed；sc.tl.score_genes 设置了 random_state=seed。JAX_PLATFORMS 固定为 cpu，Sinkhorn 是确定性求解。METHOD.md 报告 seed 0 两次运行逐字节相同。", "5 保留阶段与保留基因型：未发现问题。X3 的输入是 E8.75 和 E9.0（含 external），final 的输入按 METHOD 为 E8.5 和 E9.5，都不在 T1 的禁窗（E10.5、E12.5、9.5<E≤13.5）内；E10.5 只作为 target 的时间用于计算 dt_out，代码没有读取它的数据。没有使用 Gata4 或 β-catenin 数据。生长先验在 METHOD.md 知识来源一节注明了 Schiebinger et al. 2019 (DOI 10.1016/j.cell.2019.01.006)；事件表注明为未使用。", "6 只能在替代评测上成立：未发现问题。分支只看输入阶段个数（stage_pair 中 len(stages)<2 时走均匀照抄退路），只用时间差 dt_in 和 dt_out，没有写死阶段名或文件名。celltype 标签只在 try/except 里用于诊断输出，没有标签也不会崩。METHOD 报告伪装视图的输出逐位相同。有一个风险但不构成 fail：exact 全细胞耦合加上稠密的 (n, 32285) 矩阵（s_full、A_full、prev.X.todense）在 final 视图上的内存没有实跑验证，目前只有约 <8 GB 的预估。", "7 重新引入被合规删除的内容：未发现问题。PLAN.json 中没有 [compliance: removed] 标记或 _compliance 字段，METHOD.md 和代码里也没有保留阶段的测量值、比例或细胞数。", "附注（不影响判定）：proxy10 走单输入退路，机制只在 X3 上起作用（消融后 X3 -6.73，proxy10 +0.00）。λ=3 和 RES_LAMBDA=0.5 是用 X3 的 A 半选出来的，正式分用 B 半，可能存在选择偏差（噪声约 ±2）。最终配置没有做多 seed 复查。"]}