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

节点 n15

Gaussian VAE 潜空间转移：两输入阶段潜码 entropic-OT 耦合+kNN 局部回归得逐细胞位移，外推末阶段潜码解码，HVG 分位数匹配回末阶段边缘分布后以 γ=0.5 与真实值混合；非 HVG 保留末阶段值。

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

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

调研员的计划

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

对比：这个提交的上一版（种子程序：相对空仓库）。改动的文件：solution/EXECUTION.json +1 −0、solution/METHOD.md +30 −0、solution/run.py +299 −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..7a5a1c5--- /dev/null+++ b/solution/METHOD.md@@ -0,0 +1,30 @@+Gaussian VAE 潜空间转移：两输入阶段潜码 entropic-OT 耦合+kNN 局部回归得逐细胞位移，外推末阶段潜码解码，HVG 分位数匹配回末阶段边缘分布后以 γ=0.5 与真实值混合；非 HVG 保留末阶段值。+# METHOD+++## 方法族与机制（family: generative_latent，按 PLAN）++- HVG 2000（各输入阶段都覆盖的基因，按合并方差选）；标准化后进 VAE：enc 2000(+2 阶段 one-hot)→256→64→z(20)，dec 20→64→256→2000，Gaussian 似然，KL 权重 **0.005**（PLAN 的 0.1 使潜空间塌缩，||dz||≈0.001，机制失效；实测降到 0.005 后 ||dz||≈0.41），Adam 1e-3，batch 256，100 epochs，CPU（EXECUTION.json gpu=false）。+- 转移：Z(E8.75 类) 与 Z(E9.0 类) 之间 POT `ot.sinkhorn`（reg=0.1，代价按中位数归一，每阶段 ≤3000 细胞）；δ_i = 重心映射 − z0_i；对每个末阶段细胞取 k=15 个最近源细胞，softmax(−d²/h)（h=中位 k 距）加权平均 δ 得 dz_j。z_pred = z_last + scale·dz，scale 默认常数 1.0（auto 模式=时间差之比 clip[0,2]，X3 上为 2.0，实测略差）。+- 解码后处理（PLAN 步骤 7 的修改，实测必要）：直接把解码值写回会摧毁表达纹理（cell_state 8.9、covariation 8.4，总 28.0，远低于 copy_last）。改为 (a) 每个 HVG 基因把解码值的秩映射到末阶段该基因的经验分位数（恢复稀疏与边缘分布），(b) 与末阶段真实值按 γ=0.5 线性混合（纯解码 γ=1 时 covariation 48 但 cell_state 35；混合把 cell_state 拉到 68，代价是 covariation 降到 28.6）。clamp ≥0。+- 单输入退路：无转移，输出=该阶段 encode-decode（自编码对照），符合 PLAN 步骤 8。+- 视图无关：只用时间差与数据；seed 确定（np.random.default_rng(seed)、torch.manual_seed(seed)，纯 CPU）。++## 机制生效证据++- `[mech]` 日志：||dz|| mean=0.4097，**within-population CV=0.337 > 0.3**（PLAN 证据标准 1）：位移确实随细胞潜位置变化，不是每型常数。+- 机制开/关对照（同一程序 `--transition-scale 0` vs 默认，X3 A 半，seed 0）：+  - 关（scale=0，纯解码无 QM 时）：28.0；关+QM+γ1 未单测；+  - 开：scale=1, γ=1, QM：42.98；scale=2, γ=1, QM：42.02；**scale=1, γ=0.5, QM（默认）：48.92**。+  - scale=0 与 scale>0 输出不同（dz≠0），机制在跑；转移主要改变 covariation/de_recovery 方向（de_recovery 44.9→41.7，direction 49.8→49.1），而 QM+γ 混合主要修复 cell_state。+- 输出零值率/库大小：nnz/cell 2273 vs 输入 1974（混合使部分零被填），在格式检查内合规。++## 验证与未验证++- 已验证：X3 视图（E8.75+E9.0→E9.5）vec-check 通过，5 次 vec-score（28.0 / 42.0 / 43.0 / 48.9 / 48.9-默认复现同哈希）。默认输出与打分过的 g05 配置逐位相同。+- 未验证：final（E8.5+E9.5→E10.5）与 proxy 视图未跑分（本节点只挂 X3 尺子）；代码路径对单输入（proxy）有退路、对大细胞数有 cap（train 10k/阶段、OT 3k/阶段、kNN 分块），但耗时只在 X3（~3.5k 细胞，1–1.5 min）实测。γ、kl、scale 只做了粗网格，B 半与多 seed 未测。+- 生物学知识来源：未使用禁窗（9.5<E≤13.5）或保留阶段的任何测量信息；方法只依赖视图内输入阶段数据与通用 ML 组件，无外部先验写入。++## 下一步++γ=0.5 的线性混合牺牲了 covariation（48→28.6）；更优方向是用潜位移做"真实细胞选择/组成调制"（对 z_pred 取最近真实末阶段细胞）而非数值混合，可同时保住 cell_state 与 covariation；或把解码 dz 投回基因空间做保稀疏位移（与节点 4/5/9 的强部件组合）。diff --git a/solution/run.py b/solution/run.pynew file mode 100644index 0000000..33b0567--- /dev/null+++ b/solution/run.py@@ -0,0 +1,299 @@+"""Gaussian VAE latent transition (family: generative_latent).++Train a small Gaussian VAE on HVGs of the input stage(s), couple the last two+input stages with entropic OT in latent space, build a per-cell displacement+field by kNN-weighted local regression over the source stage's displacements,+extrapolate the last stage's latents toward the target and decode. Non-HVG+genes keep the last stage's values. Single-input fallback: autoencode the only+input stage (no transition). Everything is computed from the view; only time+DIFFERENCES are used, so the program is view/time-shift invariant.+"""+from __future__ import annotations++import argparse+import json+from pathlib import Path++import numpy as np+from scipy import sparse+++# ---------------------------------------------------------------- data view++def load_manifest(view: Path) -> dict:+    return json.loads((view / "manifest.json").read_text(encoding="utf-8"))+++def panel_genes(view: Path, manifest: dict) -> list[str]:+    text = (view / manifest["genes_file"]).read_text(encoding="utf-8")+    return [ln.strip() for ln in text.splitlines() if ln.strip()]+++def covered_mask(view: Path, entry: dict, genes: list[str]) -> np.ndarray:+    if entry.get("genes_file"):+        text = (view / entry["genes_file"]).read_text(encoding="utf-8")+        have = {ln.strip() for ln in text.splitlines() if ln.strip()}+        return np.fromiter((g in have for g in genes), dtype=bool, count=len(genes))+    return np.ones(len(genes), dtype=bool)+++def read_input(view: Path, entry: dict, genes: list[str]) -> sparse.csr_matrix:+    """CSR float32 on the panel in panel order (self-contained; no src dep)."""+    import anndata as ad++    a = ad.read_h5ad(view / entry["path"])+    names = [str(g) for g in a.var_names]+    X = a.X+    X = X.tocsr() if sparse.issparse(X) else sparse.csr_matrix(np.asarray(X))+    X = X.astype(np.float32)+    if names == genes:+        return X.tocsr()+    pos = {g: i for i, g in enumerate(names)}+    cols = np.array([pos.get(g, -1) for g in genes], dtype=np.int64)+    keep = cols >= 0+    sub = X.tocsc()[:, cols[keep]].tocsc()+    keep_idx = np.flatnonzero(keep)+    per_col = np.zeros(len(genes), dtype=np.int64)+    per_col[keep_idx] = np.diff(sub.indptr)+    indptr = np.concatenate([[0], np.cumsum(per_col)])+    full = sparse.csc_matrix((sub.data, sub.indices, indptr), shape=(X.shape[0], len(genes)))+    return full.tocsr()+++# ---------------------------------------------------------------------- VAE++def make_vae(d_total: int, d_in: int, latent: int):+    import torch.nn as nn++    class VAE(nn.Module):+        def __init__(self):+            super().__init__()+            self.enc = nn.Sequential(nn.Linear(d_total, 256), nn.ReLU(), nn.Linear(256, 64), nn.ReLU())+            self.mu = nn.Linear(64, latent)+            self.logvar = nn.Linear(64, latent)+            self.dec = nn.Sequential(nn.Linear(latent, 64), nn.ReLU(), nn.Linear(64, 256), nn.ReLU(),+                                     nn.Linear(256, d_in))++        def forward(self, x):+            h = self.enc(x)+            return self.mu(h), self.logvar(h).clamp(-8.0, 2.0)++    return VAE()+++def train_vae(Xh: np.ndarray, cov: np.ndarray | None, latent: int, kl: float,+              epochs: int, seed: int):+    import torch++    torch.manual_seed(seed)+    torch.set_num_threads(4)+    d_in = Xh.shape[1]+    d_total = d_in + (0 if cov is None else cov.shape[1])+    model = make_vae(d_total, d_in, latent)+    opt = torch.optim.Adam(model.parameters(), lr=1e-3)+    X = torch.from_numpy(np.ascontiguousarray(Xh, dtype=np.float32))+    if cov is not None:+        X = torch.cat([X, torch.from_numpy(np.ascontiguousarray(cov, dtype=np.float32))], dim=1)+    n = X.shape[0]+    g = torch.Generator().manual_seed(seed)+    bs = 256+    model.train()+    for _ in range(epochs):+        perm = torch.randperm(n, generator=g)+        for i in range(0, n, bs):+            xb = X[perm[i:i + bs]]+            mu, logvar = model(xb)+            z = mu + torch.randn_like(mu) * torch.exp(0.5 * logvar)+            rec = model.dec(z)+            mse = ((rec - xb[:, :d_in]) ** 2).mean()+            kll = (-0.5 * (1 + logvar - mu.pow(2) - logvar.exp())).sum(-1).mean()+            loss = mse + kl * kll+            opt.zero_grad()+            loss.backward()+            opt.step()+    model.eval()+    return model, d_in+++def encode_mean(model, d_in: int, Xh: np.ndarray, cov: np.ndarray | None) -> np.ndarray:+    import torch++    with torch.no_grad():+        X = torch.from_numpy(np.ascontiguousarray(Xh, dtype=np.float32))+        if cov is not None:+            X = torch.cat([X, torch.from_numpy(np.ascontiguousarray(cov, dtype=np.float32))], dim=1)+        mu, _ = model(X)+        return mu.numpy().astype(np.float64)+++def decode(model, Z: np.ndarray) -> np.ndarray:+    import torch++    with torch.no_grad():+        return model.dec(torch.from_numpy(np.ascontiguousarray(Z, dtype=np.float32))).numpy()+++# --------------------------------------------------------------------- main++def main(argv=None):+    ap = argparse.ArgumentParser()+    ap.add_argument("--data", required=True)+    ap.add_argument("--out", required=True)+    ap.add_argument("--seed", type=int, required=True)+    ap.add_argument("--hvg", type=int, default=2000)+    ap.add_argument("--latent", type=int, default=20)+    ap.add_argument("--kl", type=float, default=0.005)+    ap.add_argument("--epochs", type=int, default=100)+    ap.add_argument("--quantile-match", type=int, default=1)+    ap.add_argument("--k", type=int, default=15)+    ap.add_argument("--ot-reg", type=float, default=0.1)+    ap.add_argument("--ot-cap", type=int, default=3000)+    ap.add_argument("--train-cap", type=int, default=10000)+    ap.add_argument("--transition-scale", type=float, default=1.0,+                    help="-1: auto = clip(dt_target/dt_input, 0, 2); 0 disables the mechanism (control)")+    ap.add_argument("--no-covariate", action="store_true")+    ap.add_argument("--blend-gamma", type=float, default=0.5,+                    help="HVG output = (1-g)*last-stage real values + g*quantile-matched decode")+    args = ap.parse_args(argv)++    view = Path(args.data)+    man = load_manifest(view)+    genes = panel_genes(view, man)+    rng = np.random.default_rng(args.seed)++    inputs = sorted(man["inputs"], key=lambda e: e["time"])+    masks = [covered_mask(view, e, genes) for e in inputs]+    Xs = [read_input(view, e, genes) for e in inputs]++    n_last = Xs[-1].shape[0]+    n_out = int(np.clip(n_last, man["min_cells"], man["max_cells"]))+    rows_out = np.sort(rng.choice(n_last, size=n_out, replace=n_out > n_last))++    # HVGs among genes covered by every input stage+    cov_all = np.logical_and.reduce(masks)+    var_tot = np.zeros(len(genes), dtype=np.float64)+    for X in Xs:+        Xc = X[:, cov_all]+        m = np.asarray(Xc.mean(axis=0)).ravel()+        m2 = np.asarray(Xc.multiply(Xc).mean(axis=0)).ravel()+        var_tot[cov_all] += np.maximum(m2 - m * m, 0)+    cand = np.flatnonzero(cov_all)+    order = cand[np.argsort(-var_tot[cand], kind="stable")]+    hvg = np.sort(order[: args.hvg])++    # standardisation stats over all input cells+    S = sparse.vstack([X[:, hvg] for X in Xs]).tocsr()+    mu_g = np.asarray(S.mean(axis=0)).ravel().astype(np.float64)+    sd_g = np.sqrt(np.maximum(np.asarray(S.multiply(S).mean(axis=0)).ravel() - mu_g ** 2, 0))+    sd_g[sd_g < 1e-6] = 1.0+    del S++    A = [np.asarray(X[:, hvg].todense(), dtype=np.float32) for X in Xs]+    A = [((D - mu_g) / sd_g).astype(np.float32) for D in A]+    use_cov = not args.no_covariate+    n_st = len(A)++    def cov_of(si: int, n: int) -> np.ndarray | None:+        if not use_cov:+            return None+        c = np.zeros((n, n_st), dtype=np.float32)+        c[:, si] = 1.0+        return c++    # train (cap cells per stage)+    train_parts, cov_parts = [], []+    for si, D in enumerate(A):+        n = D.shape[0]+        idx = np.arange(n) if n <= args.train_cap else np.sort(rng.choice(n, args.train_cap, replace=False))+        train_parts.append(D[idx])+        cov_parts.append(cov_of(si, len(idx)))+    Xt = np.vstack(train_parts)+    covt = np.vstack([c for c in cov_parts if c is not None]) if use_cov else None+    model, d_in = train_vae(Xt, covt, args.latent, args.kl, args.epochs, args.seed)+    del Xt, covt, train_parts, cov_parts++    Z = [encode_mean(model, d_in, D, cov_of(si, D.shape[0])) for si, D in enumerate(A)]++    scale = args.transition_scale+    if scale < 0:+        if n_st >= 2:+            dt_in = inputs[-1]["time"] - inputs[0]["time"]+            dt_tg = man["target"]["time"] - inputs[-1]["time"]+            scale = float(np.clip(dt_tg / dt_in, 0.0, 2.0)) if dt_in > 0 else 1.0+        else:+            scale = 0.0++    # ---- latent displacement field (needs >=2 stages and scale > 0)+    dz_out = None+    if n_st >= 2 and scale > 0:+        import ot++        Z0f, Z1f = Z[-2], Z[-1]+        n0, n1 = Z0f.shape[0], Z1f.shape[0]+        i0 = np.arange(n0) if n0 <= args.ot_cap else np.sort(rng.choice(n0, args.ot_cap, replace=False))+        i1 = np.arange(n1) if n1 <= args.ot_cap else np.sort(rng.choice(n1, args.ot_cap, replace=False))+        P0, P1 = Z0f[i0], Z1f[i1]+        M = ((P0[:, None, :] - P1[None, :, :]) ** 2).sum(-1)+        M /= max(np.median(M), 1e-8)+        a = np.full(len(i0), 1.0 / len(i0))+        b = np.full(len(i1), 1.0 / len(i1))+        Pi = ot.sinkhorn(a, b, M, args.ot_reg, numItermax=300, stopThr=1e-8)+        Pi = np.nan_to_num(Pi, nan=0.0, posinf=0.0, neginf=0.0)+        rs = Pi.sum(axis=1, keepdims=True)+        rs[rs <= 0] = 1.0+        delta = (Pi / rs) @ P1 - P0                      # per-source-cell displacement+        # kNN-weighted local regression: dz for every last-stage cell (chunked)+        dz_all = np.zeros_like(Z1f)+        k = min(args.k, len(i0))+        for s in range(0, n1, 1024):+            Zq = Z1f[s:s + 1024]+            d2 = ((Zq[:, None, :] - P0[None, :, :]) ** 2).sum(-1)+            knn = np.argpartition(d2, k - 1, axis=1)[:, :k]+            dk = np.take_along_axis(d2, knn, axis=1)+            h = max(np.median(dk[:, -1]), 1e-8)+            w = np.exp(-dk / h)+            w /= w.sum(axis=1, keepdims=True)+            dz_all[s:s + 1024] = np.einsum("ij,ijk->ik", w, delta[knn])+        dz_out = dz_all[rows_out]+        nrm = np.linalg.norm(dz_all, axis=1)+        print(f"[mech] ||dz|| mean={nrm.mean():.4f} cv={nrm.std() / max(nrm.mean(), 1e-9):.4f} "+              f"scale={scale:.3f} ||delta||={np.linalg.norm(delta, axis=1).mean():.4f}", flush=True)++    # ---- predict latents for output cells, decode+    Zp = Z[-1][rows_out].copy()+    if dz_out is not None:+        Zp += scale * dz_out+    Dp = decode(model, Zp.astype(np.float32)) * sd_g + mu_g+    Dp[~np.isfinite(Dp)] = 0.0+    if args.quantile_match:+        # per-gene quantile match: map decoded ranks onto the last stage's+        # empirical marginal (restores sparsity / value distribution per gene)+        E = np.asarray(Xs[-1][:, hvg].todense(), dtype=np.float32)  # n_last x hvg+        E.sort(axis=0)+        q = (np.argsort(np.argsort(Dp, axis=0), axis=0).astype(np.float64) + 0.5) / Dp.shape[0]+        idx = np.minimum((q * E.shape[0]).astype(np.int64), E.shape[0] - 1)+        Dp = np.take_along_axis(E, idx, axis=0).astype(np.float32)+        del E+    np.maximum(Dp, 0.0, out=Dp)+    Dp[Dp < 0.01] = 0.0++    # ---- assemble full-panel output: last stage's values, HVGs blended+    Xout = Xs[-1][rows_out].toarray().astype(np.float32)+    g = float(np.clip(args.blend_gamma, 0.0, 1.0))+    if g > 0:+        Xout[:, hvg] = (1.0 - g) * Xout[:, hvg] + g * Dp+    Xout[~np.isfinite(Xout)] = 0.0++    import anndata as ad+    import pandas as pd+    adata = ad.AnnData(X=sparse.csr_matrix(Xout),+                       obs=pd.DataFrame(index=[f"c{i}" for i in range(Xout.shape[0])]),+                       var=pd.DataFrame(index=list(genes)))+    adata.X.eliminate_zeros()+    adata.write_h5ad(args.out)+    print(f"[out] {adata.shape} nnz/cell={adata.X.nnz / adata.shape[0]:.0f} "+          f"input nnz/cell={Xs[-1].nnz / n_last:.0f}", flush=True)+++if __name__ == "__main__":+    main()

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

