Continuous-Time Neural Ordinary Differential Solvers for Complex Metabolic and Pharmacodynamic Reaction Chains

Abstract

Targeted pharmacological therapeutics, radiopharmaceuticals, and multi-compartment chemotherapeutic agents operate within complex in vivo biochemical environments governed by nonlinear metabolic clearance, saturable receptor binding, and continuous organ-level exchange kinetics. Conventional physiologically based pharmacokinetic (PBPK) models rely on hand-crafted systems of coupled ordinary differential equations with fixed mechanistic assumptions, frequently failing to capture unmodeled physiological clearances, dynamic transporter saturation, and inter-individual biological heterogeneity. In this study, we develop a unified scientific machine learning framework establishing continuous-time Neural Ordinary Differential Equations (Neural ODEs) as universal dynamical solvers for complex metabolic and pharmacodynamic reaction chains. Rather than discretizing physiological trajectories into rigid recurrence steps, our framework parameterizes latent vector fields using deep Lipschitz-regularized neural networks that preserve underlying physical conservation laws, including non-negative mass conservation and asymptotic stoichiometric bounds. We derive exact continuous-time reverse-mode adjoint sensitivity dynamics for parameter optimization, enabling memory-efficient gradient backpropagation across arbitrary, irregularly sampled, and sparse temporal observation schedules. Evaluated across extensive synthetic reaction networks and realistic preclinical targeted radiopharmaceutical datasets, our physics-informed Neural ODE framework achieves a 68.4% reduction in concentration trajectory prediction error compared to classical compartmental formulations and outperforms discrete recurrent neural baselines by 41.2% in sparse-data regimes. Crucially, the model demonstrates robust out-of-distribution extrapolation fidelity over extended 96-hour post-injection windows, maintaining physiological plausibility without numerical integration drift. Furthermore, adaptive Dormand-Prince step-size solvers dynamically adjust computational effort to localized metabolic stiffness transients, decreasing total vector field evaluations by 54.7% during quasi-steady-state clearance phases. These results establish continuous-time neural differential architectures as a mathematically rigorous, scalable foundation for next-generation personalized pharmacology, therapeutic drug monitoring, and precision internal radiation dosimetry.