Research
2026
- Bayesian Bilevel Operator Learning with Low-Rank Adaptation for Efficient Uncertainty Quantification of PDE Inverse ProblemsRay Zirui Zhang, Christopher E. Miles, Xiaohui Xie, and John S. LowengrubNature Communications, Sep 2026
Uncertainty quantification in PDE inverse problems is essential in many applications. Scientific machine learning and AI enable data-driven learning of model components while preserving physical structure, and provide the scalability and adaptability needed for emerging imaging technologies and clinical insights. We develop a Bilevel Local Operator Learning framework for Bayesian inference in PDEs (B-BiLO). At the upper level, we sample parameters from the posterior via Hamiltonian Monte Carlo, while at the lower level we fine-tune a neural network via low-rank adaptation (LoRA) to approximate the solution operator locally. B-BiLO enables efficient gradient-based sampling without synthetic data or adjoint equations and avoids sampling in high-dimensional weight space, as in Bayesian neural networks, by optimizing weights deterministically. We analyze errors from approximate lower-level optimization and establish their impact on posterior accuracy. Numerical experiments across PDE models, including tumor growth, demonstrate that B-BiLO achieves accurate and efficient uncertainty quantification.
- PREDICT-GBM: A Multicenter Platform Advancing Personalized Glioblastoma Radiotherapy PlanningLucas Zimmer, Jonas Weidner, Michal Balcerak, Florian Kofler, Mara Krupa, Ivan Ezhov, Santiago Cepeda, Ray Zirui Zhang, John S. Lowengrub, Bjoern Menze, and Benedikt Wiestlernpj Digital Medicine, Sep 2026
Glioblastoma recurrence commonly extends beyond the visible tumor boundary, and standard radiotherapy relies on uniform margins that ignore individual biological and anatomical variation. PREDICT-GBM is an open-source, multi-center platform that pairs a large longitudinal patient dataset with a standardized validation pipeline to develop and benchmark personalized radiotherapy planning approaches, including biophysical tumor growth models, enabling systematic comparison of modeling strategies for improved recurrence coverage.
- Identifiability Limits of Physics-Informed Inference for Spatial Stochastic Dynamics from Static SnapshotsRujie Gu, Ray Zirui Zhang, and Christopher E. MilesJul 2026
Despite increasing scale and resolution, many biological measurements remain destructive, revealing only spatial information rather than the dynamics it encodes. By combining flexible representations with mechanistic constraints, physics-informed machine learning offers a promising route to inferring these dynamics from static snapshots. Motivated by subcellular imaging of gene expression, we ask when a static spatial pattern of molecules can identify spatially varying diffusivity, creation, destruction, and boundary exchange, and how different inference schemes perform on the task. A structural identifiability analysis shows that distributed sources are non-identifiable, whereas a point source such as a transcription site can restore identifiability. These limits are further shaped by seemingly innocuous modeling choices: the boundary conditions, the spatial regularity of the underlying dynamics, and even the stochastic calculus convention. We then adapt several physics-informed schemes, differing in how they represent the solution and enforce the governing equations, and demonstrate effective inference from a single snapshot. Physics-informed approaches can thus recover spatial heterogeneities of biological dynamics from static data, but their use should be accompanied and guided by careful identifiability analysis for meaningful interpretation of the results.
- BiLO: Bilevel Local Operator Learning for PDE Inverse ProblemsRay Zirui Zhang, Christopher E. Miles, Xiaohui Xie, and John S. LowengrubJournal of Computational Physics, Apr 2026
We propose a new neural network based method for solving inverse problems for partial differential equations (PDEs) by formulating the PDE inverse problem as a bilevel optimization problem. At the upper level, we minimize the data loss with respect to the PDE parameters. At the lower level, we train a neural network to locally approximate the PDE solution operator in the neighborhood of a given set of PDE parameters, which enables an accurate approximation of the descent direction for the upper level optimization problem. The lower level loss function includes the least-square penalty of both the residual and its derivative with respect to the PDE parameters. We apply gradient descent simultaneously on both the upper and lower level optimization problems, leading to an effective and fast algorithm. The method, which we refer to as BiLO (Bilevel Local Operator learning), is also able to efficiently infer unknown functions in the PDEs through the introduction of an auxiliary variable. We provide a theoretical analysis that justifies our approach. Through extensive experiments over multiple PDE systems, we demonstrate that our method enforces strong PDE constraints, is robust to sparse and noisy data, and eliminates the need to balance the residual and the data loss, which is inherent to the soft PDE constraints in many existing methods.
2025
- Individualizing Glioma Radiotherapy Planning by Optimization of a Data and Physics-Informed Discrete LossMichal Balcerak, Jonas Weidner, Petr Karnakov, Ivan Ezhov, Sergey Litvinov, Petros Koumoutsakos, Tamaz Amiranashvili, Ray Zirui Zhang, John S. Lowengrub, Igor Yakushev, Benedikt Wiestler, and Bjoern MenzeNature Communications, Jul 2025
Brain tumor growth is unique to each glioma patient and extends beyond what is visible in imaging scans, infiltrating surrounding brain tissue. Understanding these hidden patient-specific progressions is essential for effective therapies. Current treatment plans for brain tumors, such as radiotherapy, typically involve delineating a uniform margin around the visible tumor on pre-treatment scans to target this invisible tumor growth. This “one size fits all" approach is derived from population studies and often fails to account for the nuances of individual patient conditions. We present the Glioma Optimizing the Discrete Loss (GliODIL) framework, which infers the full spatial distribution of tumor cell concentration from available multi-modal imaging, leveraging a Fisher-Kolmogorov type physics model to describe tumor growth. This is achieved through the newly introduced method of Optimizing the Discrete Loss (ODIL), where both data and physics-based constraints are softly assimilated into the solution. Our test dataset comprises 152 glioblastoma patients with pre-treatment imaging and post-treatment follow-ups for tumor recurrence monitoring. By blending data-driven techniques with physics-based constraints, GliODIL enhances recurrence prediction in radiotherapy planning, challenging traditional uniform margins and strict adherence to the Fisher-Kolmogorov partial differential equation model, which is adapted for complex cases. Current treatment plans for brain tumors fail to account for the nuances of diversified individual patient conditions regarding the tumor infiltrating margins in brain tissues. Here this group reports a Fisher-Kolmogorov type physics model inferring the full spatial distribution of tumor cell concentration from available multi-modal imaging on glioblastoma patients for tumor recurrence monitoring.
- Personalized Predictions of Glioblastoma Infiltration: Mathematical Models, Physics-Informed Neural Networks and Multimodal ScansRay Zirui Zhang, Ivan Ezhov, Michal Balcerak, Andy Zhu, Benedikt Wiestler, Bjoern Menze, and John S. LowengrubMedical Image Analysis, Apr 2025
Predicting the infiltration of Glioblastoma (GBM) from medical MRI scans is crucial for understanding tumor growth dynamics and designing personalized radiotherapy treatment plans. Mathematical models of GBM growth can complement the data in the prediction of spatial distributions of tumor cells. However, this requires estimating patient-specific parameters of the model from clinical data, which is a challenging inverse problem due to limited temporal data and the limited time between imaging and diagnosis. This work proposes a method that uses Physics-Informed Neural Networks (PINNs) to estimate patient-specific parameters of a reaction–diffusion partial differential equation (PDE) model of GBM growth from a single 3D structural MRI snapshot. PINNs embed both the data and the PDE into a loss function, thus integrating theory and data. Key innovations include the identification and estimation of characteristic non-dimensional parameters, a pre-training step that utilizes the non-dimensional parameters and a fine-tuning step to determine the patient specific parameters. Additionally, the diffuse-domain method is employed to handle the complex brain geometry within the PINN framework. The method is validated on both synthetic and patient datasets, showing promise for personalized GBM treatment through parametric inference within clinically relevant timeframes.
2024
- A Compact Coupling Interface Method with Second-Order Gradient Approximation for Elliptic Interface ProblemsRay Zirui Zhang and Li-Tien ChengJournal of Scientific Computing, Jun 2024
We propose the Compact Coupling Interface Method, a finite difference method capable of obtaining second-order accurate approximations of not only solution values but their gradients, for elliptic complex interface problems with interfacial jump conditions. Such elliptic interface boundary value problems with interfacial jump conditions are a critical part of numerous applications in fields such as heat conduction, fluid flow, materials science, and protein docking, to name a few. A typical example involves the construction of biomolecular shapes, where such elliptic interface problems are in the form of linearized Poisson–Boltzmann equations, involving discontinuous dielectric constants across the interface, that govern electrostatic contributions. Additionally, when interface dynamics are involved, the normal velocity of the interface might be comprised of the normal derivatives of solution, which can be approximated to second-order by our method, resulting in accurate interface dynamics. Our method, which can be formulated in arbitrary spatial dimensions, combines elements of the highly-regarded Coupling Interface Method, for such elliptic interface problems, and Smereka’s second-order accurate discrete delta function. The result is a variation and hybrid with a more compact stencil than that found in the Coupling Interface Method, and with advantages, borne out in numerical experiments involving both geometric model problems and complex biomolecular surfaces, in more robust error profiles.
2023
- Explicit-Solute Implicit-Solvent Molecular Simulation with Binary Level-Set, Adaptive-Mobility, and GPUShuang Liu, Zirui Zhang, Hsiao-Bing Cheng, Li-Tien Cheng, and Bo LiJournal of Computational Physics, Jan 2023
Coarse-grained modeling and efficient computer simulations are critical to the study of complex molecular processes with many degrees of freedom and multiple spatiotemporal scales. Variational implicit-solvent model (VISM) for biomolecular solvation is such a modeling framework, and its initial success has been demonstrated consistently. In VISM, an effective free-energy functional of solute-solvent interfaces is minimized, and the surface energy is a key component of the free energy. In this work, we extend VISM to include the solute mechanical interactions, and develop fast algorithms and GPU implementation for the extended variational explicit-solute implicit-solvent (VESIS) molecular simulations to determine the underlying molecular equilibrium conformations. We employ a fast binary level-set method for minimizing the solvation free energy of solute-solvent interfaces and construct an adaptive–mobility gradient descent method for solute atomic optimization. We also implement our methods on the integrated GPU. Numerical tests and applications to several molecular systems verify the accuracy, stability, and efficiency of our methods and algorithms. It is found that our new methods and GPU implementation improve the efficiency of the molecular simulation significantly over the CPU implementation. Our fast computational techniques may enable us to simulate very large systems such as protein-protein interactions and membrane dynamics for which explicit-solvent all-atom molecular dynamics simulations can be very expensive.
2021
- Coupling Monte Carlo, Variational Implicit Solvation, and Binary Level-Set for Simulations of Biomolecular BindingZirui Zhang, Clarisse G. Ricci, Chao Fan, Li-Tien Cheng, Bo Li, and J. Andrew McCammonJournal of Chemical Theory and Computation, Mar 2021
We develop a hybrid approach that combines the Monte Carlo (MC) method, a variational implicit-solvent model (VISM), and a binary level-set method for the simulation of biomolecular binding in an aqueous solvent. The solvation free energy for the biomolecular complex is estimated by minimizing the VISM free-energy functional of all possible solute–solvent interfaces that are used as dielectric boundaries. This functional consists of the solute volumetric, solute–solvent interfacial, solute–solvent van der Waals interaction, and electrostatic free energy. A technique of shifting the dielectric boundary is used to accurately predict the electrostatic part of the solvation free energy. Minimizing such a functional in each MC move is made possible by our new and fast binary level-set method. This method is based on the approximation of surface area by the convolution of an indicator function with a compactly supported kernel and is implemented by simple flips of numerical grid cells locally around the solute–solvent interface. We apply our approach to the p53-MDM2 system for which the two molecules are approximated by rigid bodies. Our efficient approach captures some of the poses before the final bound state. All-atom molecular dynamics simulations with most of such poses quickly reach the final bound state. Our work is a new step toward realistic simulations of biomolecular interactions. With further improvement of coarse graining and MC sampling, and combined with other models, our hybrid approach can be used to study the free-energy landscape and kinetic pathways of ligand binding to proteins.
2018
- Runout Scaling and Deposit Morphology of Rapid MudflowsL. Jing, C. Y. Kwok, Y. F. Leung, Z. Zhang, and L. DaiJournal of Geophysical Research: Earth Surface, Mar 2018
Prediction of runout distance and deposit morphology is of great importance in hazard mitigation of geophysical flows, including viscoplastic mudflows. The major rheological parameters of mudflows, namely, yield stress and viscosity, are crucial factors in controlling the runout and deposition processes. However, the roles of the two parameters, especially in mudflows with high inertia, remain poorly understood and are not accounted for in runout scaling relations with source volume. Here we investigate the effects of flow rheology on runout scaling and deposit morphology using small-scale laboratory experiments and three-dimensional numerical simulations. We find that yield stress and viscosity both influence flow velocity gained during downslope propagation of mudflows, which is strongly correlated with the runout distance; the role of yield stress is more significant than viscosity. High yield stress and low viscosity lead to an elongated deposit, where longitudinal propagation is more significant than lateral spreading. In contrast, high viscosity promotes the dominance of lateral spreading of the deposit, while low yield stress and moderate viscosity produce an initial elongate deposit, followed by a secondary surge that spreads laterally near the head of the deposit. Following appropriate scaling relations for viscosity and yield stress, a general scaling function is proposed to incorporate flow properties in the well-known correlation of runout distance and source volume. Our findings regarding the inertia effects and the roles of yield stress and viscosity enhance our understanding of mudflows, muddy debris flows, and other viscoplastic geophysical flows.