PaperlayerView exact revision on arXiv

Source-linked AI summary

A hybrid CNN-adjoint optimization framework for reconstruction of viscoelastic tissue properties in magnetic resonance elastography

Anwesa Dey, Johann Rudi, Elena Cherkaev

arXiv:2610.02634v1math.OCcs.LGcs.NAmath.NA

TL;DR

Recovering complex shear modulus from displacement data is ill-posed, particularly with noisy or incomplete measurements. The paper combines a CNN initialization with adjoint-based PDE-constrained optimization and reports useful reconstructions on unseen configurations, with optimization improving accuracy in noiseless experiments.

  • Problem

    Recovering complex shear modulus from displacement measurements is ill-posed in the presence of noise, model uncertainty, and partial or corrupted data.

    Source: §1 · p. 2

  • Method

    The framework uses a CNN to map displacement data to a coefficient-field initialization for adjoint-based PDE-constrained optimization.

    Source: §1 · p. 3, §4.3 · p. 15

  • Results

    The CNN provides useful initial reconstructions across inclusion geometries and generalizes to configurations absent from training, while optimization improves accuracy in the noiseless setting.

    Source: §5 · p. 36

  • Takeaways & Limitations

    In this two-dimensional MRE model, learned reconstructions can initialize physics-based refinement rather than serve as the final reconstruction alone.

    Source: §4.3 · p. 15, §5 · p. 36

    In noisy experiments, optimization produced oscillatory artifacts away from inclusions, making regularization strength relevant.

    Source: §4.5.1 · p. 28

Abstract

from arXiv · show

Magnetic resonance elastography (MRE) is a noninvasive imaging modality for quantifying the viscoelastic properties of soft tissues from shear wave propagation. Recovering the complex-valued shear modulus from measured displacement fields leads to a severely ill-posed inverse problem, particularly in the presence of noise and limited boundary excitations. We investigate adjoint-based optimization, convolutional neural network (CNN) reconstruction, and a hybrid framework combining both approaches. The forward model is based on a scalar form of the modified stationary Stokes system with a complex shear modulus. We establish well-posedness of the forward problem, existence of minimizers, and first-order optimality conditions for the adjoint-based formulation, and implement a nonlinear conjugate-gradient method with Armijo line search. While PDE-constrained optimization can accurately refine coefficient reconstructions, its performance depends strongly on initialization. We therefore construct a two-dimensional CNN that maps complex-valued displacement measurements to spatially varying complex shear modulus fields and provides rapid, informative initial reconstructions. The proposed hybrid method uses the CNN reconstruction to initialize the adjoint-based optimization, yielding faster convergence and improved accuracy. The CNN is trained on coefficient fields containing individual perturbations and tested on both individual and previously unseen combined configurations. Numerical experiments demonstrate that the CNN generalizes to these more challenging configurations, while subsequent PDE-constrained optimization further refines the reconstructed coefficient. These results demonstrate the potential of combining data-driven initialization with physics-based optimization for efficient and accurate MRE reconstruction.

1 Introduction

MRE coefficient recovery is an ill-posed inversion, motivating a combination of learned initial guesses and physics-based optimization. The paper analyzes adjoint-based reconstruction theoretically and evaluates CNN and hybrid methods, including on unseen combined configurations and noisy data.

1 Introduction

  • MRE recovers tissue viscoelastic properties from shear-wave displacements, but noise, model uncertainty, and partial or corrupted data make the inverse problem ill-posed.MRE adds mechanical information about tissue behavior to imaging.Source: §1 · p. 2
  • Prior studies considered scalar and time-dependent Kelvin–Voigt models and modified stationary Stokes models, alongside finite-element and Levenberg–Marquardt inversion schemes.Other cited work addressed uniqueness, stability, and Landweber reconstruction [6].Source: §1 · p. 2
  • Adjoint-based methods efficiently compute objective derivatives with respect to spatially distributed unknown coefficients, while learning methods approximate inverse maps and have been applied to wave-based inverse problems.Source: §1 · pp. 2–3
  • The hybrid strategy is motivated by the complementary roles of physical models and data-driven approximations, with learned reconstructions serving to initialize or accelerate model-based optimization [17].The paper applies this perspective to MRE reconstruction.Source: §1 · p. 3
  • The paper analyzes adjoint-based reconstruction and evaluates CNN and hybrid methods, with the CNN providing an informed initial guess for optimization.The CNN is trained on individual perturbations and tested on individual and combined configurations, including combinations absent from training.Source: §1 · p. 3

2 Modified stationary Stokes system

The paper formulates MRE inversion using a modified stationary Stokes model and states the task of recovering the unknown complex shear modulus from displacement data. Under incompressibility, the model reduces to a scalar equation.

