Study of key properties behind a good undersampling pattern for quantitative estimation of tissue parameters
Riwaj Byanju1, Stefan Klein1, Alexandra Cristobal Huerta2, Juan A. Hernandez Tamames2, and Dirk H. J. Poot1

1Departments of Medical Informatics and Radiology, Erasmus MC, Rotterdam, Netherlands, 2Departments of Radiology and Nuclear Medicine, Erasmus MC, Rotterdam, Netherlands


Quantitative MR (qMRI) at present is clinically unfeasible due to long scan time. Jointly performing image reconstruction and parameters estimation is expected to allow increased acceleration. In this work, we investigate properties of undersampling patterns that are most relevant for parameter estimation using a Cramer-Rao-Lower-Bound (CRLB) based metric for such an approach. We compare key properties of undersampling patterns and conclude that one of these properties, namely low discrepancy, is most relevant for achieving time-efficient qMRI.


Quantitative MR (qMRI) has great potential, however, long scan time currently makes it unfeasible for clinical application. Conventionally in qMRI, images are reconstructed with acceleration rates limited by parallel imaging, followed by parameter quantification from those images. Jointly performing image reconstruction and parameters estimation is expected to allow increased acceleration. In this work, we investigate properties of under sampling strategies that are most relevant for parameter estimation with such approach. We use a Crame′r-Rao-Lower-Bound (CRLB) based metric for comparing different types of under sampling strategies.


Let $$$f_j(\boldsymbol\theta_{\boldsymbol{x}})$$$ predict the signal of contrast $$$j\in[1 ... J]$$$ of an acquisition scheme that acquires $$$J$$$ different contrasts with parameter vector $$${\boldsymbol\theta_{\boldsymbol{x}}}$$$. Here $$${\boldsymbol\theta_{\boldsymbol{x}}}= [M_0, T_1, T_2]_\boldsymbol{x}$$$ of voxel $$$\boldsymbol{x}$$$ where $$$\boldsymbol{x}\in\Omega$$$ is the spatial position index in image domain $$$\Omega\subset\mathcal{N}^3$$$, consisting of $$$N$$$ voxels in total, $$$\boldsymbol{k}\in\Omega$$$ is the k-space index and $$$p\in[1...P]$$$ indexes $$${\boldsymbol \theta_{\boldsymbol{x}}}$$$. Let $$$F$$$ be the Fourier transform, and $$$C_{\boldsymbol{x}, c}\in\mathcal{C}$$$ indexed as $$$c\in [1 ... C]$$$ be the coil sensitivity map. Then with an undersampling pattern, $$$U_{j,\boldsymbol{k}}\in\{0, 1\}$$$ the expected value for a k-space sample $$$\mu_{j,\boldsymbol{k},c}(\theta_{\boldsymbol{x}})$$$, is given by

$$\mu_{j,\boldsymbol{k},c}(\boldsymbol{\theta}) =U_{j,\boldsymbol{k}} {\sum_{\boldsymbol{x}\in\Omega}} F_{\boldsymbol{k},\boldsymbol{x}} C_{x,c} f_j(\boldsymbol{\theta}_{\boldsymbol{x}}) .$$

Here, $$$\boldsymbol{\theta}$$$ represents the concatenation of $$$\boldsymbol{\theta_{\boldsymbol{x}}} \forall \boldsymbol{x}\in\Omega$$$. We can derive the CRLB, which quantifies how noise on the collected k-space samples propagates to uncertainty in the estimated parameters maps $$$\boldsymbol{\theta}$$$:

$$\begin{align}\textbf{CRLB}(\boldsymbol \theta) &= I^{-1}(\boldsymbol \theta) \\I(\boldsymbol{\theta})_{(\tilde{p},\tilde{\boldsymbol{x}})(\hat{p},\hat{\boldsymbol{x}})} &= 2(\frac{\Gamma}{| \Gamma|})\sum\limits_{j=1}^{J} \sum\limits_{c=1}^{C} \sum\limits_{\boldsymbol{k}\in\Omega} \frac{\partial \mu_{j,\boldsymbol{k},c} }{\partial \theta_{\tilde{p},\tilde{\boldsymbol{x}}}} \frac{\partial \mu_{j,\boldsymbol{k},c}}{\partial\theta_{\hat{p},\hat{\boldsymbol{x}}}} \end{align}$$

where ($$$\tilde{p}, \tilde{\boldsymbol{x}})$$$ and ($$$\hat{p}$$$,$$$\hat{\boldsymbol{x}}$$$) represent row and column indices corresponding to parameters $$$p$$$ at $$$\boldsymbol{x}$$$, $$$\Gamma$$$ is the complex variance of the noise and $$$I$$$ is Fisher information matrix.

Based on the CRLB and the acquisition time $$$T_{acq}$$$, we compute for each parameter $$$p$$$ the time efficiency $$$\eta_p$$$, which will be used to compare undersampling patterns: $$$ \eta_p= \frac{1}{T_{acq} \textbf{mean}_{\boldsymbol x} \textbf{CRLB}_{(p,\boldsymbol{x})(p,\boldsymbol{x})} (\boldsymbol{\theta}) }$$$

Undersampling patterns

