Magnetic Resonance Parameter Mapping using Self-supervised Deep Learning with Model Reinforcement

Let $\mathrm{x}\in\ \mathbb{C}^{N_{x}\times N_{y}\times N_{z}\times N_{c}}$ be the measured MR image with a resolution $N_{x}\times N_{y}\times N_{z}$, where $x,y,z$ are spatial dimensions and $N_{c}$ denotes the number of coil elements in a typical multi-coil MRI reconstruction scenario. The corresponding undersampled k-space measurement ${\mathrm{f}\in\mathbb{C}}^{N_{kx}\times N_{ky}\times N_{kz}\times N_{c}}$ and $\mathrm{f}=\mathbf{U}\mathcal{F}\mathbf{C}\mathrm{x}+n$ incoproates the MRI encoding process, where $\mathbf{U}$ is an undersampling mask for accelerated acquisition, $\mathcal{F}$ is the discrete Fourier transform, $\mathbf{C}$ denotes the coil sensitivity map and $n\in\ \mathbb{C}^{N_{kx}\times N_{ky}\times N_{kz}\times N_{c}}$ is the signal noise during data acquisition, and $kx,ky,kz$ are k-space dimensions.

Supervised learning uses the fully sampled image $\mathrm{x}^{\ast}\ \in\ \mathbb{C}^{N_{x}\times N_{y}\times N_{z}\times N_{c}}$ as reference and the goal is to restore the image content from undersampled k-space data by learning network parameters $\Theta$ through minimizing the following loss function:

The loss function $\ell$ assesses the difference between the fully sampled $\mathrm{x}^{\ast}$ and the network output as the reconstructed image $\mathrm{x}\left(\mathrm{f}\middle|\Theta\right)$. Once trained with a large amount of data pairs $(\mathrm{f},\mathrm{x}^{\ast})$, the network can learn the mapping between undersampled measurements $\mathrm{f}$ and fully sampled data $\mathrm{x}^{\ast}$ for further inference.

In self-supervised learning, the network is trained to learn the mapping from undersampled k-space data $\mathrm{f}$ to itself. However, the MR physics model for image formation must be incorporated into the learning pipeline to guide the training process. The loss function in Equation (1) can be reformulated as:

Notably, as minimizing Equation (2) involves no fully sampled data and solely relies on learning the intrinsic structure and properties of the undersampled data itself, making it a powerful concept for deep learning methods with limited data [30, 31, 11].

To extend from MRI reconstruction to qMRI reconstruction, the MR parameter maps can be denoted as $\Delta=\left\{\delta_{i}\right\}_{i=1}^{N}$ where $\delta_{i}$ represents each MR parameter and $N$ is the total number of MR parameters to be estimated. Given the MR signal model $\mathbf{M}$, which is a function of $\Delta$, and the fully sampled image $\mathrm{x}^{\ast}$, the MR parameters $\Delta$ can be estimated by minimizing the following problem:

This is also referred to as model-based qMRI reconstruction. Here, the MR image $\mathrm{x}\in\mathbb{C}^{N_{x}\times N_{y}\times N_{z}\times N_{c}\times N_{k}}$ and its undersampled k-space data $\mathrm{f}\in\ \mathbb{C}^{N_{kx}\times N_{ky}\times N_{kz}\times N_{c}\times N_{k}}$ need to include $N_{k}$ measurements at varying imaging parameters dependent on individual qMRI method. The Equation (3) can be formulated as:

which incorporates the MRI encoding operators.

The above minimization problem (4) can be solved in self-supervised learning by minimizing the following loss function:

where the quantitative MR parameter $\Delta\left(\mathrm{f}\middle|\Theta\right)$ are parametrized from deep neural networks with learnable network parameters. In RELAX, minimizing Equation (5) with respect to trainable network parameters $\Theta$ is to train a deep neural network via network backpropagation for an end-to-end mapping network.

The proposed RELAX-MORE extends the RELAX framework by further considering unrolling the mapping network using an optimization algorithm. This new method is designed to solve a bi-level optimization problem, where the mapping network learning and quantitative MR parameter estimation can be jointly reformulated as:

Like RELAX, the upper level (6a) is a backward training process to optimize network parameters $\Theta$ by minimizing loss function $\ell$ via normal network backpropagation. However, unlike RELAX, the lower level (6b) can be implemented as a forward process to solve for the optimal solution of $\min\limits_{\Delta}{\phi\left(\Delta\right)}$ via a gradient descent based method. The objective function $\phi$ is a regularized least square model for quantitative MR parameter optimization. The weighted regularization $\lambda\mathcal{R}$ provides prior information pertinent to the desired MR parameter maps. Hand-crafted regularizations can be employed, such as total variation (TV) ( $\parallel\cdot\parallel_{TV}$ ) or $\ell_{1}$ norm ($\parallel\cdot\parallel_{1}$) to promote different aspects of image features. Recent studies show that regularization $\mathcal{R}_{\Theta}$ using CNNs can improve feature characterization, suppressing the hand-crafted priors [32, 33]. The weight $\lambda$ can be learned and integrated into $\mathcal{R}_{\Theta}$. The objective function $\phi_{\Theta}$ depends on the network parameters $\Theta$ that learned from regularizers, which can be parametrized as the summation of CNNs with $\ell_{2,1}$ norm:

Each ${\mathcal{D}_{\Theta}}_{i}$ learns to extract the prior features of different MR parameters $\delta_{i}$, for each $i=1,\ldots,N$. Therefore, the RELAX-MORE for qMRI reconstruction can be mathematically described as solving a CNN regularized bi-level optimization problem:

where data fidelity term $\mathcal{H}(\Delta)=\tfrac{1}{2}\parallel\mathbf{U}\mathcal{F}\mathbf{C}\mathbf{M}\left(\Delta\right)-\mathrm{f}\parallel_{2}^{2}$.

Proximal gradient descent, a widely used iterative algorithm, is implemented to unroll the lower level optimization (8b) in RELAX-MORE. For each unrolled iterative phase $t=1,\ldots,T$, the following steps is used to update MR parameter $\delta_{i}$ for each $i=1,\ldots,N$ :

where the proximal operator for regularization $\lambda\mathcal{R}(z)$ is defined as:

This step involves implementing the proximal operation for regularization $\mathcal{R}$, which is equivalent to finding the maximum-a-posteriori solution for the Gaussian denoising problem at a noise level $\sqrt{\lambda/\rho}$ [34, 35], thus the proximal operator can be interpreted as a Gaussian denoiser. However, because the proximal operator ${\rm prox}_{\mathcal{R}_{\Theta},\rho_{t}}$ in the objective function (8b) does not admit closed form solution, a CNN is used to substitute ${\rm prox}_{\mathcal{R}_{\Theta},\rho_{t}}$, where the network is constructed as a residual learning network $\mathcal{G}_{\Theta}$. Network $\mathcal{G}_{\Theta}$ is composed with CNN $\mathcal{D}_{\Theta}$, which reflects Equation (7), and its adjoint $\widetilde{\mathcal{D}}_{\Theta}$ with symmetric but reversed architecture to increase the network capability, and a soft shrinkage operator $\mathcal{S}_{\beta}$ is applied in between. All the substitutions are mathematically proven to be effective by [28, 36]. Then step (9b) becomes

where the soft thresholding operator is $\mathcal{S}_{\beta^{(t)}}={\rm prox}_{\parallel\cdot\parallel_{2,1},\beta^{(t)}}(y)=[{\rm sign}(y_{i})\max(|y_{i}|-\beta^{(t)},0)]$ for vector $y=(y_{1},\ldots,y_{n})\in\mathbb{R}^{n}$ and $\beta^{(t)}$ is a soft thresholding parameter that is updated at each phase $t=1,\ldots,T$. In summary, the forward learnable proximal gradient descent algorithm to achieve unrolled lower level optimization in RELAX-MORE can be detailed in Algorithm:

The soft thresholding parameter $\beta^{(t)}$ and step size $\alpha^{(t)}$ are learnable and updating during each phase. The final outputs of the algorithm are estimated intermediate MR images $\mathbf{x}^{(t)},\forall t=1,\ldots,T$ for all iterative phases and the reconstructed MR parameter maps $\delta_{i}^{(T)},\forall i=1,\ldots,N$ of the last phase.

Using an initial MR parameter input $\delta_{i}^{(0)}$ that is closer to the optimal solution of the Algorithm can lead to better results and eliminate the need for excessive iterations in gradient descent. To achieve a favorable input, initialization networks are used to learn $\delta_{i}^{(0)},\forall i=1,…,N$ from undersampled k-space data. First, the undersampled k-space data is converted into the image domain, followed by a Coil Combination Network to combine all coil-wise images. $N_{k}$ coil-combined images are produced to represent $N_{k}$ measurements during qMRI acquisition. Next, all images are concatenated together as input into the $\delta$-initialization Network to obtain $\delta_{i}^{(0)},\forall i=1,…,N$. The initial MR images can also be obtained from the signal model using $\mathrm{x}^{(0)}=\mathbf{M}(\{\delta_{i}^{(0)}\}_{i=1}^{N})$, which will be subsequently used in the loss function $\ell$ for training the overall network framework in upper level optimization in Equation (8a).

