Skip to content

Physics-Grounded Disentangled Flow Modeling for Brain Disease Progression Trajectory

Conference: ECCV 2026
Paper: ECCV Official
Code: https://github.com/jhuldr/PDF
Area: Medical Imaging
Keywords: Brain Disease Progression Modeling, Disentangled Flow Matching, Physics-Informed Learning, Reaction-Diffusion-Advection PDE, Longitudinal MRI

TL;DR

Addressing the physical implausibility caused by entangling structural deformation with image intensity changes in longitudinal brain disease forecasting, this paper proposes PDF, a disentangled flow matching framework regularized by Fisher-KPP PDE subspace projection, achieving state-of-the-art lesion boundary accuracy and trajectory consistency across three glioma and multiple sclerosis datasets.

Background & Motivation

Longitudinal magnetic resonance imaging (MRI) is essential in the clinical management of progressive neurological conditions, such as gliomas and multiple sclerosis (MS). It enables clinicians to track dynamic lesion evolution over time, providing vital information for early disease intervention, radiotherapy and surgical planning, and quantitative treatment response evaluation. However, modeling longitudinal medical image trajectories poses substantial real-world challenges: clinical follow-up scans suffer from irregular and variable time intervals, inter-scan protocol and scanner discrepancies, and a scarcity of densely paired patient follow-up cohorts.

Recent continuous-time generative paradigms—including Neural Ordinary Differential Equations (Neural ODEs), Schrödinger Bridge diffusion models (such as I2SB), and Trajectory Flow Matching (TFM)—have made remarkable strides in continuous trajectory fitting and image-level forecasting. Nevertheless, existing methods uniformly adopt a direct image-to-image translation formulation. This paradigm implicitly entangles complex structural deformation and subtle intensity variations into a single black-box mapping without modeling the biophysical mechanisms that drive lesion progression. Consequently, when encountering rapid tumor infiltration or dramatic ventricle distortion, these models frequently yield physically inconsistent predictions, blurred anatomical boundaries, or unnatural deformation artifacts.

The core tension lies in the fundamental nature of disease progression: anatomical progression is physically governed by cellular migration and biomechanical tissue deformation, whereas appearance changes are driven by tissue metabolic concentration and microstructural signal alterations. Learning both processes simultaneously without physics-grounded inductive biases forces models to compromise anatomical fidelity. The paper's angle of attack is to explicitly decouple the overall progression into distinct morphology evolution and intensity evolution processes, and introduce classical biophysical tumor dynamics partial differential equations (PDEs) as regularizers. Core idea: disentangle longitudinal brain disease progression into morphology flow matching constrained by Fisher-KPP reaction-diffusion-advection PDE subspace projection and intensity flow matching initialized by anatomical deformation with high-frequency texture preservation, ensuring both realistic tissue appearance and strict biophysical plausibility.

Method

Overall Architecture

The proposed PDF framework decomposes the extrapolation from a baseline image \(I_{t_1}\) and lesion mask \(M_{t_1}\) at time \(t_1\) to a future time \(t_2\) into a dual-stage flow matching pipeline. First, the morphology evolution module predicts the instantaneous velocity field of the lesion concentration and structural deformation under the guidance of a physics-grounded subspace projection regularization derived from a Fisher-KPP PDE. Second, a dense displacement field is estimated between the predicted and baseline lesion masks via variational optical flow to perform semi-Lagrangian spatial warping on the baseline anatomy. Finally, the intensity evolution module refines signal intensity changes on the warped structural scaffold while preserving intrinsic high-frequency lesion textures, generating the final predicted scan \(\hat{I}_{t_2}\).

%%{init: {'flowchart': {'rankSpacing': 24, 'nodeSpacing': 28, 'padding': 6, 'wrappingWidth': 400}}}%%
flowchart TD
    A["Baseline Input: Image $I_{t_1}$ and Mask $M_{t_1}$<br/>Time Interval $\Delta t = t_2 - t_1$"] --> B["Morphology Flow Matching<br/>SDF and Cell Concentration Modeling"]
    B --> C["Physics Span and PDE Regularization<br/>Fisher-KPP Reaction-Diffusion-Advection Projection"]
    C --> D["Dense Deformation Field Optical Flow<br/>Semi-Lagrangian Spatial Warping"]
    D --> E["Intensity Flow Matching and Texture Preservation<br/>Appearance Flow Matching with High-Frequency Guidance"]
    E --> F["Predicted Output: Future Image $\hat{I}_{t_2}$ and Mask $\hat{M}_{t_2}$"]

