ddi_static.py · Scientific Agent Skills

1#!/usr/bin/env python3

2"""Drug-drug interaction prediction: ICH M12 basic models and the mechanistic static model.

3

4ICH M12 (Step 4, 2024) sets out a stepwise risk assessment: in vitro data feed

5basic models whose cut-offs trigger further evaluation, and a

6mechanistic static or PBPK model can be used to refine a positive basic-model

7signal. The basic models are deliberately conservative — they are designed to

8over-predict under their assumptions. Uncertain inputs or unmodeled mechanisms

9can invalidate a negative screen; a positive screen needs further evaluation.

10

11 python3 ddi_static.py --basic --ki 0.5 --imax 2.0 --fu 0.05 --dose 100

12 python3 ddi_static.py --basic --tdi --ki-inact 1.2 --kinact 0.04 --imax 2.0 --fu 0.05

13 python3 ddi_static.py --msm --ki 0.5 --imax 2.0 --fu 0.05 --dose 100 --fm 0.9 --fg 0.7

14

15All concentrations are in the same molar or mass units as Ki, KI and EC50. The

16script does not know your units and cannot check them; a Ki in micromolar

17against an Imax in ng/mL is a silent, and common, error.

18"""

19

20from __future__ import annotations

21

22import argparse

23import math

24from typing import Sequence

25

26from _common import InputError, Report, add_format_argument, main_wrapper

27

28# ICH M12 basic-model cut-offs.

29CUTOFF_R1_HEPATIC = 1.02

30CUTOFF_R1_GUT = 11.0

31CUTOFF_R2_TDI = 1.25

32CUTOFF_R3_INDUCTION = 0.80

33CUTOFF_TRANSPORTER_HEPATIC_UPTAKE = 1.1

34CUTOFF_TRANSPORTER_INTESTINAL = 10.0

35CUTOFF_TRANSPORTER_RENAL = 0.1

36CUTOFF_TRANSPORTER_MATE = 0.02

37

38# Default intestinal dissolution volume used to form the nominal gut

39# concentration, and the hepatic blood flow used for the inlet concentration.

40GUT_VOLUME_ML = 250.0

41HEPATIC_BLOOD_FLOW_L_H = 97.0

42DEFAULT_KDEG_HEPATIC = 0.0005 # per minute; roughly a 23 h enzyme half-life

43

44

45def r1_reversible(i_conc: float, ki: float) -> float:

46 """Basic reversible inhibition ratio, ``1 + [I]/Ki``."""

47 if ki <= 0:

48 raise InputError("Ki must be positive")

49 return 1.0 + i_conc / ki

50

51

52def r2_time_dependent(i_conc: float, ki_inact: float, kinact: float, kdeg: float) -> float:

53 """Basic TDI ratio, ``(kobs + kdeg) / kdeg``."""

54 if ki_inact <= 0 or kdeg <= 0:

55 raise InputError("KI and kdeg must be positive")

56 kobs = kinact * i_conc / (ki_inact + i_conc)

57 return (kobs + kdeg) / kdeg

58

59

60def r3_induction(i_conc: float, emax: float, ec50: float, scaling: float = 1.0) -> float:

61 """Basic induction ratio; values at or below 0.80 flag a potential inducer."""

62 if ec50 <= 0:

63 raise InputError("EC50 must be positive")

64 return 1.0 / (1.0 + scaling * emax * i_conc / (ec50 + i_conc))

65

66

67def inlet_concentration(imax: float, fu: float, dose: float, fa: float, fg: float, ka: float, blood_ratio: float = 1.0) -> float:

68 """Maximum unbound hepatic inlet concentration.

69

70 ``Iu,inlet,max = fu * (Imax + Fa*Fg*ka*Dose / (Qh*RB))``. This is the

71 estimated inlet concentration under these assumptions, rather than

72 measured systemic Imax and is what M12 asks for in hepatic uptake-transporter

73 assessments.

74 """

75 portal = fa * fg * (ka * 60.0) * dose / (HEPATIC_BLOOD_FLOW_L_H * blood_ratio)

76 return fu * (imax + portal)

77

78