The loss function (8a) for RELAX-MORE is a summation of the weighted mean squared error (MSE) between undersampled k-space measurement $\mathrm{f}$ and the k-space representation of estimated MR images $\mathrm{x}^{(t)},\forall t=0,\ldots,T$ for all unrolled phases. Therefore, the overall network training objective in the upper level optimization is to minimize the loss function $\ell$ explicitly expressed as

where the $\gamma_{t}$ are the weighting factors for the MSE in order of the importance of the output from each phase and the initial $\mathrm{x}^{(0)}$.

The detailed framework of the proposed RELAX-MORE is illustrated in Fig. 1 The Coil combination Network and $\delta$-initialization Network applied 4 repeats of Conv($3\times 3$ kernel size, 32 kernels) followed by a ReLU activation function. The last convolution layer has a single kernel Conv($3\times 3$, 1). The CNN $\mathcal{D}$, applied 8 Conv($3\times 3$, 64) with a ReLU and the last convolution also has a single kernel Conv($3\times 3$, 1). The CNN $\widetilde{\mathcal{D}}$ applied a reversed network structure.

To investigate the feasibility of RELAX-MORE method for reconstructing accelerated qMRI, we used the widely studied $T_{1}$ mapping through the vFA method [13] as an example. Herein, the MR signal model $\mathbf{M}$ can be expressed as:

where several MR images are acquired at multiple flip angles $\eta_{k}$ for $k=1,\ldots,N_{k}$. $\mathbf{T}_{1}\in\mathbb{R}^{N_{x}\times N_{y}\times N_{z}},\mathbf{I}_{0}\in\mathbb{C}^{N_{x}\times N_{y}\times N_{z}}$ are the spin-lattice relaxation time map and proton density map, respectively, reflecting the imaged tissue relaxation properties, with $T_{1}$ value sensitive to many brain and knee joint diseases. The set of MR parameters estimated in this model are $\Delta=\{\mathbf{T}_{1},\mathbf{I}_{0}\}$. The flip angles and other imaging parameters, such as repetition time TR are pre-determined.

The experiments include in-vivo studies on the brain and knee of healthy volunteers and ex-vivo phantom studies, all of which were carried out on a Siemens 3T Prisma scanner. For the brain study, the purpose is to investigate the efficiency and performance of RELAX-MORE and compare it with other state-of-the-art reconstruction methods. The vFA on the brain of five subjects was performed in the sagittal plane using a spoiled gradient echo sequence at imaging parameters TE/TR = $12/40$ ms, FA = $5^{\circ},10^{\circ},20^{\circ},40^{\circ}$, FOV = $230\times 230\times 160$ mm and matrix size = $176\times 176\times 48$ with a dedicated 20-channel head coil. For the knee study, the purpose is to investigate the robustness of RELAX-MORE against severe noise contamination and image imperfection. Therefore, the vFA on one knee was performed using a 4-channel receiving-only flex coil at parameters TE/TR = $12/70$ ms, FA = $10^{\circ},20^{\circ},40^{\circ}$, FOV = $160\times 137\times 108$ mm and matrix size $224\times 192\times 36$. For the phantom study, the purpose is to investigate the reconstruction accuracy of RELAX-MORE. The vFA phantom data was acquired along the coronal plane with sequence parameters TE/TR = $12/80$ ms, FA = $5^{\circ},10^{\circ},20^{\circ},40^{\circ},60^{\circ}$, FOV = $170\times 170\times 60$ mm and matrix size $128\times 128\times 12$ using the 20-channel head coil. In all 3 studies, to minimize bias on $T_{1}$ estimation due to $B_{1}^{+}$ inhomogeneity and imperfect spoiling, $B_{1}^{+}$ maps were separately acquired for compensation [37], and $169^{\circ}$ linear RF phase increments between subsequent repetition cycles and strong gradient spoilers were applied to minimize the impact of imperfect spoiling [38].

