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

节点 n18

改了什么

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

来自 programs.git 93025bf760:solution/METHOD.md

调研员的计划

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

对比：这个提交的上一版（种子程序：相对空仓库）。改动的文件：solution/EXECUTION.json +1 −0、solution/METHOD.md +19 −0、solution/README.md +4 −0、solution/g37_common.py +101 −0、solution/run.py +182 −0

diff --git a/solution/EXECUTION.json b/solution/EXECUTION.jsonnew file mode 100644index 0000000..9d5125c--- /dev/null+++ b/solution/EXECUTION.json@@ -0,0 +1 @@+{"gpu": false}diff --git a/solution/METHOD.md b/solution/METHOD.mdnew file mode 100644index 0000000..47dfd59--- /dev/null+++ b/solution/METHOD.md@@ -0,0 +1,19 @@+## 改了什么+将 ot_moscot 的 moscot 耦合 + 重心位移外推替换为 manifold_ode 方向的最小实现：+1. 在联合 PCA 空间训练一个小型 MLP 向量场 v(x,t)（2 层 64 单元，输入拼接时间维），用固定步 Euler 积分。+2. 训练损失 = geomloss Sinkhorn 散度（推源到目标）+ 速度 L2 正则 + kNN 流形惩罚（轨迹点离观测数据的中位距离）。+3. 从最后输入阶段向前积分，步长按 dt_out/dt_in × 阻尼系数 0.7 缩放（官方结果：无阻尼自治场不如 copy_last，故加阻尼）。+4. 解码：PCA 线性解码差值 → 只加在非零基因上（addnz 策略，与父节点一致），夹到 ≥0。+5. 输出细胞仍按生长率 g^dt_out 重抽样（保留 WOT 生长先验）。+6. 机制对照：环境变量 MANIFOLD_CONSTRAINT=0 关闭 kNN 流形惩罚，其余不变。+7. 单输入阶段（proxy）退化为生长加权重抽样，与父节点一致。+8. 第 1 轮修复：将 `import torch` 移至模块顶层（原来在 main() 内部，导致类定义时 NameError）。++关键超参：HIDDEN=64, EPOCHS=300, INTEG_STEPS=10, TRAIN_SUBSAMPLE=2000, VEL_REG=0.01, MANIFOLD_REG=0.5, MANIFOLD_K=15, DAMPING=0.7, LR=1e-3。++## 用到的知识与出处+- k032: TrajectoryNet / MIOFlow / scNODE 文献；实践建议"PCA 30-50d + Sinkhorn + 动能惩罚 + 收缩外推"；"官方自治场不如常数位移，必须非自治或阻尼"。+- k034: torchcfm / OT-CFM 方法；"在 PCA 潜空间工作，线性解码 + 逐细胞残差"；"只有两个时间点时场弱辨识，外推需阻尼"。+- k031: geomloss SamplesLoss 用法（已安装）。+- 方法卡 T1 卡：addnz 解码策略、生长率重抽样、平移幅度收缩系数 α∈[0,1] 的建议。+- Schiebinger et al. Cell 2019（WOT 生长先验）；Tong et al. arXiv:2002.04461（TrajectoryNet）；Huguet et al. arXiv:2206.14928（MIOFlow 流形约束）。diff --git a/solution/README.md b/solution/README.mdnew file mode 100644index 0000000..dc1cbaf--- /dev/null+++ b/solution/README.md@@ -0,0 +1,4 @@+# 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。diff --git a/solution/g37_common.py b/solution/g37_common.pynew file mode 100644index 0000000..8c991b3--- /dev/null+++ b/solution/g37_common.py@@ -0,0 +1,101 @@+"""Helpers of the ot_moscot seed (copy of modeling/candidates/T1/ot_moscot/g37_common.py, G37; stage_pair reads every+input stage of the view).++- ``stage_pair``: the two latest input stages (official, or external in proxy2 / test views) and the gene mask on+  which both are really measured (external stages lack part of the panel; filled columns must not drive dynamics);+- ``embed``: HVG -> z-score (clip 10) -> PCA, fitted on the input stages only (scanpy-style preprocessing);+- ``knn_mean``: per-cell kNN average in the embedding (smooths single-cell noise before taking displacements);+- ``growth_rates``: Waddington-OT / moscot prior growth from proliferation and apoptosis gene scores.+"""++from __future__ import annotations++import numpy as np+from scipy import sparse++from src.task1_temporal.view_io import covered_mask, inputs_by_time, read_stage+++def stage_pair(view, manifest: dict, genes: list[str]):+    """(prev, last, mask, dt_in) for the two latest inputs; prev None with a single input stage."""+    stages = inputs_by_time(manifest, include_external=True)  # every input stage, whatever its source+    last_e = stages[-1]+    last = read_stage(view, last_e, genes)+    mask = covered_mask(view, last_e, genes)+    if len(stages) < 2:+        return None, last, mask, None, last_e+    prev_e = stages[-2]+    prev = read_stage(view, prev_e, genes)+    mask &= covered_mask(view, prev_e, genes)+    return prev, last, mask, float(last_e["time"] - prev_e["time"]), last_e+++def _gene_moments(X: sparse.csr_matrix):+    n = X.shape[0]+    mean = np.asarray(X.mean(axis=0)).ravel()+    sq = np.asarray(X.multiply(X).sum(axis=0)).ravel() / n+    return mean, np.maximum(sq - mean ** 2, 0.0)+++def embed(mats: list, mask: np.ndarray, n_hvg: int, n_pcs: int, seed: int):+    """PCA coordinates of each matrix in ``mats`` (same order), fitted on all of them together.++    Returns (list of Z, info) with info = {hvg, mean, std, components} so displacements can be decoded."""+    from sklearn.decomposition import PCA++    allX = sparse.vstack(mats).tocsr()+    mean, var = _gene_moments(allX)+    var = np.where(mask, var, -1.0)+    hvg = np.sort(np.argsort(var)[::-1][:n_hvg])+    sub = allX[:, hvg].toarray().astype(np.float32)+    mu, sd = mean[hvg].astype(np.float32), np.sqrt(var[hvg]).astype(np.float32)+    sd[sd == 0] = 1.0+    sub -= mu+    sub /= sd+    np.clip(sub, -10, 10, out=sub)+    pca = PCA(n_components=n_pcs, svd_solver="randomized", random_state=seed)+    Z = pca.fit_transform(sub).astype(np.float32)+    out, start = [], 0+    for M in mats:+        out.append(Z[start:start + M.shape[0]])+        start += M.shape[0]+    return out, {"hvg": hvg, "mu": mu, "sd": sd, "components": pca.components_.astype(np.float32)}+++def knn_mean(Z: np.ndarray, X: sparse.csr_matrix, rows: np.ndarray, k: int) -> np.ndarray:+    """Dense (len(rows), n_genes): mean expression of the k nearest cells (in Z, the cell itself included)."""+    from sklearn.neighbors import NearestNeighbors++    nn = NearestNeighbors(n_neighbors=k).fit(Z)+    idx = nn.kneighbors(Z[rows], return_distance=False)+    W = sparse.csr_matrix((np.full(idx.size, 1.0 / k, dtype=np.float32), idx.ravel(),+                           np.arange(0, idx.size + 1, k)), shape=(len(rows), Z.shape[0]))+    return np.asarray((W @ X).todense(), dtype=np.float32)+++# Waddington-OT prior growth (Schiebinger et al. 2019, Cell; as implemented in moscot.base.problems.birth_death):+# birth = generalised logistic of the proliferation score, death = of the apoptosis score, g = exp(birth - death) / day.+def _gen_logistic(p, sup, inf, center, width):+    return inf + (sup - inf) / (1 + np.exp(-(p - center) / width))+++def growth_rates(adata_list: list, genes: list[str], mask: np.ndarray, seed: int) -> list[np.ndarray]:+    """Per-day growth rate of every cell of each AnnData (scores computed on all of them jointly)."""+    import anndata as ad+    import scanpy as sc+    from moscot.utils.data import apoptosis_markers, proliferation_markers++    ok = set(np.asarray(genes)[mask])+    joint = ad.concat(adata_list, join="inner")+    prolif = [g for g in proliferation_markers("mouse") if g in ok]+    apopt = [g for g in apoptosis_markers("mouse") if g in ok]+    sc.tl.score_genes(joint, prolif, score_name="proliferation", random_state=seed)+    sc.tl.score_genes(joint, apopt, score_name="apoptosis", random_state=seed)+    birth = _gen_logistic(joint.obs["proliferation"].to_numpy(float), 1.7, 0.3, 0.25, 0.5)+    death = _gen_logistic(joint.obs["apoptosis"].to_numpy(float), 1.7, 0.3, 0.1, 0.2)+    g = np.exp(birth - death)+    out, start = [], 0+    for a in adata_list:+        out.append(g[start:start + a.n_obs])+        start += a.n_obs+    return outdiff --git a/solution/run.py b/solution/run.pynew file mode 100644index 0000000..e1e0f99--- /dev/null+++ b/solution/run.py@@ -0,0 +1,182 @@+#!/usr/bin/env python3+"""manifold_ode: Neural ODE in PCA latent with Sinkhorn distribution matching, velocity regularization,+and kNN manifold constraint. Extrapolate from last input stage to target time.++Two input stages:+  1. Joint PCA (HVG, z-score, 30 PCs) fitted on inputs only.+  2. Growth-weighted resampling of output cells (WOT birth-death model).+  3. Train a small MLP vector field v(x,t) with:+     - Sinkhorn divergence between ODE-pushed source and target (geomloss)+     - Velocity L2 regularization+     - kNN manifold penalty: penalize trajectory points far from observed data+  4. Integrate last-stage cells forward by dt_out with damping.+  5. Decode: linear PCA decode + per-cell residual (addnz strategy).++One input stage (proxy): growth-weighted copy only.++Mechanism control: MANIFOLD_CONSTRAINT env var (default "1" = on, "0" = off).+"""++from __future__ import annotations++import argparse+import os+import sys+from pathlib import Path++os.environ["JAX_PLATFORMS"] = "cpu"++import numpy as np+import torch++sys.path.insert(0, str(Path(__file__).resolve().parent))+from g37_common import embed, growth_rates, knn_mean, 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+K_SMOOTH = 30+GROWTH = True+N_THREADS = 8+HIDDEN = 64+N_TRAIN_EPOCHS = 300+N_INTEG_STEPS = 10+TRAIN_SUBSAMPLE = 2000+VEL_REG = 0.01+MANIFOLD_REG = 0.5+MANIFOLD_K = 15+DAMPING = 0.7+LR = 1e-3+++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))+++class VectorField(torch.nn.Module):+    def __init__(self, dim, hidden):+        super().__init__()+        self.net = torch.nn.Sequential(+            torch.nn.Linear(dim + 1, hidden),+            torch.nn.Softplus(),+            torch.nn.Linear(hidden, hidden),+            torch.nn.Softplus(),+            torch.nn.Linear(hidden, dim),+        )+        torch.nn.init.zeros_(self.net[-1].weight)+        torch.nn.init.zeros_(self.net[-1].bias)++    def forward(self, x, t):+        t_expanded = t.expand(x.shape[0], 1)+        return self.net(torch.cat([x, t_expanded], dim=1))+++def euler_integrate(vf, x0, t0, t1, n_steps):+    dt = (t1 - t0) / n_steps+    x = x0.clone()+    for i in range(n_steps):+        t = torch.full((x.shape[0], 1), t0 + i * dt, dtype=x.dtype)+        x = x + dt * vf(x, t)+    return x+++def manifold_penalty(x, data, k):+    from sklearn.neighbors import NearestNeighbors+    data_np = data.detach().numpy()+    x_np = x.detach().numpy()+    nn = NearestNeighbors(n_neighbors=k).fit(data_np)+    dists, _ = nn.kneighbors(x_np)+    return torch.tensor(dists[:, k // 2].mean(), dtype=torch.float32)+++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)+    args = parser.parse_args()++    torch.set_num_threads(N_THREADS)+    torch.manual_seed(args.seed)++    use_manifold = os.environ.get("MANIFOLD_CONSTRAINT", "1") != "0"++    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)++    if prev is None:+        (g_last,) = growth_rates([last], genes, mask, args.seed)+        rows = weighted_rows(g_last ** dt_out if GROWTH else np.ones(last.n_obs), n, rng)+        write_prediction(last.X[rows], genes, args.out, seed=args.seed)+        return++    (Zp, Zl), pca_info = embed([prev.X, last.X], mask, N_HVG, N_PCS, args.seed)+    g_prev, g_last = growth_rates([prev, last], genes, mask, args.seed)++    rows = weighted_rows(g_last ** dt_out if GROWTH else np.ones(last.n_obs), n, rng)++    Zp_t = torch.from_numpy(Zp.astype(np.float32))+    Zl_t = torch.from_numpy(Zl.astype(np.float32))++    n_prev_sub = min(TRAIN_SUBSAMPLE, Zp_t.shape[0])+    n_last_sub = min(TRAIN_SUBSAMPLE, Zl_t.shape[0])+    idx_p = torch.from_numpy(rng.choice(Zp_t.shape[0], n_prev_sub, replace=False).astype(np.int64))+    idx_l = torch.from_numpy(rng.choice(Zl_t.shape[0], n_last_sub, replace=False).astype(np.int64))+    Zp_sub = Zp_t[idx_p]+    Zl_sub = Zl_t[idx_l]++    all_data = torch.cat([Zp_t, Zl_t], dim=0)++    vf = VectorField(N_PCS, HIDDEN)+    optimizer = torch.optim.Adam(vf.parameters(), lr=LR)++    from geomloss import SamplesLoss+    sinkhorn = SamplesLoss(loss="sinkhorn", p=2, blur=0.05, truncate=5)++    for epoch in range(N_TRAIN_EPOCHS):+        optimizer.zero_grad()+        pushed = euler_integrate(vf, Zp_sub, 0.0, 1.0, N_INTEG_STEPS)+        loss_match = sinkhorn(pushed, Zl_sub)+        t_mid = torch.full((Zp_sub.shape[0], 1), 0.5)+        vel = vf(Zp_sub, t_mid)+        loss_vel = VEL_REG * vel.pow(2).mean()+        loss = loss_match + loss_vel+        if use_manifold:+            loss_man = MANIFOLD_REG * manifold_penalty(pushed, all_data, MANIFOLD_K)+            loss = loss + loss_man+        loss.backward()+        optimizer.step()++    vf.eval()+    with torch.no_grad():+        Zl_out = Zl_t[torch.from_numpy(rows.astype(np.int64))]+        ratio = dt_out / dt_in if dt_in > 0 else 1.0+        Zl_pushed = euler_integrate(vf, Zl_out, 1.0, 1.0 + ratio * DAMPING, N_INTEG_STEPS)+        Zl_pushed_np = Zl_pushed.numpy()++    components = pca_info["components"]+    hvg = pca_info["hvg"]+    mu = pca_info["mu"]+    sd = pca_info["sd"]++    decoded_step = (Zl_pushed_np - Zl[rows]) @ components+    decoded_step_full = np.zeros((n, len(genes)), dtype=np.float32)+    decoded_step_full[:, hvg] = decoded_step * sd++    X = last.X[rows].toarray()+    step = decoded_step_full.copy()+    step[:, ~mask] = 0.0+    step *= X > 0+    X += step+    np.maximum(X, 0.0, out=X)+    write_prediction(X, genes, args.out, seed=args.seed)+++if __name__ == "__main__":+    main()

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