Key Designs

1. Morphology Flow Matching: Continuous Lesion Cell Concentration Modeling and Velocity Prediction

To circumvent the non-differentiable and discrete nature of binary masks, this module first converts the discrete baseline lesion mask \(M_{t_1}\) into a smooth signed distance function (SDF) \(d(\mathbf{x})\), where negative values indicate regions inside the lesion and positive values denote regions outside. This field is subsequently mapped into a normalized continuous lesion cell concentration field \(c(\mathbf{x}) \in [0, 1]\) via a Sigmoid transformation parameterized by scale factor \(\tau\): $\(c(\mathbf{x}) = \frac{1}{1 + \exp\left(\frac{d(\mathbf{x})}{\tau}\right)}\)$ A dedicated U-Net architecture \(f_\theta\) takes the baseline mask \(M_{t_1}\), structural context image \(I_{t_1}\), and temporal interval \(\Delta t = t_2 - t_1\) as conditioning inputs to predict the instantaneous morphological velocity field \(v_\theta = f_\theta(M_{t_1}; I_{t_1}, t_2 - t_1)\). In flow matching, the empirical target velocity is defined by \(v_{\text{morph}} = \frac{M_{t_2} - M_{t_1}}{t_2 - t_1}\), enabling the network to learn smooth vector fields that handle irregularly sampled longitudinal observation times naturally.

2. Reaction-Diffusion-Advection Physics Span and PDE Regularization: Biophysical Growth Constraints

To prevent purely data-driven networks from generating biologically implausible deformations, the authors incorporate classical Fisher-KPP reaction-diffusion dynamics augmented by an advection term: $\(\frac{\partial c}{\partial t} = D \Delta c + \rho c(1 - c) - \mathbf{u} \cdot \nabla c\)$ Here, the diffusion term \(D \Delta c\) governs isotropic random cell migration and infiltrative wave expansion, the reaction term \(\rho c(1 - c)\) represents local cellular proliferation following logistic growth kinetics, and the advection term \(-\mathbf{u} \cdot \nabla c\) characterizes bulk tissue displacement caused by mass effect. To avoid solving computationally intractable patient-specific inverse PDE problems at every iteration, a two-step pseudo-inverse estimation is employed: optimal \((D, \rho)\) are first fitted via least-squares between consecutive concentration maps to generate an intermediate diffusion-reaction prediction \(\tilde{c}_{t_2}\), and the residual motion between \(\tilde{c}_{t_2}\) and observed \(c_{t_2}\) is estimated using a TV-L1 variational optical flow formulation to yield the instantaneous advection field \(\mathbf{u}(\mathbf{x}) = -\frac{\mathbf{d}(\mathbf{x})}{t_2 - t_1}\).

These components define a dynamic physics span: $\(\mathcal{S}(c_{t_1}; \mathbf{u}) = \text{span}\left\{ \Delta c_{t_1},\; c_{t_1}(1 - c_{t_1}),\; -\mathbf{u} \cdot \nabla c_{t_1} \right\}\)$ During training, the predicted velocity \(v_\theta\) is projected onto this subspace by finding optimal linear combination weights \((\alpha, \beta, \gamma) = \arg\min_{\alpha, \beta, \gamma} \|v_\theta - (\alpha \Delta c_{t_1} + \beta c_{t_1}(1 - c_{t_1}) + \gamma(-\mathbf{u} \cdot \nabla c_{t_1}))\|^2\). The PDE regularization loss \(\mathcal{L}_{\text{PDE}}\) is defined as the squared residual of this projection: $\(\mathcal{L}_{\text{PDE}} = \|v_\theta - \Pi_{\mathcal{S}}(v_\theta)\|^2\)$ This constraint penalizes any velocity components that stray from the biological diffusion-reaction-advection subspace without requiring forward PDE numerical solvers during training.

3. Intensity Flow Matching and Texture Preservation: Appearance Refinement Conditioned on Deformed Anatomy