The fully acquired MRI k-space data were undersampled retrospectively using two undersampling schemes: 1) 1D Cartesian variable density undersampling at acceleration factor AF = $3\times$ with the 16 central k-space lines fully sampled (Fig. 2(a)) and 2) 2D Poisson disk undersampling at AF = $4\times$ with the central $51\times 51$ k-space portion fully sampled (Fig. 2(b)). The undersampling patterns were varied for each flip angle, like in previous studies [18, 23]. The coil sensitivity maps were estimated from ESPIRiT [39].

The selection of network hyperparameters was non-trivial and optimized by empirically searching on a reasonably large range. More specifically, we used initial $\alpha_{i}^{(0)}=0.1,\beta_{i}^{(0)}=1e^{-5}$ for brain images; initial $\alpha_{i}^{(0)}=0.5,\beta_{i}^{(0)}=1e^{-7}$ for knee images; and $\alpha_{i}^{(0)}=0.9,\beta_{i}^{(0)}=1e^{-3}$ for phantom images. The loss function applied $\gamma_{t}=1$ at $t=T$ and $\gamma_{t}=1e^{-4}$ for other phase steps. The network parameters $\Theta$ were initialized using Xavier initialization [40]. Because RELAX-MORE trains on single-subject data, the batch size was set to include the entire image volume. The network was trained for $10,000$ epochs using the Adam optimizer [41]. The learning rate, which controls the step size taken in each iteration of optimizing the loss function, was set to $1e^{-4}$ for brain and knee images, and $5e^{-4}$ for the phantom image. It should be noted that, unlike conventional deep learning reconstruction methods where training and testing are two separate steps, in RELAX-MORE with subject-specific self-supervised learning, the reconstruction has readily concluded once the training has converged for one subject. All the programming in this study was implemented using Python language and PyTorch package, and experiments were conducted on one NVIDIA A100 80GB GPU and an Intel Xeon 6338 CPU at Centos Linux system.

The reconstruction performance of RELAX-MORE is affected by the degree of unrolling the gradient descent algorithm. Fig. 3 demonstrates the evolution of the estimated $T_{1}$ maps from fully sampled data at different total unrolling phase numbers from $T=3$ to $T=9$ of the gradient descent algorithm. The larger $T$ reflects deeper unrolling operation and thus requires more computing resources. Referring to the zoom-in views of Fig. 3, one can observe that with increasing total phase number, the algorithm starts to better distinguish the $T_{1}$ values from white matter (WM), grey matter (GM), and cerebrospinal fluid (CSF) regions. Initially blurry at a lower $T$, the tissue details start to sharpen with increasing $T$, resembling better the fully sampled $T_{1}$ map obtained using the standard pixel-wise fitting. However, the differences between phase $T=7,8$ and $9$ are negligible, indicating that performance gain can reach a plateau with a sufficient depth of unrolling steps. Therefore, to ensure consistent experiment setup and balance the trade-off between algorithm complexity and reconstruction performance, the number of $T$ was set to $8$ for all the experiments thereafter. This result illustrates the effectiveness of unrolling gradient descent algorithm in the RELAX-MORE framework.

RELAX-MORE was compared with two state-of-the-art non-deep learning qMRI reconstruction methods and one self-supervised deep learning method. These methods include 1) Locally Low Rank (LLR) [42], where image reconstruction was first performed, followed by pixel-wise parameter fitting. LLR exploits the rank-deficiency of local image regions along the acquisition parameter dimension to accelerate parameter mapping; 2) Model-TGV [43], a model-based qMRI reconstruction method that improves the piece-wise constant region restriction of total variation through a generalization of the total variation theory and 3) RELAX [23], an end-to-end self-supervised deep learning method for rapid MR parameter mapping. LLR and Model-TGV were implemented using the recommended optimization parameters in their original papers and code. Contrary to the original RELAX implementation, which was trained using many subjects, in the interest of a fair comparison, our implementation of RELAX was carried out on a single-subject training using an Adam optimizer with $10,000$ epochs.