没有记录调研来源。

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

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

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

分析员

ANALYSIS.json written

native

stream

<<<<<<< SEARCH
MANIFOLD_REG = 0.5
MANIFOLD_K = 15
DAMPING = 0.7
LR = 1e-3
=======
MANIFOLD_REG = 0.5
MANIFOLD_K = 15
DAMPING = 0.5
LR = 2e-3
BLEND_KNN = 0.3
>>>>>>> REPLACE

<<<<<<< SEARCH
N_TRAIN_EPOCHS = 300
=======
N_TRAIN_EPOCHS = 500
>>>>>>> REPLACE

<<<<<<< SEARCH
def manifold_penalty(x, data, k):
    from sklearn.neighbors import NearestNeighbors
    data_np = data.detach().numpy()
    x_np = x.detach().numpy()
    nn = NearestNeighbors(n_neighbors=k).fit(data_np)
    dists, _ = nn.kneighbors(x_np)
    return torch.tensor(dists[:, k // 2].mean(), dtype=torch.float32)
=======
def manifold_penalty(x, data, k):
    dists = torch.cdist(x.unsqueeze(0), data.unsqueeze(0)).squeeze(0)
    topk = dists.topk(k, dim=1, largest=False).values
    return topk[:, k // 2].mean()
>>>>>>> REPLACE

<<<<<<< SEARCH
    vf = VectorField(N_PCS, HIDDEN)
    optimizer = torch.optim.Adam(vf.parameters(), lr=LR)

    from geomloss import SamplesLoss
    sinkhorn = SamplesLoss(loss="sinkhorn", p=2, blur=0.05, truncate=5)