总览 · ← 返回运行 20261001-114429-search-t1-g18-continue

节点 n12

移植节点9的类内成熟度加权抽样（周期分深度残差化稳健z，w=exp(-0.5z)，E-S无放回）为主干；新增跨阶段伪批量位移分支（Qiu E9.0/两官方阶段，批校+EB收缩），实测有害，默认α=0关闭。

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

来自 programs.git 04feafd2d1:solution/METHOD.md

调研员的计划

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

对比：父节点版本 5ba0ea2aa1。改动的文件：solution/METHOD.md +49 −38、solution/maturity.py +109 −0、solution/run.py +59 −26、solution/shift.py +137 −0

diff --git a/solution/METHOD.md b/solution/METHOD.mdindex c711d1c..694e38c 100644--- a/solution/METHOD.md+++ b/solution/METHOD.md@@ -1,47 +1,58 @@-增殖驱动的组成重加权（负结果）：类型权重×(1+γ·增殖z)，γ=0 复现父本 55.97，γ=±0.25/0.35/1.0 与 5118 细胞均降分。+移植节点9的类内成熟度加权抽样（周期分深度残差化稳健z，w=exp(-0.5z)，E-S无放回）为主干；新增跨阶段伪批量位移分支（Qiu E9.0/两官方阶段，批校+EB收缩），实测有害，默认α=0关闭。  ## 做了什么 -父本（节点 2，55.97）只重加权组成、不动表达，最弱组是 de_recovery 53.06。本节点按 PLAN 试了两条只依赖输入阶段数据的新杠杆：--1. **逐类型分化轴位移（PLAN 步骤 1-3，未采用）**：在 E8.5 每个类型内，用 `prior/` 里的细胞周期基因集-   （reactome/go/msigdb 的 `*cell cycle*|*replication*|*mitotic*|*E2F*|*G2M*`，命中面板 2091 个基因）算每细胞增殖 z 分数，-   再取类型内 1000 HVG 的前 5 个 PC，选 |corr(PC, s)| 最大者作时间轴，δ_c = 晚端三分之一 − 早端三分之一。-   **实测：18 个类型无一通过 |corr| >= 0.15**（周期信号不在类型内前 5 个主成分里），全部 δ_c = 0，该分支等价于父本，故已从代码中移除。-   诊断（未查分）：直接用增殖分数三等分（不经 PCA）时，δ 的前 20 个基因是 Rpl/Rps/Tmsb10/Malat1，且 IFT-CM 的-   "晚端"下调 Myh6/Ttn/Actc1/Ankrd1 —— 与心肌成熟方向相反。说明这条轴在 E8.5 单快照上主要是核糖体/深度等技术变异，不是发育时间。-2. **增殖驱动的组成增长（已实现，默认关闭）**：`growth.py` 里 `type_growth` 给每个类型一个增殖 z 分数（对阶段中位数居中、阶段 std 归一，-   细胞数 < 30 的类型记 0），`growth_sample` 把父本的解剖权重（heart ×1.6、edge ×0.25、丢 Neural Tube）乘上-   `1 + γ·clip(g, -1, 1)`，其余抽样逻辑（类型顺序、`largest_remainder`、`take`、RNG 调用序列）与 `heart_reweight` 逐调用一致。--## 查分结果（T1:val proxy，seed 0，单次）--| 配置 | 榜分 | de_recovery | direction | cell_state | covariation |-|---|---|---|---|---|---|-| γ=0（=父本，逐位相同） | **55.97** | 53.06 | 58.56 | 56.90 | 54.97 |-| γ=+0.35 | 55.29 | 54.17 | 57.20 | 55.49 | 54.03 |-| γ=-0.25 | 54.96 | 53.06 | 58.06 | 55.12 | 53.18 |-| γ=+1.0 | 52.68 | 53.06 | 55.36 | 50.18 | 52.63 |-| N_CELLS=5118（max_cells），γ=0 | 55.12 | 53.06 | 58.38 | 55.16 | 53.57 |--γ 的两个方向都单调劣于 γ=0，说明父本手调的心脏权重已经吃掉了增殖能提供的组成信息（增殖排序本身是合理的：-NCC +0.89、aSHF +0.47、Endothelium +0.27 高，Foregut -0.50、Blood -0.68 低）。加大细胞数也降分。-所以**提交的默认配置 γ=0、N_CELLS=4000，输出与父本逐位相同**（已用 `np.array_equal` 核对）；机制保留在代码里，`VEC_GAMMA` 可复现上表。+父节点 4（55.97，=节点 2 逐位复现）最弱组 de_recovery 53.06。本节点按 PLAN 做了两件事：++1. **成熟度加权抽样（提交主干，`maturity.py`）**：重实现节点 9 的已证部件。解剖组成与父本完全一致+   （heart ×1.6、edge ×0.25、丢 Neural Tube、4000 细胞、`largest_remainder` 配额），但类内不再均匀抽样：+   每细胞算 `prior/` gmt 细胞周期基因集的 z 分数 s（全体细胞上标准化，稀疏 matvec），类型内对 `[1, log1p(库大小)]`+   闭式 OLS 取残差，median/1.4826·MAD 稳健 z，clip±2，权重 `w=exp(-0.5·z)`，Efraimidis-Spirakis+   （key=u^(1/w) 取 top-k）加权无放回抽样；池不足时有放回加权补齐。<30 细胞的类型退化为均匀。+2. **跨阶段伪批量位移（`shift.py`，默认关闭）**：视图有第二输入阶段时（proxy2 的 Qiu E9.0 外部输入 / final 的官方+   E8.5→E9.5），逐类型 δ = mean(参考阶段, 匹配类型) − mean(最后官方阶段, 类型)；跨技术批效应按 scArches+   位置漂移假设做逐基因跨类型 median 中心化；EB 收缩 δ·t²/(t²+4)（se 由两 pseudobulk 方差算）；|δ| 99 分位截断；+   周期基因与未覆盖基因 δ=0；外部阶段类型匹配用基因中位数中心化 pseudobulk 的相关（≥0.3，运行期计算，无硬编码标签映射，+   proxy2 实测匹配：Endocardial→Endothelium、FHF→SV-CM、SHF→OFT/RV-CM，r≈0.74-0.75）；应用+   x ← x + α·r·δ，r=(目标时间−最后输入时间)/(参考−最后输入)，clip≥0。α=0 时整个分支（含外部文件读取）被跳过，+   proxy/proxy2 输出逐位相同。++## 查分结果（seed 0 单次，除注明外）++| 配置 | 视图 | 榜分 | de_rec | direction | cell_state | covar |+|---|---|---|---|---|---|---|+| 提交默认（θ=0.5, α=0） | proxy | 57.76 | 53.06 | 60.92 | 60.60 | 55.41 |+| 同上 seed 1 / seed 2 | proxy | 57.61 / 57.08 | | | | |+| **3-seed 均值** | proxy | **57.48** | 52.7 | 60.8 | 60.1 | 55.4 |+| 位移 α=0.5 / α=1.0 | proxy2 | 56.78 / 55.91 | 52.0 / 52.0 | 60.86 / 60.23 | 59.5 / 58.2 | 53.57 / 51.97 |+| θ=1.0 / θ=0.25 / heart×2.0 | proxy | 57.67 / 56.94 / 57.16 | | | | |++**位移分支的负结果**：α 从 0→0.5→1.0 在 proxy2 上单调降分（57.76→56.78→55.91），四个组全部下降，de_recovery+53.06→52.00。Qiu 与官方数据的批效应不是纯基因级位置漂移，median 中心化后残留的类型×技术交互盖过了 0.5 天的+时间信号；且 δ 幅度本身很小（α=1 时受影响条目 mean|Δ|≈0.06）。按 PLAN 的预设放弃规则，提交 α=0。  ## 验证过 / 没验证 -- 验证过：`γ=0` 与父本预测逐位相同；`vec-check` 在 seed 0 和 seed 3 上均 ok；`run.py` 只按 `len(inputs_by_time(manifest))`-  读最后一个输入阶段，proxy/final 同一段代码，无阶段名分支、无硬编码统计量；prior 基因集缺失或命中 < 20 时 γ 自动置 0。-- 未验证：所有结论只基于 seed 0 单次查分（低于 PLAN 要求的 3-seed 均值），γ=+0.35 的 de_recovery +1.11 可能只是噪声；-  PLAN 步骤 4 的两阶段伪批量融合分支没有写进提交版本（单快照轴在 proxy 上就没有信号，融合无意义）；final 视图上未运行。-- 运行 ~40 s、峰值内存 ~1.5 GB（`type_growth` 只在 4000 细胞 × 2091 基因的子块上做稠密计算）。+- 验证过：默认配置在 proxy 与 proxy2 上输出逐位相同（α=0 分支跳过，`np.array_equal`）；两视图 `vec-check` ok；+  seed 0/1/2 proxy 3-seed 均值 57.48（≈节点 9 的 57.42，RNG 序列不同导致的实现差异）；θ=0.5 在 seed 0 上优于+  θ=0.25/1.0 与 heart×2.0；无阶段名分支、无硬编码统计量、`np.random.default_rng(seed)` 确定性。+- 未验证：final 视图（两官方阶段的 δ 无批混杂，理论上比 Qiu 干净，但 proxy 上无法测，α=0 提交即不启用）；+  θ/权重的 3-seed 复核（时间不够，只有 seed 0 单次）；位移分支的 k、截断分位等超参未搜（首个 α 已触发放弃规则）。  ## 下一步最值得试 -1. **表达侧必须换信号源**：单快照的增殖轴在 E8.5 上是技术轴。final 视图有 E8.5+E9.5 两个阶段，-   同名类型的两阶段伪批量差 δ2_c 是唯一有真值方向的位移信号（配 T1-02 式逐基因 t 统计收缩，k≈4），-   在 final 上直接可算；proxy 上无法验证，只能靠 final 的一次查分。-2. **组成侧改成网格搜索父本的三个手调数**（heart 权重、edge 权重、是否丢 Neural Tube），-   这是唯一被证明有效的杠杆（copy_last 49.77 → 55.97），比任何表达位移的期望收益都高；用 3-seed 均值防噪声。-3. covariation 54.97 与 mmd 是第二弱项：可试在抽出的真实细胞上做同类型内的轻微 jitter（保留协方差结构的量级），-   而不是加独立噪声。+1. **final 视图上启用官方两阶段位移**（α≈0.25-0.5）：E8.5→E9.5 同数据集 δ 无批混杂，是唯一未被否定的表达位移+   信号；`VEC_ALPHA` 环境变量即可开启，代码路径已就绪（`label_match=True` 走同名类型）。+2. de_recovery 仍是弱组（52.5）：θ=1.0 时 de_score 反降（0.111→0.056），说明"更成熟的细胞"不等于"更对的 DE"；+   试把 DE 方向信息（两阶段同名类型差异基因）直接加权到类内抽样概率，而不是平移表达。+3. 组成侧 heart/edge 权重在 θ=0.5 主干上重新网格（本节点只测了 heart×2.0 单点，57.16）。++## 来源披露++- scArches 位置漂移批校正假设：doi:10.1038/s41587-021-01001-7（用于 shift.py 的逐基因 median 中心化；实测在+  Qiu↔官方跨技术上不成立，故默认关闭）。+- Efraimidis-Spirakis 加权无放回抽样：Efraimidis & Spirakis 2006, Information Processing Letters（key=u^(1/w)）。+- 细胞周期基因集：view `prior/` 的 reactome/go/msigdb gmt 中名称匹配 `cell cycle|dna replication|mitotic|E2F|G2M|+  S phase|M phase` 的集合（运行期读取，命中约 2091 基因）。+- 解剖权重（heart×1.6、edge×0.25、丢 Neural Tube）与 4000 细胞：继承 `src.task1_temporal.reweight`（run2 冠军+  heart_jcf_peri 的手调结果），本节点未改。diff --git a/solution/maturity.py b/solution/maturity.pynew file mode 100644index 0000000..399b580--- /dev/null+++ b/solution/maturity.py@@ -0,0 +1,109 @@+"""Cell-cycle-exit (maturity) weighted sampling within each cell type.++Port of the proven node-9 component: cells of the last observed stage keep the+parent's anatomical composition, but within each type cells are drawn without+replacement with weight ``exp(-theta * z)``, where ``z`` is the robust z score+(inside the type) of the per-cell cell-cycle score after residualisation on+sequencing depth.  High-cycle (immature) cells are down-weighted, mature ones+up-weighted: over E8.5 -> E9.5 (and E9.5 -> E10.5) populations exit the cycle,+so the sampled cloud sits closer to the target state.++Everything (cycle gene set, score, depth, z) is computed from the view at run+time; no stage statistic is hard-coded.+"""++from __future__ import annotations++import numpy as np+from scipy import sparse++from growth import cell_cycle_genes++MIN_TYPOOLS = 30+Z_CLIP = 2.0+++def cycle_score(X, cc_mask: np.ndarray) -> np.ndarray:+    """Per-cell mean z score over cell-cycle genes (all cells, sparse matvec)."""+    cc_idx = np.flatnonzero(cc_mask)+    if cc_idx.size == 0:+        return np.zeros(X.shape[0], dtype=np.float64)+    A = X[:, cc_idx]+    mu = np.asarray(A.mean(axis=0)).ravel().astype(np.float64)+    sd = np.asarray(A.power(2).mean(axis=0)).ravel().astype(np.float64)+    sd = np.sqrt(np.maximum(sd - mu**2, 0.0))+    keep = sd > 1e-6+    if int(keep.sum()) < 20:+        return np.zeros(X.shape[0], dtype=np.float64)+    inv = np.zeros_like(sd)+    inv[keep] = 1.0 / sd[keep]+    s = np.asarray(A.dot(inv)).ravel() / float(keep.sum())+    s -= float((mu[keep] * inv[keep]).sum()) / float(keep.sum())+    return s+++def maturity_weights(X, labels: np.ndarray, cc_mask: np.ndarray, theta: float = 0.5) -> np.ndarray:+    """Per-cell sampling weight exp(-theta * z); z = depth-residualised robust cycle z inside type."""+    n = X.shape[0]+    s = cycle_score(X, cc_mask)+    lib = np.asarray(X.sum(axis=1)).ravel().astype(np.float64)+    d = np.log1p(np.maximum(lib, 0.0))+    w = np.ones(n, dtype=np.float64)+    for t in np.unique(labels):+        m = np.flatnonzero(labels == t)+        if m.size < MIN_TYPOOLS:+            continue+        A = np.column_stack([np.ones(m.size), d[m] - d[m].mean()])+        coef, *_ = np.linalg.lstsq(A, s[m], rcond=None)+        res = s[m] - A @ coef+        med = np.median(res)+        mad = np.median(np.abs(res - med))+        scale = 1.4826 * mad+        if not np.isfinite(scale) or scale < 1e-9:+            continue+        z = np.clip((res - med) / scale, -Z_CLIP, Z_CLIP)+        w[m] = np.exp(-theta * z)+    return w+++def weighted_take(X, pool: np.ndarray, w_pool: np.ndarray, k: int, rng: np.random.Generator):+    """Efraimidis-Spirakis weighted sampling without replacement (top-k of u^(1/w))."""+    if pool.size == 0 or k <= 0:+        return sparse.csr_matrix((0, X.shape[1]), dtype=np.float32)+    if k >= pool.size:+        head = pool+        rest = k - pool.size+        if rest > 0:+            p = w_pool / w_pool.sum()+            head = np.concatenate([pool, rng.choice(pool, size=rest, replace=True, p=p)])+        return X[head]+    keys = rng.random(pool.size) ** (1.0 / np.maximum(w_pool, 1e-12))+    top = np.argpartition(-keys, k)[:k]+    return X[pool[top]]+++def maturity_sample(X, labels: np.ndarray, n_cells: int, cc_mask: np.ndarray, seed: int = 0, theta: float = 0.5, hw: float = 1.6, ew: float = 0.25):+    """Resample with the parent's anatomical weights; within type sample by maturity."""+    from src.task1_temporal.reweight import (+        DROP_TYPES,+        largest_remainder,+        type_weights,+    )++    rng = np.random.default_rng(seed)+    w = maturity_weights(X, labels, cc_mask, theta=theta)+    types = [str(t) for t in np.unique(labels) if str(t) not in DROP_TYPES]+    counts = np.array([(labels == t).sum() for t in types], dtype=np.float64)+    alloc = largest_remainder(counts * type_weights(types, hw, ew), n_cells)+    blocks = []+    row_types: list[str] = []+    for t, k in zip(types, alloc):+        if k <= 0:+            continue+        pool = np.flatnonzero(labels == t)+        blk = weighted_take(X, pool, w[pool], int(k), rng)+        blocks.append(blk)+        row_types.extend([t] * blk.shape[0])+    out = sparse.vstack(blocks, format="csr").astype(np.float32)+    out.eliminate_zeros()+    return out, np.asarray(row_types, dtype=object)diff --git a/solution/run.py b/solution/run.pyindex 5161138..af7a42f 100644--- a/solution/run.py+++ b/solution/run.py@@ -1,17 +1,17 @@ #!/usr/bin/env python3-"""heart_jcf_peri + proliferation-driven composition growth.--Composition: the latest input stage is resampled by cell type with the parent's-anatomical weights (heart x1.6, dissection edges x0.25, neural tube dropped),-and those weights are additionally scaled by each type's own proliferation-level measured on the input stage: a type whose cells cycle harder than the-stage median gains relative mass over the interval, a slower one loses it.-The per-type growth factor is computed from the data at run time (cell cycle-gene sets come from the view's ``prior/``), never hard-coded, so the same code-runs on the proxy view (E8.5 -> E9.5) and the final view (E9.5 -> E10.5).-Expression is not shifted: cells are real cells of the last observed stage.--``gamma = 0`` reproduces the parent (heart_jcf_peri) bit for bit.+"""Anatomical reweighting + maturity-weighted sampling + cross-stage shift.++Composition: the parent's anatomical weights (heart x1.6, edges x0.25, neural+tube dropped) on the last official input stage.  Within each type, cells are+drawn without replacement with weight exp(-0.5 z), z = robust, depth-+residualised cell-cycle z (maturity.py): the cloud moves toward cycle-exit,+the state populations reach by the target stage.++Expression: when the view provides a second input stage (proxy2: external+Qiu E9.0; final: official E8.5), each type's cells are shifted along a+batch-centred, EB-shrunk pseudobulk delta times alpha * r (shift.py).  With a+single input stage (proxy) the shift branch is off and the code reduces to+composition + maturity sampling.  alpha = 0 also reduces exactly to that. """  from __future__ import annotations@@ -19,8 +19,13 @@ from __future__ import annotations import argparse import os +import numpy as np+ from src.task1_temporal.view_io import (+    covered_mask,+    external_inputs,     inputs_by_time,+    is_external,     labels_of,     load_manifest,     panel_genes,@@ -29,11 +34,15 @@ from src.task1_temporal.view_io import (     write_prediction, ) -from growth import cell_cycle_genes, growth_sample, type_growth+from growth import cell_cycle_genes+from maturity import maturity_sample+from shift import apply_shift, compute_delta  N_CELLS = 4000-# measured on proxy seed 0: 0.0 -> 55.97, +0.35 -> 55.29, -0.25 -> 54.96, +1.0 -> 52.68-GAMMA = float(os.environ.get("VEC_GAMMA", "0.0"))+ALPHA = float(os.environ.get("VEC_ALPHA", "0.0"))+THETA = float(os.environ.get("VEC_THETA", "0.5"))+HW = float(os.environ.get("VEC_HW", "1.6"))+EW = float(os.environ.get("VEC_EW", "0.25")) VERBOSE = bool(os.environ.get("VEC_VERBOSE", ""))  @@ -46,20 +55,44 @@ def main() -> None:      manifest = load_manifest(args.data)     genes = panel_genes(args.data, manifest)-    stages = inputs_by_time(manifest)-    last = read_stage(args.data, stages[-1], genes)+    official = inputs_by_time(manifest)+    ext = external_inputs(manifest)+    last = read_stage(args.data, official[-1], genes)     labels = labels_of(last)      cc_mask = cell_cycle_genes(args.data, manifest, genes)-    growth = type_growth(last.X, labels, cc_mask, rng=__import__("numpy").random.default_rng(args.seed))-    if VERBOSE:-        print(f"cell cycle genes: {int(cc_mask.sum())}")-        for t, g in sorted(growth.items(), key=lambda kv: -kv[1]):-            print(f"  {t}: g={g:+.3f} -> w_mult={max(1.0 + GAMMA * max(min(g, 1.0), -1.0), 1e-3):.3f}")--    gamma = GAMMA if int(cc_mask.sum()) >= 20 and growth else 0.0     n = target_n_cells(manifest, N_CELLS)-    X = growth_sample(last.X, labels, n_cells=n, growth=growth, gamma=gamma, seed=args.seed)+    X, row_types = maturity_sample(last.X, labels, n, cc_mask, seed=args.seed, theta=THETA, hw=HW, ew=EW)++    # --- cross-stage shift (only with a second input stage; off at alpha=0) ---+    scale = 0.0+    delta: dict[str, np.ndarray] = {}+    if ALPHA == 0.0:+        pass+    elif len(official) >= 2:+        prev = read_stage(args.data, official[-2], genes)+        prev_labels = labels_of(prev)+        dt = float(official[-1]["time"] - official[-2]["time"])+        gap = float(manifest["target"]["time"] - official[-1]["time"])+        delta, log = compute_delta(last.X, labels, prev.X, prev_labels, cc_mask,+                                   covered=np.ones(len(genes), dtype=bool), label_match=True)+        scale = ALPHA * (gap / dt if dt > 0 else 0.0)+        if VERBOSE:+            print("official two-stage matches:", log)+    elif ext:+        e = ext[-1]+        ref = read_stage(args.data, e, genes, missing="zero")+        covered = covered_mask(args.data, e, genes)+        col = e.get("celltype_col", "celltype")+        ref_labels = ref.obs[col].astype(str).to_numpy()+        dt = float(e["time"] - official[-1]["time"])+        gap = float(manifest["target"]["time"] - official[-1]["time"])+        delta, log = compute_delta(ref.X, ref_labels, last.X, labels, cc_mask, covered=covered)+        scale = ALPHA * (gap / dt if dt > 0 else 0.0)+        if VERBOSE:+            print("external matches:", log, "covered genes:", int(covered.sum()))++    X = apply_shift(X, row_types, delta, scale)     write_prediction(X, genes, args.out, seed=args.seed)  diff --git a/solution/shift.py b/solution/shift.pynew file mode 100644index 0000000..7d9c218--- /dev/null+++ b/solution/shift.py@@ -0,0 +1,137 @@+"""Cross-stage pseudobulk expression shift for sampled cells.++When the view carries a second input stage (proxy2: external Qiu E9.0 heart;+final: official E8.5 -> E9.5), each cell type of the last official stage gets a+gene-level delta = mean(reference stage, matched type) - mean(last official+stage, type), and sampled cells are moved ``alpha * r * delta`` along it+(r = time extrapolation factor), clipped at 0.++Cross-dataset batch is removed under a location-shift assumption (scArches,+doi:10.1038/s41587-021-01001-7): per gene, deltas are centred on their median+across matched types, so only between-type contrasts survive.  Deltas are then+shrunk per gene by an empirical-Bayes t statistic (delta * t^2 / (t^2 + k)) and+clipped at the 99th percentile; cell-cycle genes get delta = 0 (the maturity+sampler already moves cells along that axis).++Type matching between the reference stage and the official last stage is by+label when the vocabularies overlap, otherwise by correlation of gene-centred+pseudobulk profiles (computed at run time; no hard-coded label map).+"""++from __future__ import annotations++import numpy as np+from scipy import sparse++MIN_CELLS = 30+K_SHRINK = 4.0+CLIP_Q = 0.99+MIN_CORR = 0.3+++def _pseudobulk(X, labels, types):+    out = np.zeros((len(types), X.shape[1]), dtype=np.float64)+    for i, t in enumerate(types):+        m = labels == t+        if m.sum():+            out[i] = np.asarray(X[m].mean(axis=0)).ravel()+    return out+++def match_types(pb_ref, ref_types, pb_off, off_types, covered):+    """ref index -> official index by best correlation on covered, gene-centred profiles."""+    A = pb_ref[:, covered]+    B = pb_off[:, covered]+    A = A - np.median(A, axis=1, keepdims=True)+    B = B - np.median(B, axis=1, keepdims=True)+    A = A / (np.linalg.norm(A, axis=1, keepdims=True) + 1e-12)+    B = B / (np.linalg.norm(B, axis=1, keepdims=True) + 1e-12)+    C = A @ B.T+    pairs = []+    for i in range(len(ref_types)):+        j = int(np.argmax(C[i]))+        if C[i, j] >= MIN_CORR:+            pairs.append((i, j, float(C[i, j])))+    return pairs+++def compute_delta(+    X_ref,+    ref_labels,+    X_last,+    last_labels,+    cc_mask: np.ndarray,+    covered: np.ndarray | None = None,+    label_match: bool = False,+) -> tuple[dict[str, np.ndarray], list[tuple[str, str, float]]]:+    """Per official type delta vector (full panel, float64) + match log."""+    n_genes = X_last.shape[1]+    if covered is None:+        covered = np.ones(n_genes, dtype=bool)+    ref_types = [str(t) for t in np.unique(ref_labels) if (ref_labels == t).sum() >= MIN_CELLS]+    off_types = [str(t) for t in np.unique(last_labels) if (last_labels == t).sum() >= MIN_CELLS]+    if not ref_types or not off_types:+        return {}, []+    pb_ref = _pseudobulk(X_ref, ref_labels, ref_types)+    pb_off = _pseudobulk(X_last, last_labels, off_types)++    pairs = []+    if label_match:+        pos = {t: j for j, t in enumerate(off_types)}+        for i, t in enumerate(ref_types):+            if t in pos:+                pairs.append((i, pos[t], 1.0))+    if not pairs:+        pairs = match_types(pb_ref, ref_types, pb_off, off_types, covered)+    if not pairs:+        return {}, []++    # raw per-pair deltas on covered genes+    deltas = np.zeros((len(pairs), n_genes), dtype=np.float64)+    ses = np.zeros((len(pairs), n_genes), dtype=np.float64)+    for p, (i, j, _c) in enumerate(pairs):+        deltas[p, covered] = (pb_ref[i] - pb_off[j])[covered]+        mr = ref_labels == ref_types[i]+        mo = last_labels == off_types[j]+        nr, no = int(mr.sum()), int(mo.sum())+        vr = np.asarray(X_ref[mr][:, covered].power(2).mean(axis=0)).ravel() - pb_ref[i][covered] ** 2+        vo = np.asarray(X_last[mo][:, covered].power(2).mean(axis=0)).ravel() - pb_off[j][covered] ** 2+        ses[p, covered] = np.sqrt(np.maximum(vr, 0) / max(nr, 1) + np.maximum(vo, 0) / max(no, 1))++    # location-shift batch correction: centre per gene across matched types+    deltas[:, covered] -= np.median(deltas[:, covered], axis=0, keepdims=True)++    # EB shrinkage on per-gene t statistic+    t2 = (deltas / np.maximum(ses, 1e-9)) ** 2+    deltas *= t2 / (t2 + K_SHRINK)++    # winsorise, kill cell-cycle genes and uncovered genes+    finite = deltas[np.isfinite(deltas)]+    if finite.size:+        cap = np.quantile(np.abs(finite[finite != 0]) if (finite != 0).any() else finite, CLIP_Q)+        deltas = np.clip(deltas, -cap, cap)+    deltas[:, ~covered] = 0.0+    deltas[:, cc_mask] = 0.0++    out = {}+    log = []+    for p, (i, j, c) in enumerate(pairs):+        out[off_types[j]] = deltas[p]+        log.append((ref_types[i], off_types[j], c))+    return out, log+++def apply_shift(X, row_types, delta: dict[str, np.ndarray], scale: float):+    """X (CSR, prediction-sized) += scale * delta[type] per row, clipped at 0."""+    if scale == 0.0 or not delta:+        return X+    D = np.asarray(X.todense(), dtype=np.float32)+    for t in sorted({str(s) for s in row_types if str(s) in delta}):+        rows = np.flatnonzero(row_types == t)+        if rows.size == 0:+            continue+        D[rows] += np.float32(scale) * delta[t].astype(np.float32)[None, :]+    np.clip(D, 0.0, None, out=D)+    out = sparse.csr_matrix(D)+    out.eliminate_zeros()+    return out

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

