虚拟器官插件开发教程（18）：实战二·肝安全插件——DILIsym 风格机制模型的剂量窗口实战

虚拟器官插件开发教程（18）：实战二·肝安全插件——DILIsym 风格机制模型的剂量窗口实战

版本声明块

• 工具/软件：Python 3.11 + scipy（solve_ivp/LSODA）；机制上游=第 12 篇自制 ODE 网络；暴露上游=第 11 篇 PBPK CSV 契约（PK-Sim 开源路线）
• 商业对照：DILIsym®（Simulations Plus，商业闭源，模型 v6A/软件 DILIsym 11，2025-05；官方称已支持 38 个化合物申报、FDA 连续多年续约）——仅作方法参照，本文不臆造其任何接口或内部公式
• 本文目标：交付 DILIsym风格的自建肝安全插件 hepatox，端到端跑通"暴露→机制激活→肝细胞损失→安全窗口"，并完成官方案例的复现度声明

⚠️诚实边界（铁律 7）：调研素材库中不存在开源的 DILIsym 复现实现，本插件是自建的机制 ODE 网络；全部速率常数为合成校准值，仅复现"行为模式"，与 Simulations Plus 的任何数值结果无关。

一句话结论：hepatox 实测 0.2 秒跑完 72 h 机制仿真——APAP 式高暴露情景（合成 Cmax=200 µM、MMP IC50=80 µM、nH=2、ROS 生成项）判 DILI+（损失 1.0、ATP 谷值 1.18），治疗暴露参照判 safe（损失≈0），暴露扫描给出演示阈值loss<0.10下的安全窗口上界≈50 µM（75 µM 起损、150 µM 近饱和）；fezolinetant（CPT, doi:10.1002/cpt.70194）与 CGRP 类 vs telcagepant（PMC9237996）两个官方案例，我们只复现其"决策结构"，不复现商业内核。

〇、本篇要解决的认知问题

• 一个肝安全插件的输入输出契约到底是什么，和第 11/12 篇怎么接线？
• 闭源商业平台（DILIsym）摆在面前，插件开发者能"复现"的边界在哪里？
• fezolinetant 剂量选择、CGRP 类药 vs telcagepant 这两个官方案例，抽象成可回归的"决策结构"是什么样？
• 机制仿真输出（ATP/ROS/胆汁酸/死亡分流）怎么变成一句话能懂的"安全剂量窗口"？
• hepatox 怎么登记进平台（第 20 篇伏笔），和心脏插件共用哪套外壳？

一、机制解析：从两份上游契约到一个 ODE 网络

第 11 篇交付了暴露输入（PBPK 的 CSV 契约：时间-浓度点列 + 游离分数），第 12 篇交付了最小机制网络（ATP/ROS/MPT/胆汁酸 ODE + 自适应抗氧化）。本篇把两者拧成一个可安装插件hepatox，管线是：

PBPK CSV(第11篇) 体外机制筛查(第11篇 schema) 机制网络(第12篇) 决策输出 ┌───────────────┐ ┌──────────────────────────┐ ┌─────────────────────┐ ┌──────────────┐ │ Cmax/AUC/时程 │→→│ Bsep/Mrp2/NTCP 抑制 IC50 │→→│ MMP抑制→ATP耗竭→坏死 │→→│ 损失时程曲线 │ │ fraction_unbound│ │ MMP 抑制 IC50 / ROS 生成倾向│ │ ROS/BA 超载→凋亡 │ │ 安全剂量窗口 │ └───────────────┘ └──────────────────────────┘ └─────────────────────┘ └──────────────┘

对照 DILIsym 官方公开的机制模块（线粒体功能障碍/胆汁酸转运抑制/氧化应激/免疫介导四模块，另有 PBPK 与胆汁酸建模扩展），我们实现了前三支的最小代理；第四支（免疫介导，对应 DILIsym 11 的 CD8+ T 细胞模型）素材只给了模块名、没给可核实机制细节，本篇不实现也不臆测，留作扩展位。这正体现铁律 7 的操作方式：商业平台告诉我们"有哪几支通路"是公开信息，"每支怎么写方程"必须自己标定并声明。

复现度声明：两个官方案例我们抄什么、不抄什么。