As there are exponentially many $$$U$$$ an exhaustive search is impossible. Hence, we define patterns on the basis of key properties in order to study the features that make an optimal pattern. Key properties are:

  1. StdNK: Standard-deviation of the number of times k-space position is sampled
    $$\text{SNK:} = std_k \left(\sum_{j=1}^{J} U_{j,\boldsymbol{k}}\right)$$
  2. Discrepancy: Quantifies the deviation from regular sampling and is quantified by the L2 discrepancy.1

Evaluated undersampling patterns

We study the influence of key properties on time efficiency with following patterns, shown in figure 1.

  1. Regular: 2D regular undersampling with acceleration factors $$$R = R_m \times R_n $$$ with $$$R_m$$$ and $$$R_n$$$ ($$$m$$$ $$$\And$$$ $$$n$$$ $$$ \in [1,2,3,4,6,8,12])$$$ the acceleration in the two phase encoding dimensions. For each R with multiple $$$R_m$$$, $$$R_n$$$ pairs, the pair with the best $$$\eta_p$$$ is chosen.
  2. Translation: To modify the phase of aliasing (folding) in the echoes, each echo is translated Regularly.
  3. TS: To the Translation different sheares are applied to additionally modify the locations of aliasing.
  4. Rand: $$$N/R$$$ positions are randomly chosen for each echo independently.
  5. RUSK: To ensure each k-space position is sampled equally, a random permutation of the $$$N$$$ k-space positions is repeated $$$J/N$$$ times, selecting $$$N/R$$$ subsequent positions in each echo.
  6. RHalton: To ensure quasi-uniform distribution of sampled k-space positions over all echoes, Halton sequence2 is used to generate the pseudo-random k-space positions with a low discrepancy.

Acquisition model

A ground truth T1 and T2 map shown in figure 2 is used, coil-sensitivity map is computed by the ESPIRIT technique using BART toolbox from actual FSE acquisition with 32-channel head coil.3,4 Similar to Barral et. al. 5, we use an Inversion Recovery prepared 3D Fast Spin Echo to quantify $$$M_0$$$, $$$T_1$$$ and $$$T_2$$$ with inversion times $$$[50, 400, 1100, 2500]$$$ ms, recovery time of $$$2608$$$ ms, echo train length of 18 with $$$6$$$ ms echo spacing and constant flip angle of $$$180^\circ$$$. We model that using $$$f_j(\theta_{\boldsymbol{x}})$$$ and undersampling in phase encoding directions.


Figure 3, 4 and 5 show, respectively, the results of $$$\eta_p$$$, StdNK and discrepancy. RHalton, with lower discrepancy, performs better than others even with higher StdNK. However, low StdNK is also relevant as observed for RandT at R = (32,36); even with higher discrepancy $$$\eta_p$$$ is better due to lower StdNK.


Most patterns show similar results for lower R but diverge for higher R, indicating sampling patterns are increasingly important for higher R. Comparing StdNK and discrepancy with $$$\eta_p$$$, we observe relevant correlations indicating high influence.


This study shows that low-discrepancy is most relevant for an efficient undersampling pattern. In the future, the use of CRLB will be extended to find the optimal value for other acquisition settings.


This project has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No 764513.


  1. S Heinrich. Efficient algorithms for computing the L2-discrepancy.Mathematics of com-putation, 65(216):1621–1633, October 1996.
  2. X. Wang and F.J. Hickernell. Randomized Halton sequences.Mathematical and Com-puter Modelling, 32(7-8):887–899, October 2000
  3. Martin Uecker, Peng Lai, Mark J. Murphy, Patrick Virtue, Michael Elad, John M. Pauly,Shreyas S. Vasanawala, and Michael Lustig. ESPIRiT-an eigenvalue approach to autocal-ibrating parallel MRI: Where SENSE meets GRAPPA.Magnetic Resonance in Medicine,71(3):990–1001, March 2014.
  4. Martin Uecker, Jonathan I. Tamir, Frank Ong, Siddharth Iyer, Joseph Y. Cheng, andMichael Lustig. Bart: Version 0.3.00, January 2016
  5. Jolle K. Barral, Erik Gudmundson, Nikola Stikov, Maryam Etezadi-Amoli, Petre Stoica,and Dwight G. Nishimura. A robust methodology for in vivo T1 mapping.MagneticResonance in Medicine, 64(4):1057–1067, June 2010.


Figure 1: Different patterns (a. Reg, b. Trans, c. TS, d. Rand, e. RUSK f. RHalton) for a k-space grid of 12 x12 and R = (4,2) with each circle representing a k-space position. Size of the circle represents number of times a particular position is sampled in a complete acquisition. Each circle is further divided into diff erent sectors of di fferent colours; each sector represents an echo in which the position was sampled. Colour of the sector distinguish diff erent echoes which are illustrated in the legend (g).

Figure 2: Ground truth T1 and T2 maps each box is 2x2 voxels with the values isillustrated by the legend.

Figure 3: Time efficiency for all the evaluated patterns for all R. Time efficiency values below 2 are not shown for clarity. Only Regular had values below 2

Figure 4: StdNK as a function of acceleration factor for all patterns. Pattern Regular is omitted for clarity due to its very high values.

Figure 5: Discrepancy of different patterns; for clarity, values higher than 0.025 are not shown.

Proc. Intl. Soc. Mag. Reson. Med. 27 (2019)