High-Dimensional MR Reconstruction Integrating Subspace and Adaptive Generative Models

Objective: To develop a new method that integrates subspace and generative image models for high-dimensional MR image reconstruction. Methods: We proposed a formulation that synergizes a low-dimensional subspace model of high-dimensional images, an adaptive generative image prior serving as spatial constraints on the sequence of “contrast-weighted” images or spatial coefficients of the subspace model, and a conventional sparsity regularization. A special pretraining plus subject-specific network adaptation strategy was proposed to construct an accurate generative-network-based representation for images with varying contrasts. An iterative algorithm was introduced to jointly update the subspace coefficients and the multi-resolution latent space of the generative image model that leveraged an recently proposed intermediate layer optimization technique for network inversion. Results: We evaluated the utility of the proposed method for two high-dimensional imaging applications: accelerated MR parameter mapping and high-resolution MR spectroscopic imaging. Improved performance over state-of-the-art subspace-based methods was demonstrated in both cases. Conclusion: The proposed method provided a new way to address high-dimensional MR image reconstruction problems by incorporating an adaptive generative model as a data-driven spatial prior for constraining subspace reconstruction. Significance: Our work demonstrated the potential of integrating data-driven and adaptive generative priors with canonical low-dimensional modeling for high-dimensional imaging problems.


I. INTRODUCTION
High-dimensional imaging problems emerge in various MR imaging scenarios, including quantitative MR parameter mapping (qMR) [1], dynamic imaging [2], and MR spectroscopic imaging (MRSI) [3], [4], which introduce additional dimension(s) to encode and decode tissue properties and provide quantitative biomarkers for pathological analysis.While the introduction of additional encoding dimensions enriches the available information, it also significantly lengthens the acquisition time, a major challenge for widespread clinical applications.In recent years, reducing the number of encodings and reconstructing images from sparsely sampled data using different model constraints have been recognized as a common approach to accelerate high-dimensional imaging acquisitions.
Subspace imaging/low-rank tensor modeling is one of the most investigated approaches to enable accelerated highdimensional MRI.While many variants have been presented, the essence is to exploit the redundancy in the imaging data due to the partial separability property that exists across various high-dimensional imaging modalities [5]- [7].Partial separability can be induced by the spatial-temporal correlations in dynamic MRI [8]- [10], spatial-parametric correlations in qMR [11]- [13], and spatial-spectral correlations in MRSI [14], which lead to accurate low-rank representations of the multidimensional imaging data arrays.As a result, highdimensional image reconstruction can be transformed into the recovery of a low-rank matrix/tensor which requires a significantly fewer number of measurements.
Low-rank modeling can be realized either implicitly through low-rankness encouraging regularization of high-dimensional matrices or tensors (e.g., nuclear norm as in [10], [15]- [17]), or explicitly through a matrix factorization form that directly reduces the number of unknowns during image reconstruction (e.g., [7], [14], [18]).Explicit low-rank modeling represents the signal evolution at each voxel as a linear combination of a small number of "temporal basis" with voxel-dependent "spatial coefficients" and enables effective incorporation of other prior knowledge to facilitate image recovery.In particular, special hybrid acquisition and physics-driven subspace learning strategies tailored for different applications can be designed to predetermine the basis [19]- [21].With the predetermined basis, complementary spatial constraints can be introduced to improve the estimation of the spatial coefficients.For example, analytical sparsifying transforms have been frequently used, e.g., [12], [22], [23].Kernel-based methods have been proposed to leverage image features derived from complementary anatomical imaging modalities to re-represent the spatial coefficients, going beyond the traditional voxel basis paradigm [24], [25].However, these hand-crafted "priors" did not fully exploit the prior information about the unknown images and required carefully feature/kernel selection which can be subjective and sometimes biased.
On a parallel line of pursuit, deep generative image models have provided a new and flexible way to utilize data-driven prior for image representation [26], [27], complementary to subspace modeling.State-of-the-art generative models, such as generative adversarial networks (GANs) [28], allow for automatic learning of feature-based image representations, to generate high-quality, high-resolution images and to constrain image reconstruction.For example, the unknown image can be modeled by a pretrained GAN (using images with similar contrasts and resolutions), and the reconstruction problem can be formulated as estimating a set of latent variables input to the GAN instead of the voxel values [29], [30].Major challenges of this approach are (1) representation accuracy, i.e., even with fully sampled images as the data, the pretrained GAN may not be able to accurately represent the images; (2) optimization issue, i.e., even if the model is highly expressive, GAN inversion is solving a deep network inversion problem which is highly non-convex and sensitive to initialization.
Different models combined with various optimization strategies have been proposed to address the above mentioned issues.Asim et al. proposed to use the invertible neural networks (INNs) to mitigate the representation error [31].As INNs use latent vectors with the same dimension as the unknown image, additional regularization can be introduced to improve reconstruction [32].Image-adaptive GAN (IAGAN) based strategies that update both the lower-dimensional latent vector and weights [30], [33] can reduce data-specific representation error (sometimes referred to as "out of distribution" error) but may lead to overfiting (e.g., to noise and other artifacts).While untrained networks can also be used as an image model (e.g., deep image prior [34], [35]) for updating only weights during reconstruction, such strategies do not fully utilize data-driven prior available in existing datasets, and the effectiveness highly depends on the network structure.Style-GAN offers better representation power with flexible latent space adjustments to account for data distribution variations, which demonstrates great potential in accurately representing images and effectively constraining the reconstruction [36], [37].
In this work, we proposed a new method that integrates subspace model and an adaptive generative image prior for reconstructing high-dimensional MR images.Direct application of generative models as a spatial constraint for subspace imaging can be challenging.Besides the representation accuracy and network inversion issues, training generative models typically requires large-scale, high-resolution, fully sampled datasets, which can be rather difficult to obtain for applications like qMR, dynamic MRI and MRSI etc.To address these issues, we introduced a pretraining plus subject-specific adaptation strategy for learning an accurate StyleGAN-based image representation.This adaptive prior accounted for contrasts and geometry variations between existing database and applicationspecific data effectively.Our method does not require end-toend supervised training using high-dimensional MR data since the GAN can be pretrained on publicly available data, and then adapted to a subject-specific reference image to account for geometry variations.The geometry-adapted GAN representation was incorporated as a spatial constraint into a subspacebased reconstruction formalism.An alternating minimization algorithm was introduced to solve the optimization of jointly updating the spatial coefficients and GAN latent space.We demonstrated the effectiveness of our method using qMR and MRSI data as application examples.
The remaining of the paper is organized as follows: Section II provides some background information on subspace imaging and generative model based image reconstruction.Section III describes the proposed problem formulation and algorithm in details.Section IV presents results evaluating the utility of the proposed method in two application examples, accelerated parameter mapping and MRSI.Section V and VI provide some technical discussion and conclude the paper.

