总览 · ← 返回运行 20261002-202908-search-t1-scr-B

节点 n5

低秩投影逐簇Delta外推：双阶段KMeans簇配对，簇间伪批量差投影到目标簇自身前15个PC并乘α=2.0，仅加到非零表达项，单输入或α=0时逐比特退化为copy_last。

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

来自 programs.git a21b2adf83:solution/METHOD.md

调研员的计划

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

对比：父节点版本 75c2ae00e3。改动的文件：solution/EXECUTION.json +1 −0、solution/METHOD.md +43 −0、solution/run.py +122 −3

diff --git a/solution/EXECUTION.json b/solution/EXECUTION.jsonnew file mode 100644index 0000000..6d8012e--- /dev/null+++ b/solution/EXECUTION.json@@ -0,0 +1 @@+{"gpu": false}\ No newline at end of filediff --git a/solution/METHOD.md b/solution/METHOD.mdnew file mode 100644index 0000000..2ee9b63--- /dev/null+++ b/solution/METHOD.md@@ -0,0 +1,43 @@+低秩投影逐簇Delta外推：双阶段KMeans簇配对，簇间伪批量差投影到目标簇自身前15个PC并乘α=2.0，仅加到非零表达项，单输入或α=0时逐比特退化为copy_last。++## 方法（family: lowrank_shape，按 PLAN 实现）++1. 按时间读入最后两个输入阶段 t0、t1（`view_io.read_stage`，全基因面板）。只有一个输入阶段或 `--alpha 0` 时走与父节点完全相同的 rng/抽样路径，输出逐比特等于 copy_last（已本地验证 `np.array_equal(data/indices/indptr)` 为 True）。+2. 在两阶段合并数据上按 var/mean 离散度选 2000 个 HVG；合并中心化后 SVD 取前 25 个 PC。+3. 在 PC 空间对 t0、t1 各做 KMeans(k=8, n_init=10, random_state=0)；每个 t1 簇按最近质心匹配一个 t0 簇（X3 上 8 对匹配互不相同，配对距离 0.5–1.2，远小于簇尺度）。+4. 每对匹配簇：delta_j = mean(t1[c1]) − mean(t0[c0])，定义在全 32,285 基因上（log 空间）。+5. 核心机制：对 t1 簇 c1 的中心化表达矩阵（全基因）做经济 SVD/Gram 特征分解，取前 r=15 个右奇异向量 V；delta_lr = V(Vᵀ·delta)，即把位移限制在该簇自身共变方向张成的低秩子空间内。位移 = α·delta_lr，α=2.0。+6. 把位移加到 t1 抽样细胞（所属簇）上，`--support nz`：只加在原始非零表达项上（保留稀疏支撑）；负值截断到 0，<1e-3 归零。抽样细胞数、rows 与父节点同 rng 同路径。++## 机制生效证据（X3 视图，seed 0，stderr [diag] 输出）++- 8 个簇全部被位移，覆盖全部 652 个输出细胞（每簇 44–122 个）。+- 每簇投影保留率 ||delta_lr||/||delta|| = 0.22–0.51（低秩投影确实改变了位移方向，不是恒等）。+- ||delta|| ≈ 5.7–9.2，α=2 时 ||shift|| ≈ 1.1–3.0（log 单位，全基因 RMS 每基因 ~0.01–0.02）。++## 查分结果（vec-score --task X3_qiu_heart_early，A 半）++| 配置 | board | de_recovery | covariation | cell_state | direction |+|---|---|---|---|---|---|+| 父节点 copy_last | 47.92 | 44.09 | 48.44 | 49.68 | 49.20 |+| α=0.5, r=15, support=all | 44.08 | 43.09 | **32.43** | 48.62 | 48.95 |+| α=0.5, r=15, support=nz | 47.42 | 42.74 | 48.32 | 49.37 | 49.04 |+| α=1.0, r=15, support=nz | 47.62 | 43.80 | 48.24 | 49.26 | 48.97 |+| **α=2.0, r=15, support=nz（提交默认）** | **48.12** | **46.09** | 48.17 | 49.11 | 48.95 |++结论与教训：+- PLAN 的 support=all（全基因稠密位移）会把 covariation 从 48.4 打到 32.4——对零表达项加位移破坏稀疏支撑/基因共现结构，是本视图上最强的破坏因素；改为只作用于非零项（support=nz）后 covariation 恢复到 ≈48.2（与父节点差 <0.3，噪声内）。这验证了机制卡里"低秩投影保协变"的预期方向，但真正决定 covariation 的是支撑集而非投影。+- de_recovery 随 α 单调上升（42.7→43.8→46.1），α=2.0 时超过父节点 +2.0；covariation/cell_state 仅微降（−0.27/−0.57，噪声内）。总分 48.12 vs 47.92，提升 +0.2，**在 ±2 噪声内，不能声称确定进步**。+- 未测：α>2、rank∈{10,20,all} 的完整网格、k=12（时间预算耗尽）；de_recovery 对 α 的趋势提示 α=3–4 可能更高，留给后续节点。++## 对照（mechanism_off_control）++- `--alpha 0`：输出逐比特等于 copy_last（已验证），机制未运行时无任何差异。+- `--rank all`：跳过投影的 raw-delta 消融开关已实现（本次未查分）。++## 其他++- 未用任何保留阶段/保留基因型信息；未用 external/、prior/；无硬编码阶段名、细胞数或统计量；所有量现场从 view 计算。+- 视图无关：分支只依赖输入阶段个数与 `--seed`/超参，不读绝对时间、路径、board 字段；时间平移不影响输出。+- 确定性：rng 仅 `np.random.default_rng(seed)`；KMeans/SVD 固定 random_state 与稳定排序。+- 纯 CPU；运行时 X3 上约 90 s、峰值内存 <3 GB（final 视图更大但按 28 GB 限内估算：两阶段稠密 float32 + 每簇 Gram 分解）。diff --git a/solution/run.py b/solution/run.pyindex aa4bbf0..5b5ddfc 100644--- a/solution/run.py+++ b/solution/run.py@@ -1,11 +1,23 @@ #!/usr/bin/env python3-"""copy_last: the latest input stage, subsampled into the board's cell range."""+"""lowrank_cluster_delta.++copy_last baseline plus a per-cluster temporal displacement: for each cluster+of the last input stage, the empirical pseudobulk delta (last stage minus the+matched earlier-stage cluster) is projected onto the cluster's own top-r+principal directions of within-cluster variation (all genes), shrunk by alpha,+and added to every sampled cell of that cluster; negatives clipped to 0.++Controls: --alpha 0 (or a single input stage) degrades bit-for-bit to+copy_last; --rank all skips the projection (raw-delta ablation).+"""  from __future__ import annotations  import argparse+import sys  import numpy as np+import scipy.sparse as sp  from src.task1_temporal.view_io import (     inputs_by_time,@@ -17,20 +29,127 @@ from src.task1_temporal.view_io import (     write_prediction, ) +N_HVG = 2000+N_PC = 25+++def _dense(X) -> np.ndarray:+    if sp.issparse(X):+        return np.asarray(X.todense(), dtype=np.float32)+    return np.asarray(X, dtype=np.float32)+++def _hvg(D0: np.ndarray, D1: np.ndarray, n_hvg: int = N_HVG) -> np.ndarray:+    A = np.concatenate([D0, D1], axis=0)+    mean = A.mean(axis=0)+    var = A.var(axis=0)+    disp = np.where(mean > 1e-6, var / np.maximum(mean, 1e-6), 0.0)+    order = np.argsort(-disp, kind="stable")+    return np.sort(order[: min(n_hvg, D0.shape[1])])+++def _cluster_pcs(Xc: np.ndarray, r: int, seed: int) -> np.ndarray:+    """Top-r right singular vectors of the centered cluster matrix (n x G)."""+    n = Xc.shape[0]+    r = int(max(1, min(r, n - 1)))+    if n < 3:+        return np.zeros((0, Xc.shape[1]), dtype=np.float64)+    Xd = Xc.astype(np.float64)+    Xd -= Xd.mean(axis=0)+    if Xd.shape[0] <= Xd.shape[1]:+        K = Xd @ Xd.T+        w, V = np.linalg.eigh(K)+        take = np.argsort(-w, kind="stable")[:r]+        U = V[:, take]+        return (U.T @ Xd) / np.sqrt(np.maximum(w[take], 1e-12))[:, None]+    _, _, Vt = np.linalg.svd(Xd, full_matrices=False)+    return Vt[:r]+  def main() -> None:     parser = argparse.ArgumentParser()     parser.add_argument("--data", required=True)     parser.add_argument("--out", required=True)     parser.add_argument("--seed", type=int, default=0)+    parser.add_argument("--alpha", type=float, default=2.0)+    parser.add_argument("--rank", default="15", help="int or 'all'")+    parser.add_argument("--k", type=int, default=8)+    parser.add_argument("--support", default="nz", choices=["all", "nz"],+                        help="apply shift to all genes or only to nonzero entries")     args = parser.parse_args()      manifest = load_manifest(args.data)     genes = panel_genes(args.data, manifest)-    last = read_stage(args.data, inputs_by_time(manifest)[-1], genes)+    entries = inputs_by_time(manifest)+    last = read_stage(args.data, entries[-1], genes)     rng = np.random.default_rng(args.seed)     rows = sample_rows(last.n_obs, target_n_cells(manifest, last.n_obs), rng)-    write_prediction(last.X[rows], genes, args.out, seed=args.seed)++    if len(entries) < 2 or args.alpha == 0.0:+        write_prediction(last.X[rows], genes, args.out, seed=args.seed)+        return++    from sklearn.cluster import KMeans++    prev = read_stage(args.data, entries[-2], genes)+    D0 = _dense(prev.X)+    D1 = _dense(last.X)+    n0, n1 = D0.shape[0], D1.shape[0]++    hvg = _hvg(D0, D1)+    A = np.concatenate([D0[:, hvg], D1[:, hvg]], axis=0).astype(np.float64)+    A -= A.mean(axis=0)+    _, _, Vt = np.linalg.svd(A, full_matrices=False)+    npc = min(N_PC, Vt.shape[0] - 1)+    P = Vt[:npc]+    Z0 = D0[:, hvg].astype(np.float64) @ P.T+    Z1 = D1[:, hvg].astype(np.float64) @ P.T++    k = int(min(args.k, n0, n1))+    km0 = KMeans(k, n_init=10, random_state=0).fit(Z0)+    km1 = KMeans(k, n_init=10, random_state=0).fit(Z1)+    lab0, lab1 = km0.labels_, km1.labels_++    # match each t1 centroid to its nearest t0 centroid in PCA space+    m0 = km0.cluster_centers_+    m1 = km1.cluster_centers_+    dist = ((m1[:, None, :] - m0[None, :, :]) ** 2).sum(axis=2)+    match = dist.argmin(axis=1)++    mu0 = np.stack([D0[lab0 == j].mean(axis=0) for j in range(k)])++    rank = None if str(args.rank).lower() == "all" else int(args.rank)+    Xs = D1[rows].copy().astype(np.float64)+    lab_rows = lab1[rows]+    print(f"[diag] k={k} n0={n0} n1={n1} rows={len(rows)} alpha={args.alpha} rank={args.rank}", file=sys.stderr)+    for j in range(k):+        nj1 = int((lab1 == j).sum())+        c0 = match[j]+        delta = D1[lab1 == j].mean(axis=0).astype(np.float64) - mu0[c0].astype(np.float64)+        if rank is None:+            dlr = delta+        else:+            V = _cluster_pcs(D1[lab1 == j], rank, args.seed)+            dlr = (delta @ V.T) @ V if V.shape[0] else delta+        keep = float(np.linalg.norm(dlr) / max(np.linalg.norm(delta), 1e-12))+        sel = lab_rows == j+        nsel = int(sel.sum())+        sub = Xs[sel]+        if args.support == "nz":+            sub += args.alpha * dlr * (sub > 0)+        else:+            sub += args.alpha * dlr+        Xs[sel] = sub+        print(+            f"[diag] cluster {j}: n1={nj1} matched->t0[{c0}] dist={dist[j, c0]:.1f} "+            f"sel={nsel} ||delta||={np.linalg.norm(delta):.2f} retention={keep:.3f} "+            f"||shift||={args.alpha * np.linalg.norm(dlr):.2f}",+            file=sys.stderr,+        )++    np.maximum(Xs, 0.0, out=Xs)+    Xs[Xs < 1e-3] = 0.0+    write_prediction(Xs.astype(np.float32), genes, args.out, seed=args.seed)   if __name__ == "__main__":

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