telcagepant 是肝侧案例药物（因肝毒性退出），不要与心脏侧 CiPA 的药名混用——两个系列各走各的参考集。做肝侧插件，行业事实还有：QST（Quantitative Systems Toxicology，定量系统毒理）是 QSP 在安全性方向的分支、DILIsym 自我定位即 QST 平台；开源 ML 类（DILIPredictor/BioGLGCN-DILI/DILI Map）是分类器不是机制模型，与本篇机制网络互补而非同类（第 11 篇有定位区分表）。

二、完整代码与逐行剖析

以下三段拼接即完整hepatox.py，本文实测可跑（全程 0.21 s，比心脏插件快两个数量级——单 ODE 网络无刚性 AP 循环，这就是第 19 篇分层漏斗"肝侧天然适合大批量"的实证）。

代码 1：输入契约与 Hill（两份上游的接口固化）

# -*- coding: utf-8 -*-"""hepatox —— DILIsym 风格肝安全插件（教学演示；自制机制 ODE，非复现商业内核）"""importnumpyasnpfromdataclassesimportdataclassfromscipy.integrateimportsolve_ivp@dataclassclassHepIn:# 输入契约（铁律1）：机制命中 + 暴露，单位全部写死cmax_uM:float# 峰暴露：真实项目由第11篇 PBPK CSV 读入，此处直填演示auc_uM_h:float# 暴露量 AUC（uM·h）ic50_bsep_uM:float# Bsep 抑制（胆汁酸输出转运体）ic50_mrp2_uM:float# Mrp2 抑制（结合型胆盐排泌）ic50_ntcp_uM:float=1e9# NTCP 摄取抑制：多数药不抑 NTCP -> 默认"无穷大IC50=不命中"mmp_ic50_uM:float=1e9# 线粒体膜电位(MMP)抑制 IC50，默认不命中ros_k:float=0.0# 氧化应激生成倾向系数（0=该支路关断）fraction_unbound:float=1.0# 游离分数 fu：肝内游离浓度 = C*fu（与心脏插件同一约定）conc_unit:str="uM"# 单位进接口：nM/uM 混用是跨组交接第一大事故源defhill(C,ic50,nH=1.0):ifC<=0:return0.0returnC**nH/(ic50**nH+C**nH)

为什么缺省值全是1e9而不是None？缺省即"该机制不命中"（hill→0），新化合物只填测过的靶点、没测的自动走安全侧——这比 None 判断少一类分支 bug，也让"哪些靶点没数据"在申报档案（第 16 篇 inputs 表）里一目了然。ros_k是系数不是 IC50：素材里"ROS/氧化应激筛查"没给可核实的浓度响应形式，所以本篇把它做成生成项幅度参数并标注合成——宁可粗，不可伪。

代码 2：最小机制网络（第 12 篇方程的插件化封装）