79def mechanistic_static(

80 ih: float,

81 ig: float,

82 fm: float,

83 fg: float,

84 ki: float | None = None,

85 ki_inact: float | None = None,

86 kinact: float | None = None,

87 kdeg_h: float = DEFAULT_KDEG_HEPATIC,

88 kdeg_g: float = 0.0005,

89 ind_emax: float = 0.0,

90 ind_ec50: float | None = None,

91 ind_scaling: float = 1.0,

92) -> dict[str, float]:

93 """Mechanistic static model AUC ratio, combining reversible, TDI and induction terms.

94

95 ``AUCR = 1/(Ag*Bg*Cg*(1-Fg) + Fg) * 1/(Ah*Bh*Ch*fm + (1-fm))``

96

97 The two fractions that dominate the answer are ``fm`` (the fraction of

98 systemic clearance through the affected enzyme) and ``Fg`` (the fraction

99 escaping gut metabolism). Hepatic inhibition alone is limited

100 by ``1/(1-fm)``; including gut inhibition gives ``1/((1-fm)*Fg)``, so an fm assumed at 1.0 when it is

101 really 0.7 changes an unbounded prediction into a 3.3-fold ceiling.

102 """

103 if not all(math.isfinite(x) for x in (ih, ig, fm, fg, kdeg_h, kdeg_g, ind_emax, ind_scaling)):

104 raise InputError("MSM inputs must be finite")

105 if ih < 0 or ig < 0 or not 0 <= fm <= 1 or not 0 < fg <= 1 or min(kdeg_h, kdeg_g) <= 0 or min(ind_emax, ind_scaling) < 0:

106 raise InputError("MSM concentrations/fm/induction must be non-negative, Fg in (0,1], fm <= 1, and kdeg positive")

107 for value in (ki, ki_inact, ind_ec50):

108 if value is not None and (not math.isfinite(value) or value <= 0):

109 raise InputError("Ki, KI and EC50 must be finite and positive when supplied")

110 if (ki_inact is None) != (kinact is None) or (kinact is not None and (not math.isfinite(kinact) or kinact < 0)):

111 raise InputError("TDI requires paired KI and non-negative kinact")

112 def terms(i_conc: float, kdeg: float) -> tuple[float, float, float]:

113 a = 1.0 / (1.0 + i_conc / ki) if ki else 1.0

114 if ki_inact and kinact:

115 b = kdeg / (kdeg + kinact * i_conc / (ki_inact + i_conc))

116 else:

117 b = 1.0

118 c = 1.0 + ind_scaling * ind_emax * i_conc / (ind_ec50 + i_conc) if ind_ec50 else 1.0

119 return a, b, c

120

121 ah, bh, ch = terms(ih, kdeg_h)

122 ag, bg, cg = terms(ig, kdeg_g)

123 gut = 1.0 / (ag * bg * cg * (1.0 - fg) + fg)

124 hepatic = 1.0 / (ah * bh * ch * fm + (1.0 - fm))

125 return {

126 "Ah_reversible": ah,

127 "Bh_tdi": bh,

128 "Ch_induction": ch,

129 "Ag_reversible": ag,

130 "Bg_tdi": bg,

131 "Cg_induction": cg,

132 "gut_component": gut,

133 "hepatic_component": hepatic,

134 "auc_ratio": gut * hepatic,

135 "maximum_possible_auc_ratio": 1.0 / ((1.0 - fm) * fg) if fm < 1 else float("inf"),

136 "hepatic_only_auc_ceiling": 1.0 / (1.0 - fm) if fm < 1 else float("inf"),

137 }

138

139

140def classify(auc_ratio: float) -> str:

141 """FDA/ICH perpetrator classification bands for an AUC ratio."""

142 if auc_ratio >= 5.0:

143 return "strong inhibitor (AUCR >= 5)"

144 if auc_ratio >= 2.0:

145 return "moderate inhibitor (2 <= AUCR < 5)"

146 if auc_ratio >= 1.25:

147 return "weak inhibitor (1.25 <= AUCR < 2)"

148 if auc_ratio > 0.8:

149 return "within the 0.8-1.25 screening band; clinical relevance requires substrate-specific assessment"

150 if auc_ratio > 0.5:

151 return "weak inducer (0.5 < AUCR <= 0.8)"

152 if auc_ratio > 0.2:

153 return "moderate inducer (0.2 < AUCR <= 0.5)"