Integrating the predicted morphological velocity \(v_\theta\) over \((t_2 - t_1)\) yields the evolved lesion mask \(\hat{M}_{t_2}\). A TV-L1 optical flow algorithm is applied to compute a dense deformation field \(\mathbf{f}(\mathbf{x})\) between \(M_{t_1}\) and \(\hat{M}_{t_2}\). The baseline scan \(I_{t_1}\) is then warped via a semi-Lagrangian sampling formulation to produce \(\hat{I}^{\text{morph}}(\mathbf{x}) = I_{t_1}(\mathbf{x} + \mathbf{f}(\mathbf{x}))\).

While \(\hat{I}^{\text{morph}}\) accurately captures spatial mass shifts and ventricle displacement, it does not account for internal signal changes (e.g., central tumor necrosis, contrast enhancement variability, or peritumoral edema fluctuations). A dedicated intensity flow matching network \(g_\phi(z_\tau; M)\) takes \(\hat{I}^{\text{morph}}\) as a structural baseline and evolves the image along a linear interpolation path. To prevent generative smoothing from wiping out fine lesion heterogeneity, high-frequency textural components extracted from the baseline image \(I_{t_1}\) are incorporated as conditioning features during flow integration, producing the final appearance-consistent volume \(\hat{I}_{t_2}\).

A Worked Example

Consider an adult patient with glioblastoma undergoing post-operative surveillance. At baseline scan (\(t_1\)), a residual enhancing lesion is identified with a volume of 12.5 mL. The clinician seeks to forecast lesion status 180 days later (\(\Delta t = 0.5\) years): 1. The discrete lesion mask is converted to an SDF and mapped to a continuous concentration field \(c_{t_1}\). The morphology flow network \(f_\theta\) takes \((M_{t_1}, I_{t_1}, \Delta t)\) and predicts an unconstrained velocity field. 2. The physics span computes the local Laplacian \(\Delta c_{t_1}\) (diffusion), logistic term \(c_{t_1}(1-c_{t_1})\) (proliferation), and advection gradient \(-\mathbf{u}\cdot\nabla c_{t_1}\). Orthogonal projection onto this span filters out non-physical velocity noise, reducing the projection residual ratio from 0.20 to 0.12. 3. Integrating the regularized velocity forward produces predicted mask \(\hat{M}_{t_2}\) with an expanded volume of 21.8 mL. TV-L1 flow estimates displacement field \(\mathbf{f}(\mathbf{x})\) from the baseline to the predicted mask, and semi-Lagrangian warping realistically deforms the ipsilateral lateral ventricle. 4. The intensity network \(g_\phi\) uses the warped image as an initialization scaffold, ingests high-frequency lesion texture descriptors, and synthesizes subtle internal necrotic signal hypointensity and enhancing rim margins, delivering a calibrated T1-weighted follow-up prediction.

Loss & Training

The morphology evolution network is trained using a weighted objective combining flow matching and physics regularization: $\(\mathcal{L}_{\text{morph}} = \mathcal{L}_{\text{FM}}(v_\theta, v_{\text{morph}}) + \lambda \mathcal{L}_{\text{PDE}}\)$ where \(\mathcal{L}_{\text{FM}} = \|v_\theta - v_{\text{morph}}\|^2\), with \(\lambda = 0.01\) identified through ablation as the optimal balance between empirical data fit and physical stiffness. The intensity flow matching module is optimized via standard image-space flow matching \(\mathcal{L}_{\text{intensity}} = \|v_\phi - v_{\text{intensity}}\|^2\). Both networks use a U-Net backbone trained end-to-end on an NVIDIA RTX PRO 6000 GPU (PyTorch 2.9.1, CUDA 13.0). All datasets are strictly partitioned into 6:2:2 train/validation/test splits at the subject series level to ensure zero temporal cross-contamination.

Key Experimental Results

Main Results

Quantitative evaluations are performed across three longitudinal brain MRI cohorts: UCSF (302 adult glioma patients, T1-weighted), LUMIERE (91 glioblastoma patients, 795 post-operative longitudinal series, contrast-enhanced T1), and LMSLS (79 multiple sclerosis longitudinal cases, FLAIR). Evaluated baselines include Time-conditional UNet (T-UNet), Image-to-Image Schrödinger Bridge (I2SB), Trajectory Flow Matching (TFM), and ImageFlowNet.