defdilisym_style_ode(hin:HepIn):"""状态 y = [ATP, ROS, BA_cell, necrotic, apoptotic]，时间单位小时。 速率常数按'每小时'标定——校准要点：死亡是慢过程，72h 内不应全线饱和， 若扫描暴露后所有剂量都 loss≈1.0，一定是速率/时间尺度错了（本篇实测踩过）。"""Vmax_ATP=8.0# 线粒体 ATP 产出基准（标度值，非生理绝对值）k_ATP_cons=2.0# 基础消耗defrhs(t,y):ATP,ROS,BA,nec,apop=y C=hin.cmax_uM*np.exp(-np.log(2)*t/3.0)# 演示用一室指数衰减(半衰期3h)；Cfree=C*hin.fraction_unbound# 生产直接插值 PBPK CSV 时程(第11篇)mmp_block=hill(Cfree,hin.mmp_ic50_uM,2.0)# nH=2：线粒体毒性常呈陡的超剂量响应ntcp_in=hill(Cfree,hin.ic50_ntcp_uM)# NTCP 抑制反而降摄取=保护项(方向别反!)bsep_out=hill(Cfree,hin.ic50_bsep_uM)mrp2_out=hill(Cfree,hin.ic50_mrp2_uM)dBA=1.2*(1.0-ntcp_in)-(0.8*(1.0-bsep_out)+0.5*(1.0-mrp2_out))*max(BA,0.0)# 氧化应激：生成=药源项+线粒体泄漏项；清除力随抗氧化物耗竭而衰减(Nrf2 类比自适应,第12篇)ros_gen=hin.ros_k*hill(Cfree,60.0,2.0)+0.3*mmp_block ros_scav=0.6*max(1.0-ROS*0.25,0.1)*ROS prod=Vmax_ATP*(1.0-0.9*mmp_block)*(1.0-0.4*np.clip(ROS,0,1))# 双途径打击线粒体dATP=prod-k_ATP_cons*ATP dROS=ros_gen-ros_scav# 损伤分流：ATP 跌破阈值->坏死；ROS/BA 超阈->凋亡（系数已按小时标定）nec_rate=0.4*np.clip((3.0-ATP)/3.0,0,None)apop_rate=0.1*max(ROS-0.9,0)+0.15*max(BA-4.0,0)dNec=nec_rate*max(1.0-nec-apop,0)# 存活池收缩：死亡不会超过 100% 细胞dApo=apop_rate*max(1.0-nec-apop,0)return[dATP,dROS,dBA,dNec,dApo]returnrhsdefrun_hepatox(hin:HepIn,t_end=72.0,loss_thr=0.10):sol=solve_ivp(dilisym_style_ode(hin),[0,t_end],[5.0,0.05,1.0,0.0,0.0],method="LSODA",rtol=1e-6,atol=1e-9,max_step=0.5,dense_output=True)# LSODA：坏死阈值处 RHS 不可导，BDF/显式会抖ifnotsol.success:# 求解器失败必须上浮（第19篇落 fail_reason）returndict(ok=False,fail=sol.message)tt=np.linspace(0,t_end,300);ATP,ROS,BA,nec,apop=sol.sol(tt)loss=np.clip(nec+apop,0,1.0)# 肝细胞损失代理（ALT 是其滞后的湿测投影）returndict(ok=True,t=tt,ATP=ATP,ROS=ROS,BA=BA,loss=loss,alt_peak=float(loss[-1]),atp_nadir=float(ATP.min()),ba_auc=float(np.trapezoid(BA,tt)),safe=bool(loss[-1]<loss_thr))

踩坑实录（写进代码注释也写进正文）：初版把nec_rate系数取 1.5/h、ROS 生成阈值取 10 µM，结果暴露扫描里连 10 µM 都判 DILI+、所有剂量 loss=1.0——72 h 积分窗口下速率放大了死亡通量。修正=nH 提到 2 让剂量响应变陡 + 死亡速率降到 0.4/h 与 0.1~0.15/h。机制模型的"定性正确"一半靠速率常数的时间尺度标定，这类错误 ML 分类器不会犯（没有时间轴），恰是 QST 管线独有的调试科目。另一个方向性陷阱：NTCP 抑制减少胆汁酸摄取，在胆汁淤积叙事里是保护项，初学常把它写成毒性放大器——ntcp_in只进dBA的生成端，别的地方都不许引用它。

代码 3：三情景回归 + 安全窗口扫描（实测输出）

if__name__=="__main__":# —— APAP 式高暴露情景：合成参数模拟"对乙酰氨基酚超量"行为模式(线粒体+氧化应激主导)，# 非实测 NAPQI 动力学；DILIsym 官方案例语境(APAP 缓释过量救治指南支持)仅作叙事参照 ——apap=HepIn(cmax_uM=200.0,auc_uM_h=600.0,ic50_bsep_uM=1e9,ic50_mrp2_uM=1e9,mmp_ic50_uM=80.0,ros_k=2.0,fraction_unbound=0.7)chole=HepIn(cmax_uM=5.0,auc_uM_h=30.0,ic50_bsep_uM=1.0,ic50_mrp2_uM=2.0,mmp_ic50_uM=1e9,ros_k=0.1,fraction_unbound=0.9)# 胆汁淤积型(机制强暴露低)safe=HepIn(cmax_uM=2.0,auc_uM_h=8.0,ic50_bsep_uM=50.0,ic50_mrp2_uM=80.0,mmp_ic50_uM=1e9,ros_k=0.05,fraction_unbound=0.5)# 治疗暴露参照forname,hinin[("APAP-like",apap),("cholestatic",chole),("safe-ref",safe)]:r=run_hepatox(hin)print("%-12s loss=%.3f ATPnadir=%.2f BA_auc=%.1f safe=%s"%(name,r["alt_peak"],r["atp_nadir"],r["ba_auc"],r["safe"]))print("--- 暴露扫描 -> 安全剂量窗口(Cmax 半对数网格) ---")thr=0.10;last_safe=0.0forcmaxin[10,25,50,75,100,150,200,300]:h=HepIn(cmax_uM=float(cmax),auc_uM_h=float(cmax)*3,ic50_bsep_uM=1e9,ic50_mrp2_uM=1e9,mmp_ic50_uM=80.0,ros_k=2.0,fraction_unbound=0.7)o=run_hepatox(h)ifo["alt_peak"]<thr:last_safe=cmaxprint(" Cmax=%4d uM -> loss=%.3f %s"%(cmax,o["alt_peak"],"safe"ifo["alt_peak"]<threlse"DILI+"))print("安全窗口(演示阈值 loss<%.2f): Cmax <= %d uM"%(thr,last_safe))assertrun_hepatox(apap)["safe"]isFalseandrun_hepatox(safe)["safe"]isTrueassert0<last_safe<300,"窗口必须有界(两端都不敏感=标定失败)"print("selfcheck PASS")