用到的知识库条目

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

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

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

分析员

ANALYSIS.json written

工程师

α=2.0 gives 48.12 (de_recovery 46.09, cov 48.17) — best so far. Setting as default and finalizing.

调研员

pproach": "在 run.py 基础上（沿用 view_io：inputs_by_time/read_stage/panel_genes/sample_rows/write_prediction）：\n1) 读入两个输入阶段（按时间排序，t0、t1）；只有一个输入阶段时逐比特退化为 copy_last（与父节点同rng同抽样路径）。\n2) 选约2000个HVG（在合并两阶段上按方差/均值离散度选，通用实现，不依赖阶段名）。\n3) 对合并数据或t1做PCA，取约25维；在PCA空间对t0、t1分别KMeans聚类，k=8（搜索8,12）；每个t0簇按最近质心配到一个t1簇。\n4) 对每个匹配簇对(c0,c1)：delta_c = mean(t1[c1]) - mean(t0[c0])（在HVG上，log空间）。\n5) 核心机制：取t1簇c1内细胞的前r个主成分载荷 P（r初值15，搜索{10,20}），把delta投影到低秩子空间 delta_lr = P(P^T delta)；位移 = alpha*delta_lr，alpha初值0.5，搜索{0.3,0.5,0.8}。\n6) 把该位移向量整体加到t0该簇的每个抽样细胞的全部基因上，负值截断到0；非HVG列同样加上位移的对应分量（位移向量定义在全基因上，HVG外由PCA载荷自然决定，避免节点3'只改HVG'造成的基因子集不一致）。\n7) 抽样到 target_n_cells 后 write_prediction。\n验证流程：先2000细胞子采样跑通并用 vec-score 查分（X3），查 alpha×r 小网格（≤6次），挑 de_recovery 升且 covariation 不掉于48的组合，再全量跑并复查1次确认（差距<2分视为噪声）。总查分控制在8次内。",
 "expected_groups": ["de_recovery", "covariation"],
 "risks": "1) 若簇间delta主要落在低方差方向（与簇内PC正交），投影会把信号削没，de_recovery无改善——Engineer应先打印每个簇 ||delta_lr||/||delta|| 比值，若普遍<0.3则增大r或改用合并数据的PC；2) alpha过大仍可能轻度损伤covariation——以'前后簇内基因相关矩阵Frobenius差'做本地早筛，先于查分；3) 簇配对错误（k过大时）会引入假delta——用最近质心距离与配对唯一性检查；4) 收益可能只有2-3分，接近噪声——必须用小网格+复查确认，不要单次查分下结论；5) 30分钟时限：纯numpy+sklearn（KMeans/PCA），不用moscot/ott，避免JAX启动开销，先小样本后全量。",
 "family_id": "lowrank_shape",
 "mechanism": "把每个源-目标簇对的经验位移向量投影到目标簇自身表达变异的主分子空间（前r个PC），再做收缩后施加到该簇所有细胞：位移只能沿该群体真实存在的共变方向移动，从而在引入时间变化的同时保持基因间协方差结构。",
 "vs_constant_shift": "每型常数位移是把观测delta原样（或全局收缩后）平移整个簇，噪声方向也照搬；本方案的位移被限制在簇内协方差张成的低秩子空间，且收缩强度alpha与秩r构成对'信号子空间'的双重约束，另外簇配对是数据驱动的（跨时间点质心匹配）而非按类型名。与节点3的区别：不逐基因做独立EB收缩、不只改HVG非零项，而是整向量、全基因、沿协方差方向的位移，这正是修复协变损伤的结构性改动。",
 "mechanism_evidence": "Engineer应报告：1) 每簇 ||delta_lr||/||delta||（投影保留率）与位移范数；2) 实际改变了哪些细胞（各簇细胞数与位移分布）；3) 前后簇内基因相关矩阵的差异（应远小于节点3的对应值）；4) 四组分各自变化：预期de_recovery相对父节点+3~5且covariation不低于48；5) 与'不投影的raw-delta消融'对比，证明是投影而非单纯加位移带来协变保护。",
 "mechanism_off_control": "同一程序提供两个关闭方式：--alpha 0（位移为零）或仅一个输入阶段时，输出逐比特等于父节点 copy_last；另设 --rank all（不投影，等价于节点3式raw-delta对照）。预期：alpha=0 与父节点完全一致（机制未运行则不应有任何差异）；打开机制后 de_recovery 上升且簇内相关结构变化显著小于 rank=all 消融。",
 "sources": []}
```