Dataset Metric Ours (PDF) T-UNet [34] I2SB [23] TFM [50] ImageFlowNet [22]
UCSF (Glioma) DSC ↑ 0.752 ± 0.003 0.728 ± 0.004 0.725 ± 0.005 0.696 ± 0.006 0.724 ± 0.004
HD (mm) ↓ 25.952 ± 0.654 27.322 ± 0.812 28.839 ± 0.941 31.701 ± 1.284 28.402 ± 0.893
PSNR (dB) ↑ 31.364 ± 0.092 31.739 ± 0.103 31.393 ± 0.097 30.659 ± 0.125 31.675 ± 0.108
SSIM ↑ 0.920 ± 0.001 0.920 ± 0.002 0.920 ± 0.002 0.898 ± 0.003 0.909 ± 0.003
L1 Error ↓ 0.072 ± 0.002 0.085 ± 0.002 0.093 ± 0.002 0.091 ± 0.002 0.085 ± 0.002
LUMIERE (Glioblastoma) DSC ↑ 0.534 ± 0.013 0.520 ± 0.020 0.522 ± 0.035 0.492 ± 0.015 0.493 ± 0.033
HD (mm) ↓ 33.672 ± 1.229 38.135 ± 1.925 36.105 ± 1.720 39.211 ± 1.530 41.994 ± 1.653
PSNR (dB) ↑ 37.216 ± 0.349 36.753 ± 0.185 37.206 ± 0.348 36.317 ± 0.406 37.178 ± 0.346
SSIM ↑ 0.957 ± 0.002 0.845 ± 0.004 0.953 ± 0.003 0.894 ± 0.001 0.916 ± 0.003
L1 Error ↓ 0.071 ± 0.002 0.097 ± 0.003 0.084 ± 0.003 0.088 ± 0.002 0.078 ± 0.003
LMSLS (Multiple Sclerosis) DSC ↑ 0.696 ± 0.004 0.673 ± 0.005 0.682 ± 0.006 0.670 ± 0.007 0.687 ± 0.005
HD (mm) ↓ 20.867 ± 0.623 26.354 ± 0.932 23.617 ± 0.784 24.438 ± 0.851 24.296 ± 0.741
PSNR (dB) ↑ 29.314 ± 0.087 29.248 ± 0.118 29.312 ± 0.094 29.441 ± 0.102 28.972 ± 0.113
SSIM ↑ 0.900 ± 0.002 0.805 ± 0.003 0.893 ± 0.002 0.889 ± 0.003 0.890 ± 0.003
L1 Error ↓ 0.058 ± 0.002 0.083 ± 0.002 0.068 ± 0.003 0.062 ± 0.004 0.072 ± 0.002

Ablation Study

1. Complementary Value of Disentangled Morphology and Intensity Branches

Branch Configuration UCSF DSC ↑ UCSF HD ↓ LUMIERE DSC ↑ LUMIERE HD ↓ LMSLS DSC ↑ LMSLS HD ↓ Note
Morphology Branch Only 0.750 26.038 0.532 34.066 0.694 21.855 Excels at anatomical warping and boundary localization
Intensity Branch Only 0.749 25.806 0.529 35.143 0.690 21.749 Maintains continuous image appearance; weaker at structural deformation
Joint Full Model 0.752 25.952 0.534 33.672 0.696 20.867 Synergistic integration yields the most accurate boundary and appearance

2. Impact of PDE Regularization Loss \(\mathcal{L}_{\text{PDE}}\)

Loss Setting UCSF DSC ↑ UCSF HD ↓ LUMIERE DSC ↑ LUMIERE HD ↓ LMSLS DSC ↑ LMSLS HD ↓ Note
w/o \(\mathcal{L}_{\text{PDE}}\) 0.745 26.216 0.524 34.211 0.690 22.337 Unconstrained velocity field; residual energy ratio is 0.20
w/ \(\mathcal{L}_{\text{PDE}}\) 0.752 25.952 0.534 33.672 0.696 20.867 Subspace projection reduces non-physical error; residual drops to 0.12