2 Modified stationary Stokes system

  • The spaces C0,1(¯Ω) and Hs(Ω) denote Lipschitz-continuous functions and L2-based Sobolev spaces, respectively.The model discussion also groups multiple solutions as u = (u1, . . . , uI).Source: §2 · pp. 3–4
  • The modified stationary Stokes system is proposed as a PDE model for MRE, with the inverse problem using displacement data to recover the unknown shear modulus A.The model description identifies Ω as the region of interest, u as wave displacement, p as pressure, and ω as angular frequency.Source: §2 · p. 4

2.1 Scalar model for the modified stationary Stokes system

  • Assuming incompressibility, ∇·u = 0, the modified stationary Stokes system reduces to a scalar equation.Source: §2.1 · p. 4
  • The shear modulus A(x) is modeled as a constant baseline plus a perturbation representing deviations associated with liver abnormalities.Source: §2.1 · p. 4

3 Reconstruction of shear modulus function using adjoint-based optimization

The adjoint-based approach reconstructs the complex shear-modulus perturbation from displacement data by minimizing a data-fit objective with Tikhonov regularization. The section establishes theoretical properties for this optimization, derives equivalent optimality systems, and describes numerical optimization and CNN reconstruction concepts.

3 Reconstruction of shear modulus function using adjoint-based optimization

  • The objective fits measured displacement data while Tikhonov terms regularize the real and imaginary shear-modulus perturbations.The unknown complex modulus is represented through its real and imaginary components, B1 and B2.Source: §3 · pp. 4–5

3.1 Theoretical results

  • For given boundary data and an admissible coefficient, the forward problem has a unique solution and satisfies a stability estimate.The stated result applies for f in H1/2(∂Ω) and B in L∞(Ω).Source: §3.1 · p. 5
  • The stated uniqueness estimates require multiple boundary and interior data pairs, although one isotropic coefficient can be recovered from one measurement in any dimension.The passage cautions that the dimension-dependent pair counts are not optimal and cites [20, Section 4].Source: §3.1 · p. 5
  • The objective is sequentially weakly lower semicontinuous and Fréchet differentiable, and a minimizer exists over the admissible set.The existence result provides an admissible coefficient and associated state solutions attaining the minimum.Source: §3.1 · p. 6

3.2 Optimality conditions

  • The first-order optimality systems from the complex formulation and the separated real-imaginary PDE system are equivalent.The real-imaginary system is obtained by expanding the complex-valued PDE and separating its components.Source: §3.2 · pp. 8–9

3.3 Numerical method

  • The numerical method uses nonlinear conjugate gradients with an Armijo line search, requiring gradient evaluations rather than Hessian formation and inversion.The passage reports faster overall convergence than Newton-type schemes in the described settings.Source: §3.3 · p. 9

3.4 Reconstruction of shear modulus function using convolutional neural networks

  • CNNs combine convolutional layers and nonlinear activations to learn complex problem-specific priors from data.The passage contrasts this data-driven regularization with explicit penalties that rely on predefined assumptions such as smoothness or sparsity.Source: §3.4 · p. 10

4 Reconstruction of complex-valued shear modulus function in two dimensions

The two-dimensional reconstruction uses a scalar modified Stokes model with a complex shear modulus, represented as a constant background plus spatial perturbations in its real and imaginary parts.

4 Reconstruction of complex-valued shear modulus function in two dimensions

  • The forward model is ∇·(A(x, y)∇u(x, y)) + ρω²u(x, y) = 0, with complex displacement u and complex shear modulus A.The model is used to reconstruct spatially varying shear modulus from displacement data.Source: §4 · p. 10
  • The coefficient is written as A(x, y) = (1 + i) + B(x, y), where B has real and imaginary perturbation components.B1 and B2 perturb the real and imaginary parts of the homogeneous background coefficient.Source: §4 · p. 11
  • The experiments use Gaussian inclusions and mollified elliptical disc inclusions to represent different spatial perturbation families.The latter adds a coefficient family distinct from Gaussian inclusions.Source: §4 · p. 11

4.1 Coefficient models

The coefficient models use Gaussian perturbations and smoothly mollified elliptical discs, with independently generated real and imaginary components.

Gaussian inclusions

  • Gaussian perturbations are parameterized by a center, spatial spread, and amplitude, with real and imaginary perturbations generated independently.The Gaussian spread is controlled by σ.Source: Gaussian inclusions · p. 11

Mollified elliptical disc-shaped inclusions

  • Mollified elliptical discs use rotated coordinates and two semi-axis lengths to define an elliptical radius.The amplitude sets the perturbation magnitude, while ε controls the width of the smooth boundary transition; real and imaginary parameters are generated independently.Source: Mollified elliptical disc-shaped… · p. 11
  • The mollified discs have relatively flat interiors and localized smooth boundary transitions, unlike the Gaussian profiles in the same training set.This gives the coefficient families qualitatively different spatial structures.Source: Mollified elliptical disc-shaped… · p. 11