$T_{1}$ maps estimated from $3\times$ 1D Cartesian variable density undersampling using the different methods are presented in Fig. 4. As can be observed in Fig. 4(a), the zero-filled $T_{1}$ map obtained by pixel-wise fitting the undersampled images is noisy, displaying ripple artifacts due to aliasing. Although LLR can partially remove these artifacts, it retains the noisy signature. Model-TGV averages out the noise, producing a cleaner tissue appearance and better $T_{1}$ contrast, but it is over-smoothed, resulting in blurry maps. On the other hand, RELAX removes noises and artifacts and produces sharpened maps but remains somewhat blocky. This is hypothesized due to the difficulty of converging the end-to-end network using single-subject data. RELAX-MORE produces artifact-free $T_{1}$ maps with anticipated excellent performance for noise removal, resulting in a good performance in appearance and contrast. This is further witnessed in the zoom-in view of Fig. 4(b), where the zero-filled, LLR, and even fully sampled images cannot reliably estimate $T_{1}$ in the cervical spine where the signal is low due to insufficient head coil coverage. Those pixel-wise fitting methods are typically prone to image noise contamination. Although Model-TGV and RELAX can better estimate $T_{1}$ in these low SNR regions, the blurry Model-TGV results make it difficult to distinguish between the spinal cord and subarachnoid space, whereas RELAX shows a disconnect in the subarachnoid space which may be due to over-sharpening. However, RELAX-MORE shows remarkably good performance in maintaining contrast and tissue details at those regions, enabling clear distinction between the spinal cord and subarachnoid space regions. In Fig. 4(c) of the cerebrum zoom-in view, RELAX-MORE consistently exhibits the best reconstruction performance in correcting artifacts, removing noises, and preserving tissue contrast and details. Referring to the error maps (Fig. 4(d)), which is the absolute difference between the estimated $T_{1}$ maps in Fig. 4(c) and the fully sampled $T_{1}$ map, overall, zero-filled shows the most significant error followed by LLR. Model-TGV and RELAX exhibit similar error maps, whereas RELAX-MORE produces the least error.

$T_{1}$ maps estimated from $4\times$ 2D Poisson disk undersampling are shown in Fig. 5. Comparing the images in Fig. 5(a), the overall signature differences are like the $3\times$ 1D undersampling case except for the undersampling artifacts being noise-like due to 2D undersampling. Referring to the zoom-in view and comparing the estimated $T_{1}$ maps again show similar differences to the $3\times$ 1D undersampling case. However, compared with other methods, RELAX-MORE can clearly distinguish between WM and GM $T_{1}$ of the cerebellum in Fig. 5(b), particularly in the posterior part. The absolute difference maps taken between the estimated $T_{1}$ maps in Fig. 5(c) and the fully sampled map, shown in Fig. 5(d), exhibit similar results to the $3\times$ 1D undersampling case with RELAX-MORE showing the least error. This is likely achieved through integrating deep learning capability and unrolling gradient descent algorithm with proximal operator prioritizing noise suppression without compromising the fidelity and clarity of the underlying tissue structure.

Further performance comparison was carried out using peak-signal-to-noise-ratio (PSNR), structural similarity index measure (SSIM), and normalized mean squared error (NMSE) as evaluation metrics for $T_{1}$ maps. PSNR is a measbure of the quality of the reconstruction, while SSIM is a measure of the similarity between two maps. PSNR, SSIM, and NMSE are defined in the following:

where $\mathrm{v},\mathrm{v}^{*}$ represent the estimated map and reference map. $\mu_{\mathrm{v}},\mu_{\mathrm{v}^{*}}$ are local means of pixel intensity, $\sigma_{\mathrm{v}},\sigma_{\mathrm{v}^{*}}$ are the standard deviations and $\sigma_{\mathrm{v}\mathrm{v}^{*}}$ is covariance between $\mathrm{v}$ and $\mathrm{v}^{*}$, $C_{1}=(K_{1}L)^{2},C_{2}=(K_{2}L)^{2}$ are two constants that avoid denominator to be zero, and $K_{1}=0.01,K_{2}=0.03$. $L$ is the largest pixel value of the magnitude of image.

$T_{1}$ maps obtained from fully sampled data were used for $\mathrm{v}^{*}$, as a reference to compare the performance among different methods. The metric calculations were carried out on results from $3\times$ 1D Cartesian variable density undersampling, using the brain regions from all five subjects.

The analysis results are presented in Fig. 6, where the average values of PSNR, SSIM and NMSE along with corresponding standard deviations among five subjects are presented using box-whisker plots. For both PSNR and SSIM, RELAX-MORE produces the highest mean value, followed by RELAX, Model-TGV and LLR. Both PSNR and SSIM of RELAX-MORE show the highest consistency, as evidenced by its small standard deviation compared to other methods. Zero-filled shows the lowest mean values for PSNR and SSIM with the largest standard deviation, which can be attributed to aliasing artifacts. RELAX-MORE shows the lowest mean NMSE, followed by RELAX, Model-TGV and LLR. Zero-filled shows the highest mean NMSE. Overall, in agreement with the qualitative observation in Figs. 4 and 5, RELAX-MORE performs superiorly in terms of reconstruction fidelity, structure and texture preservation, and noise suppression, outperforming other state-of-the-art methods.