II. BACKGROUND
The acquisition process (after proper discretization) for many high-dimensional MRI problems can be defined as: where ρ ∈ C N ×T corresponds to the discretized representation of the underlying high-dimensional image function with N being the number of voxels and T the number of "time" points to be determined, y ∈ C M ×T ′ denotes the measured data (T ′ may not necessarily be the same as T ), A and n represent the encoding operator and measurement noise, respectively.The goal is to recover ρ from measurements y which can be incomplete and/or noisy.

A. Subspace Reconstruction
While direct reconstruction of the high-dimensional ρ can be rather challenging and sensitive to noise perturbation, lowdimensional subspace models offer a way to significantly reduce the dimensionality of the problem by exploiting "correlations" among different dimensions, enabling image recovery from highly sparsely sampled and/or noisy data.More specifically, if ρ is a Casorati matrix [5] with spatial dimensions cascaded along the first dimension and a temporal/parametric dimension along the second, the subspace model can be realized as a matrix factorization form: ρ ≈ UV, where U ∈ C N ×R and V ∈ C R×T are low-rank matrices.U is often referred to as the spatial coefficient matrix and V as the temporal basis matrix.The model order/rank R, is typically much smaller than the ambient image dimensions (N and T ), which means that the number of unknowns/degree-of-freedom (DOF) is significantly reduced.
With the low-rank model, the image reconstruction problem can be formulated as where R(•) is a spatial/spatiotemporal regularization that incorporates prior information such as sparsity [7], [10], piecewise smoothness [38], [39] or structured low rankness [40] etc.This regularization, often including spatial constraints, plays a critical role for the success of subspace reconstruction, as it provides a complementary prior and effectively addresses the instability/ill-conditionness issues related to subspace fitting.Furthermore, for many imaging applications, the temporal basis V can be pre-determined using a problem-specific hybrid data acquisition strategy [7], [9], [12], [19]- [21].As a result, the problem in Eq. ( 2) can be simplified to only estimating U, circumventing the nonconvex problem of jointly solving for U and V.