154 return "strong inducer (AUCR <= 0.2)"

155

156

157# --------------------------------------------------------------------- CLI

158

159

160def build_parser() -> argparse.ArgumentParser:

161 parser = argparse.ArgumentParser(

162 description="ICH M12 basic DDI models and the mechanistic static model.",

163 formatter_class=argparse.RawDescriptionHelpFormatter,

164 epilog=__doc__,

165 )

166 parser.add_argument("--basic", action="store_true", help="run the basic models and report against their cut-offs")

167 parser.add_argument("--msm", action="store_true", help="run the mechanistic static model")

168

169 parser.add_argument("--imax", type=float, help="maximum total systemic plasma concentration of the perpetrator")

170 parser.add_argument("--fu-validated", action="store_true", help="accuracy and precision of fu < 0.01 demonstrated")

171 parser.add_argument("--fu", type=float, default=1.0, help="unbound fraction in plasma (default: 1.0)")

172 parser.add_argument("--dose", type=float, help="perpetrator amount; use micromoles when concentrations are micromolar (amount/L)")

173 parser.add_argument("--ka", type=float, default=0.1, help="perpetrator absorption rate constant, per min (default: 0.1)")

174 parser.add_argument("--fa", type=float, default=1.0, help="fraction absorbed (default: 1.0)")

175 parser.add_argument("--fg-perpetrator", type=float, default=1.0, help="perpetrator gut availability (default: 1.0)")

176 parser.add_argument("--blood-ratio", type=float, default=1.0, help="blood-to-plasma ratio (default: 1.0)")

177

178 parser.add_argument("--ki", type=float, help="reversible inhibition constant")

179 parser.add_argument("--tdi", action="store_true", help="evaluate time-dependent inhibition")

180 parser.add_argument("--ki-inact", type=float, help="KI for time-dependent inactivation")

181 parser.add_argument("--kinact", type=float, help="maximum inactivation rate constant, per min")

182 parser.add_argument("--kdeg", type=float, default=DEFAULT_KDEG_HEPATIC, help="hepatic enzyme degradation rate, per min")

183 parser.add_argument("--ind-emax", type=float, default=0.0, help="maximum fold induction minus one")

184 parser.add_argument("--ind-ec50", type=float, help="induction EC50")

185 parser.add_argument("--ind-scaling", type=float, default=1.0, help="induction scaling factor d (default: 1.0)")

186

187 parser.add_argument("--fm", type=float, help="fraction of victim clearance via the affected enzyme")

188 parser.add_argument("--fg", type=float, default=1.0, help="victim fraction escaping gut metabolism (default: 1.0)")

189 parser.add_argument("--transporter", choices=("hepatic-uptake", "intestinal", "renal", "mate", "systemic-efflux"), help="also run the transporter basic model")

190 parser.add_argument("--transporter-ki", type=float)

191 add_format_argument(parser)

192 return parser

193

194

195def run(argv: Sequence[str] | None = None) -> int:

196 args = build_parser().parse_args(argv)

197 from _common import validate_numeric_args

198 validate_numeric_args(args)

199 if args.basic == args.msm:

200 raise InputError("choose exactly one of --basic or --msm")

201 if args.imax is None:

202 raise InputError("--imax is required")

203 if not 0 < args.fu <= 1:

204 raise InputError("--fu must be in (0, 1]")

205

206 report = Report()

207 for name in ("imax", "dose", "kinact", "ind_emax", "ind_scaling"):

208 value = getattr(args, name)

209 if value is not None and value < 0:

210 raise InputError(f"--{name.replace('_', '-')} must be non-negative")

211 for name in ("ki", "ki_inact", "ind_ec50", "transporter_ki", "ka", "kdeg", "blood_ratio"):

212 value = getattr(args, name)

213 if value is not None and value <= 0:

214 raise InputError(f"--{name.replace('_', '-')} must be positive")

215 if not 0 <= args.fa <= 1 or not 0 <= args.fg_perpetrator <= 1:

216 raise InputError("--fa and --fg-perpetrator must be in [0, 1]")

217 if args.fu < 0.01 and not args.fu_validated:

218 report.finding("fu below 0.01 without demonstrated assay reliability: using M12's 0.01 floor")

219 args.fu = 0.01

