从零构建虚拟细胞与生物标度律：一份硬核工程师的落地指南

# 我用的 Python 3.11，别用 3.12，某些生物信息库还没跟上
conda create -n vcell python=3.11 -y
conda activate vcell

pip install torch==2.4.0 --index-url https://download.pytorch.org/whl/cu121
pip install scanpy==1.10.3 anndata==0.10.8
pip install hydra-core==1.3.2 omegaconf==2.3.0  # 配置管理，后面全靠它
pip install wandb  # 实验追踪，不装你会后悔

# configs/physics_constraints.yaml
scaling_laws:
  cell_size:
    # 表面积/体积比约束：细胞越大，单位体积可用的膜面积越小
    surface_to_volume:
      min_ratio: 0.15
      max_ratio: 8.5
    # 分子拥挤约束：内部大分子体积占比上限
    molecular_crowding:
      max_volume_fraction: 0.40
    # 能量预算：维持一个细胞的最小 ATP 通量
    energy_budget:
      min_atp_flux_mol_per_s: 1.0e-17

  metabolism:
    # Kleiber 律的幂指数，真核生物普遍落在 0.75 附近
    mass_exponent: 0.75
    exponent_tolerance: 0.05

  growth:
    # 生长速率不能超过能量供应允许的上限
    max_specific_growth_rate_per_h: 2.5

validation:
  on_input: true      # 数据进模型前校验
  on_output: true     # 模型输出后裁剪
  strict_mode: false  # false = 报警但不中断，生产环境建议 false
  log_violations: true

# vcell/physics/validator.py
import numpy as np
from dataclasses import dataclass
from omegaconf import DictConfig
import logging

log = logging.getLogger(__name__)

@dataclass
class ConstraintViolation:
    rule: str
    value: float
    bound: float
    severity: str  # "warn" | "error"

class ScalingLawValidator:
    def __init__(self, cfg: DictConfig):
        self.cfg = cfg
        self.violations: list[ConstraintViolation] = []

    def check_cell_size(self, radius_um: float) -> bool:
        """表面积/体积比校验。球体近似，够用了。"""
        if radius_um <= 0:
            self._record("cell_size.positive_radius", radius_um, 0.0, "error")
            return False

        volume = (4/3) * np.pi * radius_um**3
        surface = 4 * np.pi * radius_um**2
        ratio = surface / volume  # 单位: 1/um

        lo = self.cfg.scaling_laws.cell_size.surface_to_volume.min_ratio
        hi = self.cfg.scaling_laws.cell_size.surface_to_volume.max_ratio

        if not (lo <= ratio <= hi):
            self._record("cell_size.surface_to_volume", ratio,
                         lo if ratio < lo else hi, "warn")
            return False
        return True

    def check_metabolism(self, mass_pg: float, metabolic_rate: float) -> bool:
        """Kleiber 律校验：代谢率正比于质量^0.75。"""
        exp = self.cfg.scaling_laws.metabolism.mass_exponent
        tol = self.cfg.scaling_laws.metabolism.exponent_tolerance
        expected = mass_pg ** exp
        # 允许一个宽泛的容差带，别卡太死
        if metabolic_rate <= 0:
            self._record("metabolism.positive_rate", metabolic_rate, 0.0, "error")
            return False
        rel_err = abs(np.log(metabolic_rate) - np.log(expected))
        if rel_err > tol * 10:
            self._record("metabolism.mass_exponent", rel_err, tol*10, "warn")
            return False
        return True

    def _record(self, rule, value, bound, severity):
        v = ConstraintViolation(rule, float(value), float(bound), severity)
        self.violations.append(v)
        if severity == "error":
            log.error("约束违反 [%s]: value=%.4e bound=%.4e", rule, value, bound)
        else:
            log.warning("约束违反 [%s]: value=%.4e bound=%.4e", rule, value, bound)

    def report(self) -> dict:
        return {
            "total": len(self.violations),
            "errors": sum(1 for v in self.violations if v.severity == "error"),
            "warnings": sum(1 for v in self.violations if v.severity == "warn"),
        }

# configs/rules/protein_interaction.yaml
rules:
  - id: "A_promotes_B"
    description: "蛋白 A 促进蛋白 B 的表达（Hill 激活）"
    type: hill_activation
    reactant: "protein_A"
    target: "protein_B"
    params:
      K_half: 0.5      # 半饱和浓度 uM
      hill_n: 2.0      # 协同系数
      v_max: 1.2       # 最大速率 uM/s
    constraints:
      - "v_max <= energy_budget.max_atp_flux_mol_per_s * 1e4"

# train.py 核心片段
import hydra
from vcell.physics.validator import ScalingLawValidator

@hydra.main(config_path="configs", config_name="train")
def main(cfg):
    validator = ScalingLawValidator(cfg.physics_constraints)
    # ... 数据加载、批次校正 ...
    for batch in dataloader:
        # 进模型前校验
        if cfg.physics_constraints.validation.on_input:
            for cell in batch["cells"]:
                validator.check_cell_size(cell["radius_um"])
        preds = model(batch)
        # 出模型后校验
        if cfg.physics_constraints.validation.on_output:
            for p in preds:
                validator.check_metabolism(p["mass_pg"], p["metabolic_rate"])
    print(validator.report())  # {'total': 137, 'errors': 3, 'warnings': 134}