B. Generative Model Constrained Image Reconstruction
Recent advances in deep generative models provide a new way to incorporate domain/data-specific prior into image reconstruction problems.With their state-of-the-art high-quality image generation performance, image representations using convolutional networks, GAN and extensions have been used in inverse problems in imaging.Specifically, a common approach is to represent the unknown image of interest x as x = G θ (w), where G is a generative network (can be pretrained using available domain-specific images or untrained), θ contains the network parameters, and w is a set of latent variables (often sampled from i.i.d distributions) that typically has a lower dimensionality than x.With this explicit generative model-based image function, the reconstruction problem can be formulated as which seeks to recover the image by determining the latent w and/or θ instead of solving voxel values in the conventional paradigm.This approach assumes that the unknown image lies in the range space of G θ , which may yield substantial modeling errors, depending on the network structure and capacity.This representation accuracy issue has been observed in the early works of compressed sensing using generative model (CSGM) [29], where a pretrained Deep Convolutional GAN (DCGAN) was used and the reconstruction process only updated w.IAGAN methods [30], [33] addressed this issue by simultaneously updating w and the network parameters θ to fit the data in Eq. ( 3).Early stopping is required to avoid overfitting (e.g., to the noisy measurements), particularly when the network has a high capacity.But determining the stopping criterion is heuristic, can be problem dependent and rather challenging.
The development of advanced generative models, such as StyleGAN2, has demonstrated superior image representation power than early generations of GANs, particularly for highresolution images [36].This makes it a powerful tool in image generation and representation.Many works have been proposed to utilize this representation power for image processing and reconstruction.For example, Kelkar et al. leveraged the multi-resolution latents of StyleGAN2 by constraining the reconstructed images to have similar features with a reference image [37].In the meantime, the high-dimensional, nonconvex network inversion problem associated with image reconstruction using GAN-based representation remains challenging and computationally expensive.There are also application specific challenges for high-dimensional MR imaging problems that are in need of formulation and algorithmic innovation.

A. Problem Formulation
We propose to integrate a generative image model as a spatial constraint and the subspace model for high-dimensional image reconstruction via the following formulation: The first term in Eq. ( 4) is the well-established subspace/lowrank model with the final reconstruction formed as ρ = Û V and V being the predetermined "temporal "basis.The second term introduces a GAN-based image representation as a spatial regularization for images at different time points (t = 1, 2, ..., N t ).The key assumption is that using a properly trained and adapted StyleGAN with fixed parameters in θ, the contrast variations in the image sequence can be accurately accounted for by updating only a set of low-dimensional latent space variables {w t } (referred to as the latents below; Note that w t denotes the entire multi-resolution latent variable for the tth image not a particular segment of the latents as described in [37]).This effectively serves as a data-adaptive reference image to constrain the spatial variations.An updating strategy to explore the multi-resolution latent structure and impose similarity for {w t } at neighboring time points will be specified in the Algorithm section.The final term R(.) is a hand-crafted regularization term, e.g., the commonly used sparsity constraint [7], [12].The subspace, sparsity and GAN prior play complementary roles to enhance the reconstruction performance.
The key challenges to integrate a GAN-based image prior are: 1) The representation accuracy of the GAN needs to be validated in application-specific context; 2) For many highdimensional imaging applications, there may be a lack of high-resolution, high-SNR images for training GAN; and 3) Even with a sufficiently accurate representation, solving the non-convex GAN inversion problem can be quite challenging.To address these issues, we proposed a pretraining plus subject-specific adaptation strategy to construct a StyleGAN prior G θ (•) adaptive to different imaging contexts, leveraging reference images readily available in a specific experiment.We also used an intermediate layer optimization (ILO) algorithm to address the GAN inversion challenge [41].More details are provided in the subsequent sections.

B. Adaptive Generative Image Model
While it is challenging to directly train a generative model that can accurately capture the geometry and contrast variations in every high-dimensional imaging problems (in many cases such large datasets are not available for training generative networks), a transfer learning plus application-specific adaptation strategy can be developed.Specifically, for a specific organ of interest, we can first train a StyleGAN2 model on large, publicly available data, e.g., the Human Connectomme Project (HCP) database for brain imaging applications [42].We chose StyleGAN2 due to its state-of-the-art performance in generating high-quality, high-resolution images, compared to older generations of DCGANs, and its superior image representation capability due to higher dimensional multiresolution latents.
With the pretrained model, parameterized by θ p , we proposed a subject-specific adaptation strategy to construct the generative image representation for a specific reconstruction task.The key assumptions are 1) there exists a high-quality reference image for network geometry adaptation, and 2) the adapted network can accurately represent images with various contrasts (often acquired in high-dimensional imaging problems) by tuning only the latent space.Given a subject-specific high-resolution reference image x p (e.g. from a reference T1w anatomical scan), the adaptation step can be formulated as: where both the latents w p and network parameters θ are updated to minimize the representation error of the network for x p .The second term is introduced to penalize large deviation of the adapted network from the pretrained one.To reduce the optimization error and fully exploit the latent-space structures, we proposed to solve Eq. ( 5) with a two-step algorithm, which alternates between: and The latents ŵp were solved by the ILO algorithm with fixed network parameters θ p (details provided below), and then the parameters were updated with fixed latents ŵp .This pretraining plus adaptation strategy produces an effective adapted GAN prior for which updating only the low-dimensional latents can accurately account for contrast variations, as demonstrated in Fig. 1 using multicontrast brain MR images.The adaptive GAN can serve as a prior in many high-dimensional MR imaging applications, i.e., in this case introduced as a spatial constraint to enhance the subspace reconstruction.