用到的知识库条目

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

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

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

分析员

ANALYSIS.json written

工程师

完成。draft 了 generative_latent 家族：Gaussian VAE 潜空间 + entropic-OT 耦合 + kNN 局部回归逐细胞位移（||dz|| CV=0.337，机制非恒定），解码后按基因分位数匹配回真实边缘分布，再以 γ=0.5 与末阶段真实值混合。

查分（X3 A 半，5/20 次）：机制关 28.0 → 纯解码+QM 42.0/43.0 → 默认（γ=0.5, scale=1, kl=0.005）**48.92**（cell_state 68.1, direction 49.2, de_recovery 42.1, covariation 28.6），vec-check 通过，默认输出与打分文件同哈希。

下一步最值得试：用潜位移做真实细胞选择/组成调制（对 z_pred 取最近真实细胞）替代数值混合，可同时保住 cell_state 与 covariation；或把潜位移投回基因空间做保稀疏位移，与节点 4/5/9 的强部件组合。

调研员

lls × 2000 genes should be <3 min GPU; set hard 5-min timeout. 4) Dense HVG output destroys sparsity pattern: mitigated by only replacing HVG genes. 5) With only 2 snapshots, transition direction is weakly identified (k030); keep scale conservative. 6) POT sinkhorn may be slow for large types: cap at 2000 cells/type.",
  "family_id": "generative_latent",
  "mechanism": "Gaussian VAE learns a 20-dim latent representation of input cells; a kNN-weighted local regression in latent space produces per-cell displacement vectors that depend on each cell's latent position (not just its type label), then the decoder maps displaced latents back to correlated gene expression.",
  "vs_constant_shift": "Per-type constant shift adds the same vector to every cell of a type. Here, dz varies within a type based on the cell's latent coordinates and its local neighborhood's observed E8.5→E9.5 displacements. Two cells of the same type at different latent positions receive different dz vectors. The decoder also produces correlated multi-gene changes rather than independent per-gene offsets.",
  "mechanism_evidence": "1) Report within-type coefficient of variation of ||dz||: if >0.3, displacement is genuinely cell-state-dependent. 2) Report per-group scores separately for zero-transition vs full-transition: covariation and de_recovery should improve with transition if the decoder captures gene correlations. 3) Report zero-value rate and library-size distribution of output vs input. 4) Report reconstruction MSE of zero-transition control vs copy_last to verify decoder quality.",
  "mechanism_off_control": "Set --transition-scale 0 (or equivalently skip the displacement step): the program encodes E9.5 cells to latent z and immediately decodes back without any shift. Expected: output ≈ input (autoencoding), scores should be within 2 points of copy_last. If zero-transition output is identical to full-transition output (scale>0), the transition mechanism is not running.",
  "sources": []
}
```