实测输出：

APAP-like loss=1.000 ATPnadir=1.18 BA_auc=66.5 safe=False cholestatic loss=0.000 ATPnadir=3.98 BA_auc=78.0 safe=True safe-ref loss=0.000 ATPnadir=3.98 BA_auc=66.6 safe=True --- 暴露扫描 -> 安全剂量窗口(Cmax 半对数网格) --- Cmax= 10 uM -> loss=0.000 safe Cmax= 25 uM -> loss=0.000 safe Cmax= 50 uM -> loss=0.000 safe Cmax= 75 uM -> loss=0.269 DILI+ Cmax= 100 uM -> loss=0.592 DILI+ Cmax= 150 uM -> loss=0.948 DILI+ Cmax= 200 uM -> loss=1.000 DILI+ Cmax= 300 uM -> loss=1.000 DILI+ 安全窗口(演示阈值 loss<0.10): Cmax <= 50 uM selfcheck PASS

三情景的形状差异比结论重要：APAP 型由线粒体+ROS 双通路驱动、ATP 谷值 1.18 直接触发坏死分流；胆汁淤积型BA_auc=78三组最高却判 safe——胆汁酸负荷升了但没越过凋亡率阈值，这正是"决策结构"层面教给注册审评的话：机制激活 ≠ 损伤，要有阈值-时程联读。窗口扫描给出 50/75 µM 之间的拐点，对应 fezolinetant 案例那种"在疗效暴露与预测损伤拐点之间选剂量"的论证骨架。

注册进平台（第 20 篇伏笔）：hepatox 对平台的可见面只有三样——HepIn的 JSON schema、run_hepatox()入口、{alt_peak, atp_nadir, ba_auc, safe, fail}输出字典，与 cardiotox 共用同一 Task 外壳；source字段一路带到审计表。

三、常见报错与排查

• 现象：所有剂量 loss≈1.0 或全 0，扫描无拐点。根因：速率常数时间尺度与积分窗不匹配（本篇实测踩过 1.5/h 版死亡项）。解法：把t_end缩到 12 h 重跑，若形态不变则网络太刚性；拐点至少出现在网格中间 1/3 才算标定成功。
• 现象：solve_ivp报 stiffness 或结果锯齿。根因：坏死阈值clip((3-ATP)/3,0,..)处 RHS 不光滑。解法：LSODA/BDF +max_step限步（本篇 0.5 h），或把硬阈值换 sigmoid 软化。
• 现象：胆汁酸组dBA发散。根因：Bsep 与 Mrp2 同时 100% 抑制时ba_out→0，dBA 退化为常速率线性积累。解法：给BA加生理上限项或在档案里标注"完全双抑制=非生理外推区"。
• 现象：NTCP 抑制越抑越毒。根因：把ntcp_in同时写进了摄取和外排两端。解法：NTCP 只进生成端（保护向），Bsep/Mrp2 只进出排端，画一张转运体-方程项对照表贴 README。
• 现象：文档里出现"DILIsym 的 API/公式"。根因：把商业闭源当开源写。解法：全文只允许"DILIsym 风格的自建模型"表述；其内部实现以厂商文档为准/需向厂商索取（铁律 7）。