C. Algorithm
With the adapted G θ (•), we proposed to solve the problem in Eq. ( 4) using an alternating minimization algorithm.Specifically, the original reconstruction problem is decoupled into two subproblems, i.e., one StyleGAN inversion problem and the other GAN-constrained spatial coefficient update problem, which can be mathematically expressed as: Subproblem (I): Update latent w by solving and Subproblem(II): Update spatial coefficients U by solving where i denotes the iteration number.Since the StyleGAN prior has already been adapted to the subject-specific reference image, for Subproblem (I) in Eq. ( 8), we minimized the image representation error with respect to the latents only.To address the optimization challenge for GAN inversion, we adapted the ILO algorithm considering the unique latent space structure in StyleGAN (similar to solving Eq. ( 6)).More specifically, instead of directly solving Eq. ( 8) w.r.t.{w t }, the loss was first minimized w.r.t. the latents in an intermediate layer and then projected back to the range space of the previous layers.
Our procedure is illustrated in Fig. 2. Mathematically, we can denote the network as a composition of G θ = G θ2 •G θ1 , with θ1 and θ2 containing parameters for different partitions of the network (Fig. 2).The multiresolution latents in StyleGAN, w t , can then also be divided Subsequently, the update for the lower-resolution latents w 1,t was obtained by fitting This step essentially "projects" the intermediate layer output ŵint,t to the range space of G θ1 by representing it with the lower resolution latents { ŵi+1 1,t }.The resulting { ŵi+1 1,t }, { ŵi+1 2,t } can then be combined into an initialization for the problem in Eq. ( 8) to further minimize the loss.In this work, to fully explore the multi-resolution latent structure in StyleGAN2, the network was divided into multiple partitions (with multiple intermediate latents) to break down the original challenging one step deep network inversion problem into multiple steps.Specifically, the minimization starts from a certain intermediate layer depending on the application and then backpropagated to the earlier layers (see the Experiments and Results section).
Additionally, we can assume that images sharing similar features, such as geometry and contrast, are likely to have similar latents (often demonstrated via "style mixing" experiments in the computer vision literature).Exploiting this assumption, we can mitigate potential overfitting to noise/artifacts by introducing a similarity constraint in the latent space across images with different contrasts.For example, we can introduce an l 1 ball constraint on the latents for adjacent "time" points, which can be mathematically formulated as: Here ŵi+1 t−1 ⊕ B t (l t ) denotes an l 1 ball with radius l t centered at ŵi+1 t−1 from an adjacent "time" point (e.g., an image corresponding to a previous TE).The radius was determined empirically as hyperparameters (see the Discussion and Conclusion Section).Another strategy is to enforce latents at a certain resolution to be consistent when representing images across various contrasts (application dependent).
For Subproblem (II), with the updated GAN representation, different regularization choices can be made for R(U V) in Eq. ( 9) in an application-specific context.This is essentially the well-studied regularized subspace reconstruction, but incorporating additional, reference-adapted network constraints.The initialization of alternating minization algorithm can be application specific, which will be discussed in the Experiments and Results section.

D. Training and Other Implementation Details
For the generative network, we used the StyleGAN2 architecture proposed in [36] with a latent space dimensionality of 512.Considering the complex-value nature of MRI data, we trained a magnitude image network and a phase image network separately.The magnitude network was trained on HCP database [42] with 120,000 T 1 -weighted & 120,000 T 2weighted images.The phase network was trained on NYU fastMRI database [43] with 70,000 phase images after coil combination [44].The network training, subject-specific adaptation and StyleGAN2 inversion were performed on a Linux server with NVIDIA A40 GPU and implemented in PyTorch 1.12.1.We used Adam optimizer [45] with batch size 80 for training while other hyperparameters remained unchanged as described in [36].The subspace reconstruction and other methods for comparison were implemented in Matlab R2020b.