220 if args.tdi and (args.ki_inact is None or args.kinact is None):

221 raise InputError("--tdi needs --ki-inact and --kinact")

222 if (args.ind_emax > 0) != (args.ind_ec50 is not None):

223 raise InputError("supply --ind-emax > 0 together with --ind-ec50")

224 unbound = args.imax * args.fu

225 gut = (args.dose / (GUT_VOLUME_ML / 1000.0)) if args.dose is not None else None

226 report.scalar("imax_total", args.imax)

227 report.scalar("imax_unbound", unbound)

228 if gut is not None:

229 report.scalar("gut_concentration", gut)

230 report.note(

231 "Units are not checked. Ki, KI, EC50, Imax and dose must all be expressed on the same molar "

232 "basis; the gut concentration is dose divided by 250 mL."

233 )

234

235 if args.basic:

236 rows = []

237 if args.ki:

238 r1 = r1_reversible(unbound, args.ki)

239 rows.append(

240 {

241 "model": "reversible inhibition, hepatic",

242 "value": r1,

243 "cutoff": f">= {CUTOFF_R1_HEPATIC}",

244 "triggers_study": r1 >= CUTOFF_R1_HEPATIC,

245 "basis": "1 + Imax,u/Ki",

246 }

247 )

248 if gut is not None:

249 r1g = r1_reversible(gut, args.ki)

250 rows.append(

251 {

252 "model": "reversible inhibition, intestinal (CYP3A4)",

253 "value": r1g,

254 "cutoff": f">= {CUTOFF_R1_GUT}",

255 "triggers_study": r1g >= CUTOFF_R1_GUT,

256 "basis": "1 + Igut/Ki, Igut = dose/250 mL",

257 }

258 )

259 if args.tdi:

260 if args.ki_inact is None or args.kinact is None:

261 raise InputError("--tdi needs --ki-inact and --kinact")

262 r2 = r2_time_dependent(unbound * 5.0, args.ki_inact, args.kinact, args.kdeg)

263 rows.append(

264 {

265 "model": "time-dependent inhibition, hepatic",

266 "value": r2,

267 "cutoff": f">= {CUTOFF_R2_TDI}",

268 "triggers_study": r2 >= CUTOFF_R2_TDI,

269 "basis": "(kobs + kdeg)/kdeg at 5 x Imax,u",

270 }

271 )

272 if args.ind_ec50 and args.ind_emax:

273 r3 = r3_induction(unbound * 10.0, args.ind_emax, args.ind_ec50, args.ind_scaling)

274 rows.append(

275 {

276 "model": "induction",

277 "value": r3,

278 "cutoff": f"<= {CUTOFF_R3_INDUCTION}",

279 "triggers_study": r3 <= CUTOFF_R3_INDUCTION,

280 "basis": "1/(1 + d*Emax*I/(EC50+I)) at 10 x Imax,u",

281 }

282 )

283 if args.transporter:

284 if not args.transporter_ki:

285 raise InputError("--transporter needs --transporter-ki")

286 if args.transporter == "hepatic-uptake":

287 if args.dose is None:

288 raise InputError("hepatic uptake needs --dose for the inlet concentration")

289 conc = inlet_concentration(args.imax, args.fu, args.dose, args.fa, args.fg_perpetrator, args.ka, args.blood_ratio)

290 value = 1.0 + conc / args.transporter_ki

291 cutoff, triggers = CUTOFF_TRANSPORTER_HEPATIC_UPTAKE, value >= CUTOFF_TRANSPORTER_HEPATIC_UPTAKE

292 basis = "1 + Iu,inlet,max/Ki,u (OATP1B1/1B3)"

293 elif args.transporter == "intestinal":

294 if gut is None:

295 raise InputError("intestinal transporter assessment needs --dose")

296 value = gut / args.transporter_ki

297 cutoff, triggers = CUTOFF_TRANSPORTER_INTESTINAL, value >= CUTOFF_TRANSPORTER_INTESTINAL

298 basis = "Igut/IC50 (P-gp, BCRP)"

299 else:

300 value = unbound / args.transporter_ki

301 cutoff = CUTOFF_TRANSPORTER_MATE if args.transporter in {"mate", "systemic-efflux"} else CUTOFF_TRANSPORTER_RENAL