四、动手练习

1. 把chole的ic50_bsep_uM从 1.0 改到 0.2（更强抑制）重跑。判定标准：ba_auc显著上升（>78.0 基线），loss 可能突破 0.10 判 DILI+——观察胆汁酸通路何时从"负荷"翻成"损伤"。
2. 将扫描网格auc_uM_h改为float(cmax)*3→float(cmax)*1.5（半衰期缩短）。判定标准：窗口上界右移（同 loss 阈值下允许更高 Cmax），理解"暴露形态与效力谁主导窗口"。
3. 给safe情景加"ros_k": 0.9。判定标准：ATPnadir变化不大但 ROS 时程抬升、loss 仍 <0.10——复现"生物标志物动了、损伤没动"的亚临床形态。

五、小结与下一篇预告

本篇交付了肝侧端到端插件 hepatox：两份上游契约（PBPK 暴露 + 机制筛查 IC50）→ 自制五状态 ODE 网络 → 损失时程与安全窗口，三情景回归 + 有界性断言全部实测通过；对 fezolinetant 与 CGRP/telcagepant 两个官方案例，明确了"复现决策结构、不复现商业内核"的边界声明——这句话本身就是申报材料里限制声明一节（第 16 篇）的正文。心脏与肝两个插件已具备同一外壳，但 100 倍速差（40 s vs 0.2 s）暴露出管线设计的真问题：单细胞电生理怎么扩到万级化合物？GPU 和 MPI 在哪一层介入？审计怎么统一？这些是下一篇（第 19 篇）生产化主题，也是第 20 篇平台化的地基。

本篇认知问题回显（FAQ）

Q1：肝安全插件 hepatox 的输入输出契约是什么？
A：输入=机制命中参数（Bsep/Mrp2/NTCP 抑制 IC50、线粒体 MMP 抑制 IC50、ROS 生成倾向，均对游离浓度、单位 µM）+ PBPK 暴露（Cmax/AUC/时程与游离分数 fu，第 11 篇 CSV 契约）；输出=肝细胞损失时程（坏死+凋亡分数）、ATP 谷值、胆汁酸 AUC、safe/DILI+ 判定与安全剂量窗口，另带 fail 原因字段。

Q2：DILIsym 是闭源的，肝插件"复现"的边界在哪里？
A：公开可核实的只有事实与结构——归属 Simulations Plus（2017 收购 DILIsym Services）、模型 v6A、DILIsym 11（2025-05，儿科+CD8+ T 细胞）、四大机制模块、38 化合物申报与 FDA 续约。内部方程与校准集不公开且无开源复现（调研未查到），所以教程写"DILIsym 风格自建机制模型"并全参数标注合成；细节以厂商文档为准/需向厂商索取。

Q3：fezolinetant 与 CGRP vs telcagepant 案例在本教程里被复现到什么程度？
A：只复现决策结构：前者=多剂量档暴露→预测肝细胞损失→对照疗效选剂量（CPT doi:10.1002/cpt.70194 的论证骨架）；后者=同类化合物批量过机制模型、用已知肝毒性药（telcagepant，PMC9237996 语境）当阳性对照支持类内换药。不复现商业内核参数与原报告数值，本教程输出与厂商结果无对应关系。

Q4：安全剂量窗口是怎么从机制仿真里算出来的？
A：对同一化合物把暴露 Cmax 在半对数网格上扫描，逐点积分 72 h 机制 ODE 得终末损失分数，取"损失<演示阈值（如 10%）"的最大暴露为窗口上界；实测合成情景给出 50 µM 拐点（75 µM 起损、150 µM 近饱和），并要求窗口有界（两端全 safe 或全 DILI 都判标定失败）。

Q5：肝插件与心脏插件如何共用同一平台外壳？
A：两者都对平台暴露三样东西：输入 dataclass 的 JSON schema（HepIn/DrugProfile）、单一入口函数（run_hepatox()/assess()）、统一输出字典（含 safe 或 tdP_tier 与 fail 字段），并共用 Task 哈希键（化合物|浓度|模型版本|种子）作为幂等与审计主键；批量引擎按 organ 字段分发插件入口，许可元数据随注册表行携带，外层无需为不同器官改代码。