An ablation study was carried out further to evaluate the parameter estimation accuracy of the proposed method. To analyze RELAX-MORE in discrete tissue environment, it was applied to undersampled data from a phantom composed of five 20 mL vials consisting of peanut oil and $2\%,4\%,8\%$ Agar dissolved in DI water and boiled egg white, all immersed in a 200 $\mu$M $\rm{MnCl_{2}}$ water bath (Fig. 7(a)). From left to right, Fig. 7(b) shows the $T_{1}$ maps obtained from zero-filled data, RELAX-MORE and fully sampled data. Zero-filled $T_{1}$ map exhibits ripple artifacts due to the undersampling, whereas these artifacts are eliminated in RELAX-MORE. Comparing the line profile indicated by the line in Fig. 7(b) and shown in Fig. 7(c), zero-filled $T_{1}$ profile resembles amplitude modulation, where the periodic amplitude fluctuation stems from aliasing effects and the higher frequency fluctuation reflects noise. In good agreement with the fully sampled $T_{1}$ profile, RELAX-MORE’s $T_{1}$ profile not only removes the periodic amplitude fluctuation but also smoothens the oscillations, demonstrating its ability to average out noise, remove aliasing artifacts while maintaining high accuracy for parameter estimation in a wide range.

The ablation study was also carried out to evaluate the parameter estimation robustness of RELAX-MORE against severe imaging imperfection and the generalizability of RELAX-MORE on different anatomies such as knee joints. Fig. 8 shows the $T_{1}$ maps from a representative sagittal slice of the knee obtained from zero-filled, RELAX-MORE, and fully sampled data. The MRI data was purposely acquired using a 4-channel receiving-only flex coil, which provides insufficient coil coverage and a poor signal-to-noise ratio for knee MRI. As evident in Fig. 8, the fully sampled $T_{1}$ map presents enormous noises at bone, cartilage, ligament, meniscus, and muscle and the fine structure, such as the epiphyseal plate (highlighted by a white arrow) is contaminated by noises. The zero-filled $T_{1}$ map from $3\times$ 1D Cartesian acceleration also presents severe undersampling artifacts combined with the noise rendering erroneous $T_{1}$ quantification. However, RELAX-MORE successfully removes all the artifacts, suppresses the unwanted image noises in $T_{1}$ quantification, and provides a surprisingly favorable quantification of different knee joint structures. RELAX-MORE demonstrates its high robustness and generalizability for reconstructing the in-vivo brain and in-vivo knee joint, making it a widely applicable method for different anatomies.

RELAX-MORE is a subject-specific method that utilizes self-supervised learning for efficient qMRI reconstruction. While RELAX-MORE can perform well using single-subject data, as shown in all other experiments, the reconstruction process (e.g., network training) needs to be conducted for each subject data. In this experiment, we investigated transfer learning to improve the computing time efficiency for the training/reconstruction of RELAX-MORE. Using transfer learning, the network weights after training on one brain data were then applied as the starting point to train the brain data of another subject.

In Fig. 9, the performance of RELAX-MORE with and without transfer learning is presented by plotting the SSIMs for $T_{1}$ estimations as a function of epoch number at $3\times$ 1D Cartesian acceleration. Comparing the two SSIM curves in Fig. 9(a), it is observed that transfer learning starts at a higher initial SSIM value during the initial training phase and increases faster compared to the case without transfer learning. This implies that transfer learning benefits the model to converge more efficiently due to commencing training with pre-trained weights. This is further confirmed by comparing the $T_{1}$ maps generated with (Fig. 9(b) top row) and without (Fig. 9(b) bottom row) transfer learning, where at low epoch values (500, 1000), transfer learning generated $T_{1}$ maps show better $T_{1}$ clarity and contrast between different tissue types compared to its non-transfer learning counterpart. As the training epoch increases, both SSIM curves flatten, indicating that the model has reached stable performance. The results suggest that transfer learning can improve the computing time efficiency of RELAX-MORE with reduced training/reconstruction time. At approximately 0.5 seconds per epoch using our GPU device for training the 3D brain data, transfer learning can reduce the reconstruction time to less than 10 min to reach good performance. With the advance of new techniques [44], transfer learning can be a practical approach to further improve the reconstruction timing efficiency for RELAX-MORE.