302 triggers = value >= cutoff

303 basis = "Imax,u/IC50,u (MATE1/2-K or systemic P-gp/BCRP)" if cutoff == CUTOFF_TRANSPORTER_MATE else "Imax,u/IC50,u (OAT1/3, OCT2)"

304 rows.append(

305 {

306 "model": f"transporter inhibition, {args.transporter}",

307 "value": value,

308 "cutoff": f">= {cutoff}",

309 "triggers_study": triggers,

310 "basis": basis,

311 }

312 )

313 if not rows:

314 raise InputError("--basic needs at least one of --ki, --tdi, or --ind-ec50 with --ind-emax")

315 report.table("ICH M12 basic models", rows)

316 for row in rows:

317 if row["triggers_study"]:

318 report.finding(

319 f"{row['model']}: {row['value']:.4g} meets the cut-off {row['cutoff']} - a clinical "

320 "DDI study or a refined mechanistic/PBPK assessment is indicated"

321 )

322 report.note(

323 "Basic models are intentionally conservative. A result below the cut-off supports not "

324 "studying the interaction; a result above it is a trigger for further evaluation, not a "

325 "prediction of clinical magnitude."

326 )

327

328 else:

329 if args.fm is None:

330 raise InputError("--msm needs --fm")

331 if not 0 < args.fm <= 1:

332 raise InputError("--fm must be in (0, 1]")

333 if not 0 < args.fg <= 1:

334 raise InputError("--fg must be in (0, 1]")

335 if args.dose is None:

336 raise InputError("--msm needs --dose (use 0 for a systemic-only scenario)")

337 ih = inlet_concentration(args.imax, args.fu, args.dose, args.fa, args.fg_perpetrator, args.ka, args.blood_ratio)

338 ig = args.fa * args.ka * 60.0 * args.dose / 18.0

339 report.scalar("msm_hepatic_inlet_unbound", ih)

340 report.scalar("msm_enterocyte_concentration", ig)

341 result = mechanistic_static(

342 ih=ih,

343 ig=ig,

344 fm=args.fm,

345 fg=args.fg,

346 ki=args.ki,

347 ki_inact=args.ki_inact if args.tdi else None,

348 kinact=args.kinact if args.tdi else None,

349 kdeg_h=args.kdeg,

350 ind_emax=args.ind_emax,

351 ind_ec50=args.ind_ec50,

352 ind_scaling=args.ind_scaling,

353 )

354 report.table("mechanistic static model", [{"term": k, "value": v} for k, v in result.items()])

355 report.scalar("predicted_auc_ratio", result["auc_ratio"])

356 report.scalar("classification", classify(result["auc_ratio"]))

357 report.note(

358 f"With fm = {args.fm} and Fg = {args.fg}, the combined hepatic/gut inhibition ceiling is "

359 f"{result['maximum_possible_auc_ratio']:.2f}-fold. If the prediction approaches that "

360 "ceiling, fm is doing more work than the inhibition constants."

361 )

362 if gut is None:

363 report.note("no --dose given, so the intestinal component was evaluated at zero inhibitor concentration")

364 if args.fm > 0.95:

365 report.finding(

366 f"fm = {args.fm} implies almost all clearance goes through one enzyme. This is rarely "

367 "established and it drives the prediction; state the evidence for it or run a sensitivity "

368 "analysis across a plausible range."

369 )

370 if result["auc_ratio"] >= 2.0:

371 report.finding(

372 f"predicted AUC ratio {result['auc_ratio']:.2f} - {classify(result['auc_ratio'])}; "

373 "this screening prediction needs confirmation before clinical or labelling conclusions"

374 )

375 report.note(

376 "The mechanistic static model assumes a single constant perpetrator concentration and no "

377 "time course. It is a screening refinement; where the interaction is decision-relevant, "

378 "ICH M12 points to PBPK with a verified perpetrator model."

379 )

380

381 if args.msm and args.ind_emax and (args.ki is not None or args.tdi):

382 report.finding("Run inhibition-only and induction-only scenarios separately as well: combined effects can mask either mechanism (ICH M12 7.5.1.1)")

383 return report.emit(args.format)

384

385

386if __name__ == "__main__":

387 raise SystemExit(main_wrapper(run))