用到的知识库条目

计划里引用的来源

• DOI 10.1038/s41587-021-01001-7 — scArches：跨数据集参考映射中，基因级位置（location-shift）批校正可在去批的同时保留生物状态差异，支撑步骤 2 对 Qiu/官方跨技术 δ 做逐基因 median 中心化的做法
• PMID 30451828 — Jia 2018（Nkx2-5/Isl1 心脏祖细胞 E7.5–E9.5 单细胞图谱）：Isl1/Nkx2-5 标记心脏祖细胞、Nkx2-5 表达指向单向心肌命运，用作 Qiu 心脏谱系与官方类型标签对不上时的保守标记映射依据

文献检索

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

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

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

分析员

ANALYSIS.json written

工程师

完成。总结：

**做了什么**：在父节点 4（55.97）上移植节点 9 的类内成熟度加权抽样（周期分深度残差化稳健 z，w=exp(-0.5z)，E-S 无放回，组成不变）作为主干；新实现跨阶段伪批量位移分支（Qiu E9.0 外部输入，scArches 式逐基因 median 批校 + EB 收缩 + 相关匹配类型），α 环境变量门控。

**查分结果**（用了 8/20）：默认配置 proxy 3-seed 均值 **57.48**（57.76/57.61/57.08），proxy2 与 proxy 逐位相同；位移分支负结果——α=0.5/1.0 在 proxy2 单调降分（56.78/55.91 vs 57.76），四组全降，按 PLAN 预设规则放弃，提交 α=0；θ=1.0/0.25、heart×2.0 均不优于 θ=0.5。两视图 vec-check ok，运行 ~2.4s。