4.2 Forward problem and numerical discretization

The forward problem generates complex displacement fields from spatially varying complex shear moduli; the inverse problem seeks the coefficient from displacement measurements despite ill-posedness.

4.2 Forward problem and numerical discretization

  • The computational domain is [0, 2] × [0, 2], discretized on a 48 × 48 grid, with density 1 kg/m3 and angular frequency 1 Hz.The PDE is discretized by central finite differences, with coefficient values arithmetically averaged at cell interfaces, and the resulting sparse system yields the complex displacement field.Source: §4.2 · p. 12
  • Complex-valued Dirichlet data are imposed on all four boundaries to model actuators vibrating tissue edges and generating shear waves.The experiments include Gaussian and mollified elliptical disc coefficient fields and their corresponding forward-PDE displacement fields.Source: §4.2 · p. 12, Figure 1 · p. 12, Figure 2 · p. 13

4.2.1 The inverse problem

  • The inverse problem reconstructs the spatially varying complex shear modulus from internal displacement measurements, but small data perturbations can cause significant coefficient errors.Accurate and stable recovery therefore requires a reconstruction strategy and regularization.Source: §4.2.1 · pp. 13–14
  • The hybrid strategy uses a CNN to produce an approximate coefficient field, then refines it with PDE-constrained optimization and regularization.The CNN provides an informative initialization because optimization performance can depend strongly on the initial guess.Source: §4.2.1 · p. 14

4.3 Hybrid reconstruction framework based on convolutional neural network and optimization

The framework uses a two-dimensional CNN to map complex displacement data to an initial coefficient field, then applies PDE-constrained optimization to refine that estimate.

4.3 Hybrid reconstruction framework based on convolutional neural network and optimization

  • The hybrid method uses the CNN prediction as the optimization initial guess, so the network need only provide a sufficiently close approximation for PDE-based refinement.This reduces reliance on arbitrary initial guesses while retaining the governing physical model in the final reconstruction.Source: §4.3 · pp. 14–15

4.3.1 CNN-based initialization

  • Unlike a parameter-based model, the two-dimensional CNN predicts the complete spatial coefficient field rather than inclusion parameters such as center, amplitude, spread, or orientation.Real and imaginary displacement components form two input channels, and the network outputs separate real and imaginary coefficient fields.Source: §4.3.1 · p. 15, §4.3.3 · p. 16

4.3.2 Training dataset

  • The 4,000-pair training set contains equal numbers of Gaussian and mollified elliptical disc perturbations, combining smoothly varying profiles with flatter interiors and localized boundaries.The two families were selected to expose the CNN to qualitatively different spatial features.Source: §4.3.2 · p. 15

4.3.3 CNN architecture and training

  • The encoder–decoder reduces spatial resolution from 48×48 to 12×12 and restores it, with separate output branches for the real and imaginary coefficients.The model has 444,226 trainable parameters and is trained for 300 epochs using Adam and mean-squared-error loss; training uses an 80/20 split of 3,200 and 800 samples.Source: §4.3.3 · pp. 16–17, Table 1 · p. 17

4.4 Hybrid CNN–optimization reconstruction

The hybrid method uses the CNN prediction to initialize PDE-constrained optimization, which refines the coefficient against measured displacement data with regularization. A discrete adjoint supplies gradients for nonlinear conjugate-gradient updates, with Armijo backtracking selecting step lengths.

4.4 Hybrid CNN–optimization reconstruction

  • The CNN prediction initializes minimization of displacement misfit with a zero-order Tikhonov penalty controlling perturbation magnitude.The regularization parameter β balances data fitting against coefficient magnitude; the CNN prediction is used only as the starting point, not in the penalty.Source: §4.4 · pp. 17–18

4.4.1 Adjoint-based gradient computation

  • Each iteration solves the forward problem and an adjoint problem driven by the displacement residual, then computes gradients for both coefficient components.The real and imaginary parts are optimized simultaneously to update the complete complex-valued coefficient.Source: §4.4.1 · p. 18

4.4.2 Nonlinear conjugate-gradient iteration

  • Nonlinear conjugate gradient starts with the negative gradient and then combines the current negative gradient with the previous search direction.The coefficient is updated along the resulting direction using a step length.Source: §4.4.2 · p. 18

4.4.3 Armijo line search

  • Armijo backtracking reduces trial step lengths until sufficient decrease holds, limiting oversized updates while permitting larger corrections when the objective drops enough.Optimization stops when its prescribed criterion is met or the maximum iteration count is reached.Source: §4.4.3 · p. 19