IV. EXPERIMENTS AND RESULTS
We evaluated the utility of the proposed method in two high-dimensional MR imaging applications: quantitative MR parameter mapping and MR spectroscopic imaging.Specific GAN adaptation considering the unique acquisition design in each imaging scenario, as well as modification to the overall formulation will be discussed.Improved reconstruction performance over state-of-the-art subspace reconstruction methods in each case will be demonstrated.All in vivo studies were performed with local IRB approval.

A. Accelerated MR Parameter Mapping
qMR provides quantitative tissue properties beyond traditional subjective and qualitative contrast-weighted images for diagnosis and disease characterization.One major challenge for clinical applications of qMR is the prolonged acquisition due to the need to acquire multiple images for parameter estimation.In this section, we demonstrate the effectiveness of our method for accelerating qMR (i.e., T 2 mapping) experiments.
In vivo T 2 mapping data were acquired using a multispin-echo sequence on a 3T scanner (Siemens Trio) using a 12-channel head coil.Specific acquisition parameters are: 16 echoes with TE 1 = 8.8 ms and echo spacing ∆TE = 8.8 ms, TR = 4000 ms, slice thickness = 3 mm, matrix size = 192 × 192 and a FOV of 192×192 mm 2 .The acquired fully-sampled data were retrospectively undersampled using a 1D random phase encoding pattern at different acceleration factors (AFs).The central 12 k-space lines were acquired at all TEs for subspace determination [12], [22], and the central k-space was fully sampled at the first TE with coverage dependent on AFs, e.g., 48 for AF = 4, for coil sensitivity estimation by ESPIRiT [46].The data and reconstruction workflow for this application is further illustrated in Fig. S1 in the supplementary information.
To obtain a reference image for the proposed subjectspecific GAN adaptation (Eq.( 5)), we used a sum-of-squares (SoS) combination of images across all TEs from an initial subspace reconstruction.Based on our observation, the SoS image exhibits negligible aliasing even at high AFs.While additional high-resolution reference images can be acquired and used, e.g., a T 1 w image typically acquired in neuroimaging experiments, our strategy might alleviate the need of additional acquisitions.The coil sensitivity maps were integrated into the forward encoding model for reconstruction from multichannel data and the specific reconstruction formulation becomes: where Ω, F and S c denote the sampling operator, Fourier encoding matrix and sensitivity encoding matrix respectively.The additional regularization chosen here was a joint sparsity constraint with D a finite difference operator [12].We used magnitude network only and the phase {Φ t } came from a subspace reconstruction with joint sparsity constraint (Assuming the smooth phase from subspace reconstruction was good enough).The parameters {λ 1,t } can be TE-dependent: we used 0.2 for the last five TEs and 0.1 for the other TEs in the reconstruction using subspace and adaptive GAN constraint alone.When using all three constraints, we simply chose 0.04 for all λ 1,t and achieved the best reconstruction results.For the joint sparsity contraint, we used a mild choice of λ 2 = 2e −7 (λ 2 = 1e −6 when using joint sparsity constraint alone).The algorithm converged after five iterations.More discussion on parameter selection can be found in the Discussion section.
When applying the proposed reconstruction, a similarity constraint in latent space was introduced for images at adjacent TEs in Subproblem (I) (as shown in Eq. ( 12)).For initialization, we started the alternating minimization by first updating the latents (Eq.( 8)).More specifically, we used a Û0 from a subspace reconstruction with a mild sparsity constraint, and then further updated the latents (before moving to Subproblem II) by solving the following problem each TE using k-space data directly to ensure data consistency: A set of multicontrast images reconstructed from data sampled with AF=4 are shown in Fig. 3. Introducing either the sparsity constraint or the adapted GAN prior improved upon the subspace reconstruction, as expected, with the GAN prior producing a slightly lower image-domain error.The proposed method integrating both the constraints yielded the best reconstruction, as shown in the error images (right panel) and the overall relative ℓ 2 errors.The T 2 estimation results corresponding to the same data (AF=4) are compared in Fig. 4, demonstrating effectiveness of the proposed method.More quantitative comparisons across different AFs are shown in Fig. 5.The proposed method consistently achieved lower errors at different AFs.