3. Sensitivity to Regularization Weight \(\lambda\) - With \(\lambda = 1.0\): UCSF DSC is 0.745, LUMIERE DSC is 0.526, LMSLS DSC is 0.691 (excessive physical rigidity suppresses the neural network's empirical flexibility). - With \(\lambda = 0.1\): UCSF DSC is 0.746, LUMIERE DSC is 0.529, LMSLS DSC is 0.693. - With \(\lambda = 0.01\) (default): UCSF DSC reaches 0.752, LUMIERE DSC is 0.534, and LMSLS DSC is 0.696, establishing an optimal balance between empirical data fit and physical regularization.

Key Findings

  • Lesion-specific metrics exhibit superior clinical sensitivity: On whole-brain metrics (PSNR, SSIM), baseline models yield scores close to PDF because normal brain parenchyma dominates the total volume. In contrast, on lesion-specific metrics (DSC, HD), PDF delivers striking improvements—reducing Hausdorff Distance on LUMIERE from 38.135 mm to 33.672 mm.
  • Robustness in aggressive progression edge cases: In high-grade tumor cases where aggressive lesion expansion compresses and distorts the ipsilateral ventricle, all baseline methods struggle (DSC ranging between 0.63 and 0.77). In contrast, PDF accurately forecasts the mass effect and ventricle shift, achieving a DSC of 0.85 and reducing HD to 16.4 mm.

Highlights & Insights

  • Reformulating biophysical dynamics as lightweight subspace projection: Traditional physics-informed modeling depends on costly iterative parameter inversion. PDF recasts the three components of the Fisher-KPP equation (diffusion, reaction, advection) as an empirical span and constrains predicted velocities via ordinary least squares projection, combining generative flexibility with biophysical validity at negligible computational cost.
  • Disentangled morphology and intensity flow matching: Decoupling macroscopic tissue deformation from microscopic signal intensity variations cleanly resolves the blurriness and distortion common to single-network end-to-end models.
  • Broad cross-pathology and cross-modality generalization: The model demonstrates robust generalizability across diverse clinical presentations—ranging from mass-expanding gliomas on T1 and contrast-enhanced T1 scans to subtle, non-mass-effect multiple sclerosis plaque evolution on FLAIR MRI (where HD drops from 26.35 mm to 20.87 mm).

Limitations & Future Work

  • Dependence on baseline segmentation accuracy: The physics span relies on initial lesion segmentations to fit baseline diffusion \(D\), reaction rate \(\rho\), and advection fields \(\mathbf{d}\). Severe segmentation errors or missing focal lesions at baseline may degrade the quality of the physics span.
  • Isotropic diffusion assumption vs. white matter tract anisotropy: The model currently assumes isotropic diffusion across the entire parenchymal volume, overlooking the directed axonal migration pathways along white matter tracts. Integrating diffusion tensor imaging (DTI) priors could further improve directional infiltration accuracy.
  • Absence of clinical treatment intervention modeling: Patients frequently undergo chemotherapy, surgical resection, or radiotherapy between imaging visits. Explicitly conditioning the reaction coefficient \(\rho\) and flow fields on treatment variables represents a compelling direction for counterfactual trajectory simulation.
  • vs ImageFlowNet [22]: ImageFlowNet trains Neural ODEs directly in multi-scale latent image feature space. While continuous in time, it blends structure and intensity into an entangled latent representation without mechanistic physical priors. PDF decouples morphology and intensity and enforces PDE projection, outperforming ImageFlowNet in boundary overlap (DSC).
  • vs I2SB (Image-to-Image Schrödinger Bridge) [23]: I2SB models transitions between endpoint distributions via diffusion bridges. Although it produces sharp images, it tends to misinterpret spatial boundary shifts as smooth gray-level fades, introducing diffusive artifacts. PDF preserves structural integrity through explicit optical flow warping.
  • vs Classical Physics Inversion Methods (e.g., Learn-Morph-Infer [9], GliODIL [3]): These techniques optimize cell density distributions via iterative PDE solvers. In contrast, PDF embeds PDE dynamics as a geometric regularizer within a forward flow-matching architecture, enabling real-time generative forecasting during inference.

Rating

  • Novelty: ⭐⭐⭐⭐⭐ [Pioneers the use of Fisher-KPP reaction-diffusion-advection PDE subspace projection to regularize generative flow matching for longitudinal medical imaging]
  • Experimental Thoroughness: ⭐⭐⭐⭐⭐ [Extensively validated across three diverse longitudinal brain MRI cohorts with comprehensive lesion-level and whole-brain metric evaluations and detailed ablations]
  • Writing Quality: ⭐⭐⭐⭐⭐ [Rigorous mathematical formulation, clear narrative architecture, and compelling comparative figures]
  • Value: ⭐⭐⭐⭐⭐ [Provides an impactful paradigm for combining biophysical domain knowledge with modern continuous generative models in computational neuro-oncology]