**下一步最值得试**：final 视图上启用官方两阶段（E8.5→E9.5）位移——同数据集无批混杂，代码路径已就绪（α≈0.25-0.5）；以及把 DE 方向信息加权进类内抽样而非平移表达（θ=1.0 时 de_score 反降，说明"更成熟"≠"更对的 DE"）。

调研员

al Tube/4000 细胞）。查分 seed 0 验证 ≈57.92、p1≈p2，并用 α=0 时 np.array_equal 锁定基线。\n步骤 1（发现第二输入）：在 proxy2 数据目录里找 Qiu E9.0——先查 manifest 的 inputs_by_time，再查 external/ 目录（任务书明示 external/ 属于视图数据）；读入其表达矩阵、细胞类型标签、阶段时间。取与官方面板的共有基因掩码（预期 ~27,883）；按标签精确匹配两阶段的共有类型，只保留两侧 n≥30 细胞的类型（Qiu 仅心脏谱系，预期只有 CM/心 progenitor 等少数类型能配对；若共有类型 <2 或标签完全对不上，用 Isl1/Nkx2-5/Tnnt2 等心脏标记的定性对应做保守映射，仍不行则放弃位移分支、只交步骤 0）。\n步骤 2（位移估计，全部在运行期由视图数据计算，无硬编码统计量）：对每个共有类型 t、共有基因 g：raw δ_g(t)=mean_Qiu,E9.0(t,g) − mean_official,E8.5(t,g)。跨技术批效应按基因级位置漂移去除（scArches 式 location-shift 假设）：δ_g(t) ← δ_g(t) − median_over_shared_types(δ_g(t))，即只保留类型间对比。EB 收缩：per-gene t 统计 t_g=δ_g/se_g，se_g=sqrt(v1/n1+v2/n2)，δ_shrunk=δ·t²/(t²+k)，k 初值 4（不搜，样本小）。细胞周期基因集内 δ 强制置 0；|δ| 在 99 分位截断；非共有类型与缺失的 ~4,402 个基因 δ=0。\n步骤 3（应用）：对抽样出的官方最后阶段细胞，属于共有类型者 x ← x + α·r·δ_shrunk(t)，负值 clip 到 0。r=Δt(最后输入→目标)/Δt(两输入间隔)，由阶段时间算出：proxy2 为 0.5/0.5=1，final 为 1/1=1（final 输入是官方 E8.5+E9.5，同数据集，δ 无批混杂，类型间 median 中心化仍照做、无害），代码零分支。门控：仅当视图里存在第二输入阶段时启用；proxy（单输入）恒等，用一次查分确认 p1 不变。α 用环境变量 VEC_ALPHA，初值 0.5，搜索 {0, 0.25, 0.5, 1.0}。\n步骤 4（查分预算，共 ≤14/20）：α=0 基线 1 次；α∈{0.25,0.5,1.0} 各 seed 0 一次看 p2/de_recovery 走向；选最优 2 个 α（含 0）各跑 3 seeds（6 次）；剩余预算做诊断复跑。出货规则：3-seed 均分最高者；若最优 α>0 与 α=0 差 <1（噪声内），保守出 α=0（=节点 9 克隆，下界 57.92）。",
  "expected_groups": ["de_recovery", "direction"],
  "risks": "1) Qiu 批效应不是纯基因级位置漂移（类型×技术交互），median 中心化去不干净 → α>0 时 p2 的 de_recovery/direction 反降；α=0.25 的首次查分即可暴露，de_recovery 相对 α=0 降 >1 就放弃位移分支只交步骤 0。2) external/ 里 Qiu 的格式/标签体系与预期不符，或共有类型 <2 → 步骤 1 用 10 行脚本先探查再写主代码，探查失败即触发只交步骤 0 的退路（下界仍 57.92，高于父本）。3) 位移破坏协方差或产生负表达 → clip 0、α 上限 1.0、细胞周期基因置 0；3-seed 检查 covariation 降幅 >2 则回退。4) 噪声：T1 约 2 分，0.5·Δp2 需 Δp2≥4 才在榜分上稳超噪声，单次差一律不作数，只认 3-seed 均值。5) 30 分钟时限：步骤 0 完成后即有一个可交版本，位移分支任何环节超时都直接交步骤 0。",
  "sources": ["DOI 10.1038/s41587-021-01001-7 — scArches：跨数据集参考映射中，基因级位置（location-shift）批校正可在去批的同时保留生物状态差异，支撑步骤 2 对 Qiu/官方跨技术 δ 做逐基因 median 中心化的做法", "PMID 30451828 — Jia 2018（Nkx2-5/Isl1 心脏祖细胞 E7.5–E9.5 单细胞图谱）：Isl1/Nkx2-5 标记心脏祖细胞、Nkx2-5 表达指向单向心肌命运，用作 Qiu 心脏谱系与官方类型标签对不上时的保守标记映射依据"]
}
```