B. High-Resolution MRSI
In this section, we demonstrated the effectiveness of the proposed method on another high-dimensional imaging problem, i.e., MRSI.The goal here is to evaluate our method's utility for SNR-enhancing reconstruction from noisy, high-resolution MRSI acquisitions, benchmarking against the state-of-the-art subspace reconstruction (SPICE).
As it is rather difficult to obtain "reference" high-SNR data for high-resolution MRSI, we performed both numerical simulations (with ground truth) and in vivo experiments for evaluation.In vivo MRSI data were acquired from healthy volunteers using a fast sequence described in Ref. [47], with TR/TE = 1100/30 ms, matrix size = 64 × 64 × 12, FOV = 220 × 220 × 64 mm 3 , spectral bandwidth (BW) of 1.67kHz and 320 echoes (spectral encodings).Our sequence also generated a set of interleaved high-resolution, high-SNR water spectroscopic data, which is an excellent data for simulation-based evaluation, particularly for the adaptive GAN representation for different contrasts.For simulations, we first interpolated the water spectroscopic images to a grid of 128 × 128 × 12 × 150 (where 150 is the number of samples along the FID dimension) and used them as the ground truth to simulate noisy water MRSI data (Fig. 6).The noisy metabolic MRSI data from the same sequence were used to evaluate the proposed method for metabolite spatiospectral reconstruction.
In MRSI experiments, a high-resolution, anatomical scan is typically acquired and can be used as a reference for GAN adaptation.Here, a 3D T 1 w image (from an MPRAGE scan) was acquired.
The StyleGAN2 network was pretrained using the same HCP data but resized to 128 × 128 (as MRSI data have a lower resolution than T 2 mapping) and skull stripped.The pretrained network was first adapted to the T 1 w reference image to generate G θ (•) for the proposed reconstruction.The phase network G θϕ (•) pretrained on NYU fastMRI data was introduced to account for the complex-valued data.Considering that a) the MRSI data have a significantly larger number of "time points" than parameter mapping and b) the images Fig. 6.
Simulation results for water MRSI reconstruction: (Row 1) Ground truth images simulated from interpolated water spectroscopic data (128 × 128 × 12 × 150 spatiospectral encodings); (Row 2) Noisy images generated by adding Gaussian noise to ground truth; (Rows 3 and 4) Reconstruction results from the SPICE subspace reconstruction and proposed method (Proposed).The proposed method produced better contrast, a higher SNR and a lower reconstruction error.Fig. 7.In vivo MRSI reconstruction: (Row 1) T 1 -weighted anatomical images (T 1 w, serving as the reference for network adaptation); (Rows 2-4) Metabolite maps for the corresponding slices from noisy data (left three columns), SPICE (middle three columns) and the proposed method (right three columns), respectively.Both SPICE and the proposed method produced significant SNR enhancement for this high-resolution acquisition, revealing tissue specific metabolite variations, but better contrast and less artifacts can be observed in the proposed reconstruction.The last two rows show two voxel spectra from respective methods (voxel locations marked in the T 1 w images).Similar spectral quality was achieved by the SPICE and proposed methods (not surprising as both methods used the same subspace).
at individual time points (especially later echoes) are much noiser, instead of representing images at different echoes using the adapted StyleGAN2, we proposed to use it to represent the spatial coefficient maps (for subspace reconstruction) as they can be treated as "images" with different contrasts.The specific formulation can be written as: where U r is the r th column in U (spatial map of each coefficient) and R is the number of columns (aka the rank), B models the field inhomogeneity phase effects and D w denotes an edge-weighted finite difference operator.This led to a slightly different Subproblem I than the previous application example, where the latents for magnitude network ({w r }) and phase network ({w ϕ,r }) were determined alternatively.We performed a SPICE-based subspace reconstruction with a mild edge-preserving regularization as the initialization for Û0 .The high-resolution latents for the StyleGAN2 (inputs for the last layer) were kept unchanged (after adaptation in Eq. ( 7)) during the GAN inversion subproblem to alleviate overfiting to noise considering the low SNR.The effects of changing latents at different resolution are further illustrated in the Discussion.We chose λ 1,r = 1.6 and λ 2 = 0.3 for proposed method (λ 2 = 0.6 for the SPICE reconstruction with a similar data consistency level).The ILO started from the 4th layer of the network.
Reconstruction from the simulated noisy water images supported the idea of using adaptive GAN to represent the spatial coefficient maps.The subspace for water reconstruction was obtained by applying an SVD to the "ground truth" water spectroscopic images.As shown in Fig. 6, the proposed method achieved a significantly lower reconstruction error and noticeably better "contrast recovery" than the SPICE reconstruction.Impressive SNR enhancement was achieved by The effects of iteration for our proposed algorithm: Reconstruction errors of the multi-TE images for a T 2 mapping data set w.r.t. the number of iteration (AF=4).The error changes substantially for the first 2-3 iterations, but minimally after 5, with and without integrating the adaptive GAN and joint sparsity constraints.
both.Fig. 7 shows a set of metabolite spatiospectral reconstruction from an in vivo MRSI acquisition.The metabolite maps from the proposed reconstruction exhibit a higher SNR and clearer contrast compared to the noisy data and subspace reconstruction (SPICE).Note that increasing the spatial regularization parameter for the SPICE reconstruction can achieve a similar level of SNR enhancement but at the expense of oversmoothing.The regularization parameters were tuned to a similar data consistency level, matching noise energy for both methods.Therefore, the proposed method provided a better trade-off in SNR and resolution by incorporating an adaptive GAN-based spatial constraint.The subspace used for the MRSI reconstruction was pre-estimated and adapted to this specific data (using the higher-SNR portion) as described in [48]- [50].