4.5 Numerical experiments

The numerical experiments evaluate CNN reconstruction on individual inclusions and test hybrid refinement, robustness to noisy displacement data, and generalization to unseen combined structures.

4.5.1 Experiment I: reconstruction of individual shape families

  • Across 500 individual-shape test samples, mean L2 errors were 6.4693×10−2 for the real component, 6.6093×10−2 for the imaginary component, and 9.4681 × 10−2 for the full complex coefficient.The corresponding medians were slightly below the means.Source: Table 2 · p. 20, §4.5.1 · p. 20
  • Error histograms across the 500 individual-shape test samples characterize variability in CNN accuracy beyond a single representative reconstruction.The test set contains both Gaussian and mollified-ellipse fields.Source: Figure 5 · p. 20, §4.5.1 · p. 20
  • For both shape families, the CNN recovered principal perturbation locations and geometries, and PDE-constrained optimization refined its initialization using the governing PDE and displacement measurements.Smoothing and small spatial artifacts could remain in CNN predictions; optimization substantially improved a representative mollified-ellipse reconstruction.Source: §4.5.1 · pp. 23, 25, Figure 11 · p. 24, Figure 12 · p. 24
  • With 1% displacement noise, the pretrained CNN retained reasonable accuracy despite noiseless-only training, but optimization produced oscillatory artifacts that underscored the importance of regularization and stopping.The noisy-data results indicate that continued reduction of displacement misfit may introduce spatial variations associated with measurement perturbations.Source: §4.5.1 · pp. 25, 26, 28

4.5.2 Experiment II: reconstruction of combined Gaussian and mollified-ellipse structures

  • Experiment II evaluates whether a CNN trained on separate shape families generalizes to 500 noiseless coefficients combining Gaussian and mollified-ellipse perturbations, configurations absent from training.The experiment tests whether features learned separately transfer to their superposition; the supplied passages describe its error assessment but do not report the results.Source: §4.5.2 · pp. 28–30
  • Despite training only on individual inclusion types, the CNN identifies principal locations and approximate geometries in combined cases, though amplitudes, shapes, and transitions remain discrepant.The prediction provides an informative approximation across the spatial domain and is used to initialize PDE-constrained optimization.Source: §4.5.2 · p. 31, Figure 18 · p. 31
  • For a spatially separated combined-shape sample, optimization substantially improves overall coefficient recovery despite some localized artifacts.The comparison evaluates the true coefficient, the CNN initialization, and the optimized reconstruction.Source: Figure 19 · p. 32, Table 12 · p. 33, §4.5.2 · p. 33
  • PDE-constrained refinement improves CNN initial reconstructions for combined structures, including overlapping inclusions, by enforcing consistency with the forward model and measured displacement.The overlapping case tests transfer to spatial configurations more complex than the individual inclusion types used for training.Source: Figure 20 · p. 33, §4.5.2 · p. 34
  • The CNN was applied without retraining to combined-shape test coefficients using displacement data contaminated with 1% relative noise, then supplied as the optimization initial guess.This tests the CNN under both measurement noise and coefficient configurations absent from its training data.Source: §4.5.2 · p. 34, Figure 21 · p. 35
  • With 1% measurement noise, optimization improves the CNN initialization, sharpening separated inclusions in the real component and distinguishing partially overlapping structures in the imaginary component.The CNN was trained on noiseless samples containing individual inclusion types.Source: Figure 22 · p. 35, Table 14 · p. 36, §4.5.2 · p. 36

5 Conclusion

The two-dimensional hybrid framework uses CNN reconstructions to initialize adjoint-based PDE-constrained optimization for spatially varying complex shear modulus fields. Experiments show generalization to unseen inclusion configurations, with optimization refining noiseless reconstructions but becoming more noise-sensitive.

5 Conclusion

  • The CNN initializes nonlinear conjugate-gradient optimization with Armijo line search, which refines complex shear-modulus reconstructions and enforces consistency with the governing physical model.The framework combines data-driven initialization with adjoint-based PDE-constrained optimization.Source: §5 · p. 36
  • CNN initial reconstructions generalize across inclusion geometries, including configurations not encountered during training.Source: §5 · p. 36
  • In the noiseless setting, PDE-constrained optimization improves reconstruction accuracy; under moderate noise, the CNN retains dominant spatial features while optimization may introduce artifacts.The authors identify noise-aware regularization and stopping criteria as priorities for future work.Source: §5 · p. 36
  • The framework is evaluated in a two-dimensional MRE model, with future work targeting more heterogeneous geometries and realistic three-dimensional models.The authors also propose improved learning and optimization strategies.Source: §5 · p. 36
Loading 2610.02634v1…