V. DISCUSSION
We presented a new adaptive generative model based spatial constraint for subspace reconstruction of high-dimensional imaging data.One unique advantage of our proposed method is that the subject-specific adaptation step allows our GAN model to be flexibly adapted to different imaging contexts, even without high-quality, application-specific training data.We validated that the subject-specific, adapted GAN (to an application-specific reference image) can be accurate representation of images from the same subject with different contrasts (a common feature in high-dimensional imaging) by updating latent variables alone.This circumvented the overfitting issue in methods updating both latents and network parameters during reconstruction, such as IAGAN.We demonstrated the effectiveness of our method in two application examples.
There are several points worth discussing.The first is the selection of regularization parameter(s) (and other hyperparameter(s)), particularly when our adaptive GAN prior is combined with other regularization functional.To alleviate possible overregularization/oversmoothing introduced by other handcrafted constraints, we chose a small λ 2 in this work for a mild additional regularization effect, and our proposed GAN constraint played a more important role.For {λ 1,t }, selection was based on a combination of discrepancy principle and visual inspection.We also observed that in our T 2 mapping experiments, similar {λ 1,t } worked well across different AFs, and a relatively large range of {λ 1,t } can result in similar reconstruction results.One advantage is that {λ 1,t } can be "time" specific.More specifically, when using adaptive GAN only to constrain the subspace reconstruction, a 2× larger regularization parameter for the last five TEs (due to relatively lower SNRs) can slightly improve the T 2 mapping reconstruction.Different combinations of {λ 1,t } may be further explored.Similarly, the determination of the radius l t that imposes the l 1 ball constraint (Eq.( 12)) can be "time" specific and further explored.The second is the initialization step.We observed that the reconstruction can be substantially affected by initialization which should be carefully chosen for different applications.For our T 2 mapping experiments, starting the iterations using Û0 from a subspace reconstruction (with a mild joint sparsity constraint) may introduce potential bias.Therefore, we proposed to further update the latents using the data consistency loss described in Eq. ( 14) which led to a better performance.For MRSI, we found that initializing the algorithm with Û0 from a initial subspace reconstruction worked well.
Another issue is the reconstruction speed.Currently, we used five outer iterations (Fig. 8) that achieved empirical convergence and satisfactory reconstruction results.Each iteration needs about 40 minutes (can be application dependent).The most computationally expensive step is Subproblem (I), where multiple gradient descent updates are required to explore the intermediate latent space.Currently, we solved the latents TE by TE for MR parameter mapping, for which more efficient implementations should be considered.Other methods that better exploit parallel computation to simultaneously update the latents for all TEs and even circumvent network inversion will be investigated.
While our work and previous work [37] using GAN for image reconstruction have demonstrated that imposing constraints on latents can help better utilize prior information from reference images and reduce overfitting to artifacts, an important issue that presents unique research opportunities is the properties of latents (feature) space.Here, we further investigated the latent space structure through a style-mixing experiment on the adapted StyleGAN2 (after pre-training and network parameter adaptation to experimentally acquired T 1 w image).Fig. 9 demonstrated that manipulating low-resolution latents can introduce more global and significant changes in structure and contrast of the generated images than changing the high-resolution latents.This is consistent to us keeping the high-resolution latents unchanged for our MRSI reconstruction (lower-resolution data and avoiding overfitting to noise).In the meantime, it can be observed that for the current StyleGAN representation, different types of image features (e.g.contrast and geometry) are still somewhat entangled, and it is challenging to control certain semantic features by changing latents at a specific scale.Future work can focus on developing disentangled representation to enable more explicit control of different image features, which will be desirable for high-dimensional MR imaging problems, since in many situations, only parts of image features are changing across specific dimensions.
Recent development of generative models may also be helpful to improve the accuracy and efficiency of image representation, and thus improve the reconstruction performance.Transformer-based GAN has been proposed [51] to better utilize the global correspondence of image features through self-attention mechanism.Korkmaz et al. [52] demonstrated that generative vision transformers were more robust to artifacts and achieved better performance than CNN-based generative models when worked as a deep image prior (without pretraining) [34].Diffusion model is an emerging approach to generate high-quality images.Diffusion model training is more stable than GAN and can obtain better sample quality over the state-of-the-art GANs [53], [54].While diffusion models have demonstrated its potential for image reconstruction [55], [56], many opportunities remain on how to incorporate diffusion models into high-dimensional MR imaging problems [57].

VI. CONCLUSION
We proposed a novel reconstruction method for highdimensional MR imaging that integrated subspace modeling and an adaptive generative-network-based image prior.The proposed adaptive generative model served as an accurate representation of images from the same subject with different contrasts and an effective spatial constraint for subspace reconstruction.The complementary powers of subspace, generative image constraint, and sparsity regularization produced improved performance over state-of-the-art subspace methods in two application examples.We believe our work offers a new perspective on high-dimensional image reconstruction via integrating learning-based spatial priors and low-dimensional modeling.

Fig. 1 .
Fig. 1.Validation of the proposed adaptive StyleGAN representation: (Row 1) The pretrained StyleGAN2 can generate high-quality brain images with different geometries and contrasts; (Row 2) Representative images with different contrast weightings from an independent experiment (not from the database).The first image served as a reference for the GAN adaptation (see texts); (Row 3) Images from the range space of the pretrained StyleGAN2 that were closest to the images in Row 2 in terms of l 2 error, from updating the latents of the StyleGAN2.Noticeable representation errors can be observed; (Row 4) The adapted GAN achieved significantly better representation accuracy.Specifically, updating both latents w and network parameters θ (left) yielded an image almost the same as the reference image.With the adapted subject-specific network, updating latents only accurately accounted for contrast variations (right three figures in Row 4).These results support using the adaptive StyleGAN2 as an accurate and flexible image representation in a constrained reconstruction.

Fig. 2 .
Fig.2.Illustration of the network partition in the ILO algorithm: The original StyleGAN G θ can be divided into a cascade of two sub-networks, G θ1 (lighter yellow) and G θ2 (darker yellow), with corresponding latents w 1,t and w 2,t .The output of G θ1 , w int,t (originally hidden in G θ ), can be viewed as intermediate latents and the input to G θ2 , which will be first updated in the GAN inversion process (see texts for details).

Fig. 3 .
Fig. 3. Reconstructed multi-TE images (Left panel) by different methods from an in vivo data set and the corresponding error maps (Right panel) at AF=4.Images for different TEs were shown in different columns while results for different methods in respective rows.The adaptive GAN representation worked as an effective spatial constraint that improved the subspace reconstruction (Row 4, Subspace+GAN constraint).The proposed method integrating all components outperformed the state-of-the-art subspace+sparsity method and achieved the lowest error (Last row), better visualized in the error images.

Fig. 4 .
Fig. 4. Estimated T 2 maps (Row 1) from the data in Fig. 3 and the corresponding error maps (Row 2) for different reconstruction methods (AF=4).The overall T 2 estimation errors are shown in the upper left corners of the images, which were calculated for the brain region only.As can be seen, reconstructions using joint sparsity constraint and GAN constraint yielded similar reconstruction errors.Integrating both sparsity and GAN constraints produced the best result.

Fig. 5 .
Fig. 5. Reconstruction errors at different AFs for the T 2 mapping application: The proposed method shows lower errors across all AFs in both the reconstructed contrast-weighted images (left) as well as estimated T 2 maps (right).
was trained with 192 × 192 images and 7 layers in total to match the reconstruction resolution of the T 2 mapping data.The ILO algorithm started from the 5th layer.

Fig. 8 .
Fig. 8.The effects of iteration for our proposed algorithm: Reconstruction errors of the multi-TE images for a T 2 mapping data set w.r.t. the number of iteration (AF=4).The error changes substantially for the first 2-3 iterations, but minimally after 5, with and without integrating the adaptive GAN and joint sparsity constraints.

Fig. 9 .
Fig.9.A style mixing experiment for our adapted StyleGAN2: A single image (Source A) and a set of images (Source B) were generated from the network; The rest of the images were generated by replacing a specific subset of latents in Source A by corresponding latents from Source B. Specifically, for each row, images in different columns were obtained by taking the same level of latents from different images in Source B. From the top to bottom rows, different levels of latents from Source B were used (for the same image A), from low-resolution to highresolution subsets, respectively.As can be seen, lower-resolution latents from Source B are more related to global features such as shape and contrast, while higher-resolution latents lead to more subtle changes in fine details and less/negligible contrast and geometry variations.