Compressed Sensing Plus Motion (CS+M): A New Perspective for Improving Undersampled MR Image Reconstruction

10/25/2018 ∙ by Angelica I. Aviles-Rivero, et al. ∙ University of Cambridge 0

Purpose: To obtain high-quality reconstructions from highly undersampled dynamic MRI data with the goal of reducing the acquisition time and towards improving physicians' outcome in clinical practice in a range of clinical applications. Theory and Methods: In dynamic MRI scans, the interaction between the target structure and the physical motion affects the acquired measurements. We exploit the strong repercussion of motion in MRI by proposing a variational framework - called Compressed Sensing Plus Motion (CS+M) - that links in a single model, simultaneously and explicitly, the computation of the algorithmic MRI reconstruction and the physical motion. Most precisely, we recast the image reconstruction and motion estimation problems as a single optimisation problem that is solved, iteratively, by breaking it up into two more computationally tractable problems. The potentials and generalisation capabilities of our approach are demonstrated in different clinical applications including cardiac cine, cardiac perfusion and brain perfusion imaging. Results: The proposed scheme reduces blurring artefacts and preserves the target shape and fine details whilst observing the lowest reconstruction error under highly undersampling up to 12x. This results in lower residual aliasing artefacts than the compared reconstructions algorithms. Overall, the results coming from our scheme exhibit more stable behaviour and generate a reconstruction closer to the gold-standard. Conclusion: We show that incorporating physical motion to the CS computation yields a significant improvement of the MR image reconstruction, that in fact, is closer to the gold-standard. This translates to higher reconstruction quality whilst requiring less measurements.

READ FULL TEXT VIEW PDF
POST COMMENT

Comments

There are no comments yet.

Authors

page 6

page 12

page 13

page 14

page 16

page 24

page 25

page 26

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

Abstract


Purpose: To obtain high-quality reconstructions from highly undersampled dynamic MRI data with the goal of reducing the acquisition time and towards improving physicians’ outcome in clinical practice in a range of clinical applications.


Theory and Methods: In dynamic MRI scans, the interaction between the target structure and the physical motion affects the acquired measurements. We exploit the strong repercussion of motion in MRI by proposing a variational framework - called Compressed Sensing Plus Motion (CS+M) - that links in a single model, simultaneously and explicitly, the computation of the algorithmic MRI reconstruction and the physical motion. Most precisely, we recast the image reconstruction and motion estimation problems as a single optimisation problem that is solved, iteratively, by breaking it up into two more computationally tractable problems. The potentials and generalisation capabilities of our approach are demonstrated in different clinical applications including cardiac cine, cardiac perfusion and brain perfusion imaging.


Results: The proposed scheme reduces blurring artefacts and preserves the target shape and fine details whilst observing the lowest reconstruction error under highly undersampling up to 12x. This results in lower residual aliasing artefacts than the compared reconstructions algorithms. Overall, the results coming from our scheme exhibit more stable behaviour and generate a reconstruction closer to the gold-standard.


Conclusion: We show that incorporating physical motion to the CS computation yields a significant improvement of the MR image reconstruction, that in fact, is closer to the gold-standard. This translates to higher reconstruction quality whilst requiring less measurements.


Keywords: Image Reconstruction; Compressed Sensing; Physical Motion Estimation; Variational Methods; Dynamic MRI.

1 Introduction

A central limitation of Magnetic Resonance Imaging (MRI) is the linear relation between the number of necessary measurements to form an image and the acquisition time. This constraint causes negative effects including [1, 2, 3]: (i) sensitivity to motion causing image degradation, (ii) reduced clinical throughput and (iii) patient non-compliance which introduces more artefacts during the image formation. Thus, a central challenge in MRI is to respond to the question How to decrease the long data-acquisition time?

There have been different attempts to answer this question. These attempts can rely on two categories: hardware-based or software-based approaches. The approaches from the first category aim to redesign the internal mechanisms of the scanner. However, despite of the continuous development of technologies, at the present, there is no evidence that the MRI data-acquisition time problem can be solved using a hardware-based approach [1, 4]. This is primarily due to both the physical limitations of gradient hardware, and the physiological limits imposed for safe pulse sequences.

Unlike hardware-based approaches, the second category has shown significant potential for effectively reducing acquisition times. This is particularly due to the advent of Compressed Sensing (CS) in MRI [5], which is motivated by the solid mathematical foundations of CS [6, 7]. The key idea of CS in MRI is to form an image represented in an appropriate transform basis from a considerably reduced finite-dimensional subset of k-space data acquired in an incoherent manner.

The potential and benefits of using CS in MRI have been demonstrated in different works starting with the pioneering work reported by Lustig et al. in [5], and then following by different approaches such as k-t SPARSE [8], k-t SPARSE-SENSE [9] and L1-SPIRiT [10], in which CS implications have been applied alone or in combination with parallel imaging.

A different direction based on low-rank matrix completion, which extends the idea of CS, has gained the attention of other scientists that have explored the effects of applying either local or global low-rank constraints. Liang in [11]

proposed to compute temporal basis functions, using singular value decomposition, as a new form to achieve MR image reconstruction from undersampled

k-

space data. This work set the basis and the motivation in the use of principal component analysis (PCA) for reconstructing a small subset of

k-space data (e.g. [12, 13, 14]). The idea of computing locally low-rank (LLR) constraint, in the context of undersampled MR image reconstruction, was reported in [15, 16, 17], in which spatiotemporal correlations were analysed in small regions. Although LLR has been suggested to decrease computational load [15], in comparison to global low-rank constraint, it comes at the expense of dealing with block artefacts [18].

The concept of combining sparse and low-rank constraints was also successfully reported by the MRI community. The authors in [19, 20, 21] aimed to find a solution that is both sparse and low-rank. A different approach based on decomposing the acquired measurements as a linear combination of low-rank (L) and sparse (S) components known as Robust Principal Component Analysis (RPCA) [22] or L+S decomposition was investigated in [23, 24, 25, 26, 27]. Particularly Otazo et al. demonstrated the feasibility of L+S model in [26] for various clinical applications in MRI including cardiac cine and abdominal perfusion. Even though the L+S model has proved to be effective for MR image reconstruction, the separation of background and dynamic components, in fact, is not always possible since the incoherence required by this model, i.e. for L and S, is not always fulfilled.

CONTRIBUTIONS. Notwithstanding that the body of literature has evidenced powerful results when using CS and low-rank methods for MR image reconstruction, in what follows we show that there still is significant room for improvement in terms of reconstruction quality in the context of spatio-temporal MRI reconstruction. Namely, in this work we introduce a mathematical framework - which follows the philosophy of [28]- to improve undersampled MR image reconstruction in a dynamic setting which we call Compressed Sensing plus Motion (CS+M). This approach was first discussed in our ISMRM abstract [29].

The key idea of our CS+M model is to link the computation in a single model of the MRI reconstruction to the complex motion patterns derived either from physiological or involuntary motion. More precisely, we recast the image formation problem as an optimisation problem that is solved, iteratively, by breaking it up into two more computationally tractable problems. Whilst this is an important part of our solution, our main contributions are:

  1. [label=0)]

  2. We introduce a computationally tractable variational model that allows establishing, explicitly and simultaneously, the connection between the MR image reconstruction and the inherent physical motion captured during the data acquisition. We show how the resulting optimisation problem can be solved efficiently by an alternating minimisation technique.

  3. We show that incorporating physical motion into the CS computation improves the MR image reconstruction from undersampled data, and produces images closer to fully sampled data whilst requiring less measurements. Furthermore, we show that image quality persists for higher acceleration factors than can be tolerated if motion is not included in the model.

  4. We provide evidence of the feasibility and generalisation capabilities of our model with several clinical applications including cardiac cine, cardiac perfusion and brain perfusion data.

2 Theory

In this section, we introduce our CS+M model for improving undersampled MRI reconstruction in dynamic settings. We first start describing the MR image formation problem and the current challenges in this context. We then introduce our CS+M model, which shall mitigate a major drawback of the existing solutions.

2.1 Undersampled MR Image Reconstruction

Consider a MRI setup, in which signal measurements are collected in spatial-frequency space (i.e. space) composed of time samples expressed as , where and denote a bounded image domain and temporal location respectively, and is the inherent noise during the acquisition. The matrix form of the measured data can be expressed as where is the target object to be reconstructed, refers to the space measurements and is the Fourier operator for a multiple receiver coil,

encodes coils sensitivities and the Fourier transform.

A current research direction in MRI is focused on reconstructing a small number of measurements with the aim of decreasing the acquisition time. This reconstruction task can be achieved by applying the notion of sparsity (i.e. CS implications), which is predominately approximated using the -norm. Formally, the problem of reconstructing undersampled MRI data can be cast as an optimisation problem, which reads:

(2.1)

where is the operator promoting sparse representation common transforms include a Wavelet Transform, temporal Fourier-Transformation (tFT), Total Variation (TV) and Total Generalized Variation (TGV) - e.g. [47, 43]. In what follows, we choose TV as sparsifying transformation. For tractability purposes, the algorithmic approach from (2.1) can be relaxed to the following regularised least squares functional:

(2.2)

where is the regularisation parameter controlling the importance of the two terms. Although the body of literature for reconstructing undersampled MRI data has provided promising results - for generalised CS approaches such the scheme described in (2.2) and extended ideas e.g. [19, 20, 26]- MRI reconstruction is still a challenging and open problem and there is plenty of room for further improvements. In particular, one seeks to cope with How to reconstruct high-quality MR images with less measurements? We respond to this question with our CS+M approach, which unlike existing approaches it considers, explicitly and simultaneously, the computation of the inherent complex motion patterns derived from ph ysiological or involuntary motion.

Figure 1: (From left-to-right) Reconstructed cine cardiac samples, selected region for displaying the motion estimation, and temporal evolution of the estimated physical motion.

2.2 CS+M: A New Perspective for MRI Reconstruction

In this section, we detail the proposed method and how it can be solved in a tractable computational manner. We first start describing how to compute the physical motion in a temporal scene. We then connect this computation, via a single model, to the notion of sparsity yielding to the CS+M model.

The key idea of the CS+M model is to link, in a single model, the algorithmic computation of the MRI reconstruction with the inherent physical motion in the scene. In the body of literature, the problem of computing two tasks simultaneously has been investigated in the community; examples are [30, 31], in which different constraints and assumptions of the scenes are applied. More recently, the philosophy of simultaneously computing the reconstruction of a given image sequence along with the motion estimation from a variational perspective has been proved in [28], which was later used in the medical domain for photoacoustic imaging [32] and X-ray tomography [33]. In this work, we follow the philosophy of [28] to build our CS+M model. An illustration of the temporal evolution of the estimated motion is displayed in Figure 1 along with a reconstructed sample. A key motivation for exploring the effects of incorporating motion to the MRI reconstruction comes with the response to the following question.

WHY IS MOTION IMPORTANT IN MRI? MRI is well-known to be highly sensitive to motion since its early development [34, 35]. In dynamic MRI scans, in particular, the interaction between the target structure and the physical motion affects the acquired measurements, meaning that, image degradation is produced which compromises clinical interpretations. Therefore, our hypothesis is that incorporating physical motion to a given algorithmic MRI reconstruction approach results in a faster acquisition and higher resolution in time, and hence, higher image quality.

We now turn to define the CS+M model. To do this, we address two questions: (i) How to compute the inherent physical motion in the scene? and (ii) How to plug in the algorithmic MRI reconstruction and the physical motion in a single model? Both questions are addressed next.

MODELLING PHYSICAL MOTION. Emphasising on our assumption of a dynamic setting, a key part of our solution is the explicit integration of an estimate of the scene’s motion to the MRI reconstruction. The problem of estimating physical motion from a set of images has been extensively investigated in the community (however it is still a challenging problem), in which solutions are roughly categorised into direct and indirect methods [36, 37]. In particular, Optical Flow is one of the most well-established methods, in which the predominant estimation scheme is based on variational approaches, starting with the pioneering work of Horn and Schunck [38].

Optical flow is based on the classical brightness constancy assumption, in which relative motion is approximated from the velocities of patterns’ brightness as where denotes a reconstructed image sequence and refers to time parameter. Linearised form of the constancy assumption yields to the optical flow constraint equation, in which the velocity is associated to the space-time image derivatives at location , expressed as:

(2.3)

where the path derivative is the velocity field, and is the gradient image. Using  (2.3) for recovering pointwise is not possible since there are two unknowns in one equation (termed aperture problem). However, one can deal with the underdetermined nature of  (2.3) by introducing further constraints such as a smoothness assumption on , and therefore, get a good approximation of the motion via a variational model [38]. Meaning that, one seeks to minimise an energy functional composed of two terms: a data fidelity term which models the brightness constancy deviations, and a regularisation term which penalises high variations in . Following this direction, we estimate the physical motion using the following energy functional:

(2.4)

where weights both terms. Our motivation to use the combination of Total Variation (TV) regularisation and

data fidelity terms is twofold: firstly, we aim to gain robustness in terms of outliers, this is an important factor since in a MRI setting inherent noise during the acquisition is assumed. Secondly, this combination allows for discontinuities in the flow field. Overall, the TV-

optical flow yields to a robust estimation of the physical motion (e.g. [39, 40, 41]).

CS+M MODEL WHEN CS MEETS MOTION. We now turn to describe how (2.4) is incorporating into the algorithmic MRI reconstruction. Generally speaking, this problem can be cast as an optimisation problem, in which one seeks to optimise, simultaneously, over the MRI image formation and the estimated motion as:

(2.5)

where are a positive weighting parameters, and we consider the choice of the sparsifying transformation as TV (see Section 3.3 for discussion). More precisely, plugging-in the algorithmic MRI reconstruction approach and the physical motion defined in  (2.2) and  (2.4) respectively into (2.5), the CS+M model seeks to solve the following optimisation problem:

(2.6)

For tractability purposes, we define a time-discrete version of (2.6) by performing the discretisation of the time-interval into -steps, this shall result in reconstructed images and vector fields. Therefore, the scheme from (2.6) is now defined as:

(2.7)

1

Algorithm 1 CS+M Algorithm
  Start from and ;
  Set ;
  while   do
     ;
     ;
     Optimise over , for the current value of :
           solve(primalDualOf(2.8)) using [42];
      solve(absdiff(, ));
     Optimise over , for the current value of :
           solve(primalDualOf(2.9)) using [42];
      solve(absdiff(, ));
     ;
  end while
  return .

From a computational point of view, the model from (2.7) presents some numerical challenges it is non-convex, and non-differentiable (due to the norm). Therefore, to achieve more computational tractability we solve (2.6) by alternating minimisation, that is, we break up the (2.6) into two more tractable sub-problems. Following [28], we alternate between optimising over and whilst keeping the other fixed. Details on the solution of each sub-problem, at each iteration, is as follows.

Sub-problem 1: Optimisation over . The optimisation problem over while keeping fixed is defined as:

(2.8)

Sub-problem 2: Optimisation over . For a fixed , the optimisation problem in is given by:

(2.9)

the splitting of (2.6) allows now for a tractable and straightforward computation, in which approximate solutions can be obtained by applying a primal-dual algorithm [42] which, generally speaking, considers the saddle-point version of (2.8) and (2.9) in the form . The overall procedure is listed in Algorithm 1.

3 Methods

This section describes in detail the experimentations that we conducted to validate our CS+M reconstruction technique.

3.1 Data Description

We evaluated our approach extensively and prove its generalisation by using data coming from: cine cardiac, cardiac perfusion and brain perfusion MR imaging. The datasets were saved as fully sampled raw data and then were retrospectively undersampled using a variable-density random sampling pattern as suggested by Lustig in [5] and using cartesian sampling. Datasets characteristics are as follows:

  • Dataset I A cine cardiac dataset which was acquired from a healthy volunteer, from [43]. Measurements were collected using a 3T Siemens scanner. Matrix size, heart phases, coils, FOV = mm and TA=16s.

  • Datasets II and II Realistic cardiac cine simulation generated using the MRXCAT phantom framework [44]. Whilst the Dataset II was simulated during breath-holding, the III was set with free respiratory motion. The duration of the heart beat and respiratory cycles were set as 1 sec and 5 sec respectively. Both were generated with Matrix size, heart phases, coils.

  • Datasets IV and V Realistic cardiac perfusion simulations using MRXCAT[44] framework. Simulations were set with Matrix size, Frames , coils. The difference, for analyses purposes, relies on the fact that the Dataset IV was simulated with breath-holding whilst the dataset V during free-breathing with the respiratory cycle set as 5 sec.

  • Dataset VI A cardiac perfusion dataset acquired during free-breathing. Available from the computational biomedical imaging group at the University of Iowa 111https://research.engineering.uiowa.edu/cbig/ . It was acquired using FLASH sequence on a 3T Siemens scanner (TR/TE=2.5/1.5ms) with a matrix size of .

  • Dataset VII A single-coil brain perfusion dataset acquired from multi slice 2-D dynamic susceptibility contrast (DSC) from with frame separation of TR=2s, and with a matrix size of .

All the measurements and reconstructions in this section were taken from these datasets. All tests and comparisons were run under the same conditions in a CPU-based implementation. Sensitivity estimation along with computational details are discussed in Section 3 of the supplementary material.

Figure 2: Sample frames of the clinical datasets used for evaluating our approach from several clinical applications (from left-to-right): cine cardiac, cardiac perfusion and brain perfusion.

3.2 Evaluation Scheme

We validate our approach based on both qualitative and quantitative analyses. Whilst the qualitative analysis is derived from the observations and interpretations of physicians, for the quantitative analysis, we rely on computing a set of metrics between the gold-standard and the reconstructed MR image. From now, we refer to the reconstruction obtained from computing the sum-of-squares (SoS) using fully sampled data as gold-standard. More precisely, we base our findings on the following evaluation scheme.

  1. [label=(

    We address the quantitative analysis relying on two well-established image-quality metrics: (i) the Structural Similarity (SSIM) Index [45], which calculates the similarity of two image reconstructions from the contrast, luminance and structure, and (ii) the inverted Localised Mean Squared Error (sLMSE) [46] that computes the local similarity based on local patches. The computation of the sLMSE is a normalised and inverted measure such that closer to 1 means higher quality reconstruction. The practicalities of the metrics’ computation are discussed in Section 4 of the supplementary material.

    3.3 Parameter Selection

    Generally speaking, parameters of our approach and the ones we compare with were set individually chosen from a range of values with respect to the best SSIM and sLMSE, for each type of data. In particular for the case of the L+S approach, we define the range of values as suggested in [26] along with the code provided on the authors’ website 222http://cai2r.net/resources/software.

    For our approach, there are three parameters to set . These parameters were tested in a range of values and for each clinical application type as follows. For the parameter , we found the best outcome by setting it to whilst for the other two parameters, , we tested them in the range . It is to be noticed that the parameters of the CS+M model can be easily fine-tuned for a particular clinical application and kept fixed for other datasets with similar dynamic information. Moreover, we used TV [47] as an operator for promoting sparse representation in our and compared schemes. The objective of this work is to have a proof of concept to show the potentials and generalisation capabilities of our approach for MRI reconstruction. Therefore, a detailed regularisers investigation shall be tackled in future work.

    3.4 Results and Discussion

    We describe our findings following the scheme described in Section 3.2. We start by evaluating classic CS MRI scheme against our CS+M model. In Figure 3, we display three reconstructed samples with an acceleration factor of 8x, and for three different datasets. Visual assessment of the reconstructed samples agrees with our initial hypothesis that the incorporation of motion in the reconstruction model benefits the output quality. Most notably, it can be observed in a comparison with the gold-standard, the right and left ventricular endocardial borders (see outlined green and red regions), that CS+M offers better reconstructions in terms of contrast and shape than the CS reconstructions (see yellow arrows). Moreover, in a closer inspection at the zoomed-in views, it is to be noticed the loss of fine details and blurring effect at the papillary muscles.

    This is further reflected in the signal intensity profile at right side of Figure 3, where, and for all the displayed cases, the CS+M approach (red line) is closer to the gold-standard (blue line). Although CS based reconstruction (green line) offers a good approximation to the gold-standard, it fails to eliminate all perturbations such as blurring artefacts (see yellow arrows), which our approach is able to prevent. This is reflected in the behaviour of signal intensity profile in which significant oscillations are observed; as is visible on the zoomed-in views at right side of Figure 3. In particular, it can be observed at the signal intensity profile of Dataset III which was acquired during free-breathing (last row of Figure 3) that signals generated from CS described strong oscillations yielding to an unstable behaviour (see black squares).

    For a more detailed analyses, we evaluate our approach by comparing its performance against three different reconstruction schemes: zero-filling, CS and L+S. Figure 4 displays reconstructed samples, taken from datasets I and II and undersampled at 8x, of the chosen schemes and our proposed one along with the gold-standard. By visual evaluation, we observed that the compared schemes tends to produce blurring artefacts and loss of fine details (see red arrows at (a),(b) of Figure 4). This is further reflected and better appreciated in the displayed SSIM maps at (a.1),(b.2) from the same figure, which offer meaningful comparison of local image quality over space. A closer inspection shows that our approach produces a reconstruction with higher similarity metric to the gold-standard. These results are compatible with the signal intensity profiles plotted in Figure 4-(a.2),(b.2), in which the compared approaches exhibit a fluctuation behaviour compared to the gold-standard whilst our approach (red line) is closer to and agrees better with the gold-standard (blue line). We also quantified the temporal LV blood pool area (green region) which is displayed in Figure 4-(a.3), (b.3), in which our approach show a better LV function compared to the gold-standard.

    Next, we investigate - how both CS+M and the compared schemes perform in the temporal domain?

    Figure 3: Displayed reconstruction of three samples. (From left-to-right and top-to-bottom) Reconstructed full-sampled data (i.e. ground truth), reconstructed samples using CS and our CS+M. Zoom-in views in the middle part show more details in which the left (red) and right (green) ventricular endocardial borders have been outlined. Yellow arrows show artefacts created during the CS reconstruction. Signal intensity profiles retrieved from the reconstructed sample, blue lines in the plots show ideal signal intensity green and red lines represent the CS and CS+M approaches. The results show that measurements generated from our approach are closer to the gold-standard.
    Figure 4: Comparison performance of our approach and three reconstruction schemes along with the ground truth. (a) and (b) reconstructed samples from datasets I and II in which loss of detail and contrast, and introduction of blurring artefacts are pointed out with the red arrows. (a.1) and (b.1) corresponding SSIM maps, in which values as closer to 1 as better reconstruction quality. (a.2) and (b.2) numerical visualisation of the signal intensity profile. (a.3) and (b.3) temporal LV blood pool area.
    Figure 5: Comparison performance on the temporal domain of our approach against three from the body of literature. (a) and (b) display a reconstructed sample for each evaluated scheme along with the gold-standard using Datasets V and VI. (a.1) and (b.1) show the corresponding temporal signal intensity profiles and error plots of three regions of interest at myocardium.
    4x 6x 8x
    Exp. Dataset
    Reconstruction
    Scheme
    sLMSE SSIM sLMSE SSIM sLMSE SSIM
    1
    Dataset I
    Cardiac Cine
    Zero-Filling 84.65 82.15 83.56 79.76 83.34 76.51
    CS 98.81 90.90 98.73 89.56 96.13 85.71
    L+S 99.34 92.07 98.99 89.67 97.09 87.21
    CS+M 99.87 94.62 99.71 92.40 98.30 89.90
    2
    Dataset II
    Breath-holding
    Cardiac Cine
    Zero-Filling 83.85 84.69 83.05 81.70 82.91 78.16
    CS 98.91 91.41 97.74 88.07 96.69 84.28
    L+S 98.86 92.33 98.34 89.68 97.19 87.55
    CS+M 99.87 93.36 99.70 91.07 98.54 89.81
    3
    Dataset III
    Free-breathing
    Cardiac Cine
    Zero-Filling 84.88 80.10 84.46 76.51 83.86 73.74
    CS 98.84 90.07 97.56 86.69 95.82 83.14
    L+S 99.32 91.85 98.41 88.02 97.03 85.60
    CS+M 99.81 93.09 99.56 90.87 98.70 88.14
    4
    Dataset IV
    Breath-holding
    Myocardial Perfusion
    Zero-Filling 85.63 82.80 85.34 80.47 84.65 78.07
    CS 98.95 91.53 98.89 88.62 97.74 84.15
    L+S 99.50 92.45 99.21 89.13 98.86 86.11
    CS+M 99.81 92.95 99.86 91.52 99.08 88.67
    5
    Dataset V
    Free-breathing
    Myocardial Perfusion
    Zero-Filling 85.51 82.11 84.73 80.45 83.64 76.37
    CS 98.70 89.75 96.22 84.45 94.88 81.11
    L+S 99.10 91.37 98.41 87.22 96.88 83.67
    CS+M 99.72 92.05 99.15 89.16 98.30 86.25
    6
    Dataset VI
    Free-breathing
    Myocardial Perfusion
    Zero-Filling 84.72 82.57 84.50 79.57 83.29 75.55
    CS 98.03 89.99 96.85 86.80 93.98 82.81
    L+S 98.26 91.75 97.05 88.03 95.61 85.54
    CS+M 98.81 93.96 98.09 90.01 97.45 88.68
    7
    Dataset VII
    Brain Perfusion
    Zero-Filling 83.85 81.77 83.05 78.03 82.91 76.32
    CS 98.02 92.91 96.86 89.46 94.60 85.59
    L+S 98.57 93.57 97.88 90.90 96.52 86.27
    CS+M 99.33 95.04 98.21 92.75 97.89 88.36
    Table 1: Global performance analyses from our approach and three other reconstruction schemes. Displayed measures averaged over all the corresponding dataset denoted in .

    We executed another experiment using the whole cardiac perfusion datasets, we selected this application since the dynamic information contained in it differs from the cardiac cine information. In Figure 5, we show four reconstructed samples resulted from the compared schemes and our approach along with the gold-standard, (A) at 7.5x and (B) at 8x acceleration. By visual inspection, we observed that the effects discussed above from the cardiac cine datasets also prevail for the perfusion cardiac case, meaning that, the compared approaches reflect perturbations in the reconstructions such as blurring effects, and loss of contrast and details. To do this, we computed the reconstructions over the whole corresponding datasets, and extracted the temporal intensity for three regions of interest in the myocardium. The results are displayed at Figure 5 - (a.1),(b.1), a closer inspection show higher temporal fluctuations which translate in higher residual artefacts for the compared schemes, whilst our approach reflected a more stable signal which is in fact closer to the gold-standard in all cases. This is further reflected by the error plots (Root-mean-square error RMSE), in which our approach had the lowest error value for all reconstructed samples. From the compared schemes the one that performs better is L+S, however, we observed on the temporal intensity signals and the resulted reconstructions (Figures 45), contrast variations and oscillations in the signal intensity. By contrast, our approach displayed more stable signals intensities with less fluctuations and blurring artefacts, and better preservation of the shape. An accompanying comparison at the LV cavity is provided in Section 2 of the supplementary material.

    But Is there a significant difference in reconstruction quality between our approach and the compared reconstruction schemes? To respond to this question and to further support Figures 4 and 5, we computed the reconstruction of all datasets under high undersampling factors up to 12x. The plot of Figure 6 shows the SSIM and sLMSE curves of both approaches where we can observe that the CS+M outperforms CS and L+S at all acceleration counts. For example, with an acceleration of 8x, the reconstruction quality obtained with our approach can only be achieved with an acceleration of 6x for CS for both metrics. Meaning that, despite reducing the measurement samples, our approach is still able to generate higher-quality images. A similar behaviour, but with a less numerical difference, is observed in a comparison between our approach and L+S.

    Figure 6: Plots comparing CS and CS+M using two reconstruction quality metrics, SSIM and sLMSE, over the whole cardiac perfusion datasets.
    Figure 7: (From left-to-right) Reconstructed samples from the DSC brain perfusion- dataset VII. (a.1) Temporal signal intensity profiles, in which significant fluctuations are produced by the CS approach along with (a.2) the error plot.

    For quantitative evaluation of generalisation capabilities of our approach and for a global analysis performance, we report the results in Table 1. The reported numbers are the average of the selected image metrics across the entire corresponding dataset. It is to be noticed, and as we previously mentioned, both metrics as closer to 1 as higher reconstruction quality. From the results, we observe that our approach outperforms the compared reconstruction schemes with respect to both the sLMSE and SSIM metrics. We also noticed that in the case where there are different types of dynamics - for example the organ’s physiological motion along with the breathing - the reconstruction benefits from this knowledge reflecting in a significant improvement of the reconstruction (for example see Exps. 2-3, and Exps. 4-5).

    To visualise the motion incorporated to the image formation, we show in Figure 7 a reconstructed sample and the temporal evolution of an estimate of the physical motion. Right side of Figure 7 displays the temporal signal intensity profiles, in which a closer inspection shows that CS reconstructions (green line) exhibit significant oscillations which our approach (red line) is able to prevent. This is further reflected in terms of image quality at Exps. 7 and 8 of Table 1, in which our approach produced higher values than CS for both quality metrics

    Overall, a global inspection of the qualitative and quantitative analyses from the compared schemes shows the following drawbacks, on which our CS+M approach improves:

    • Introduction of Temporal Artefacts. We notice that the compared schemes exhibit significant oscillations in the temporal signal intensity, which is translated in higher residual aliasing artefacts. This effect is prevailed for all datasets. By contrast, the CS+M approach displayed temporal stable signals which are closer to the gold-standard yielding to higher reconstructions quality.

    • Fine Details and Shape Preservation. From the results, we observe that compared approaches exhibit loss in details and blurring effect, for example the papillary muscles of the heart or in the brain, and even more significantly noticed during complex dynamic transitions such as the contraction-expansion of the heart. Conversely, our approach offered a better spatio-temporal fidelity for all applications.

    • Stability under Physiological Motion. We notice that the compared schemes tend to have significant fluctuations during changes in the dynamic information such as systolic phases (in the case of the heart) or when addition dynamic information appears as in the case of free-breathing datasets. During these events, the results coming from our scheme exhibit a more stable behaviour resulting in a reconstruction closer to the gold-standard.

    These points comes to highlight the benefit of incorporating motion into the algorithmic MRI reconstruction whilst supporting our initial hypothesis.

    4 Conclusion, Limitations and Future Work

    In this paper we addressed a central question in MRI which is How to get high quality MRI reconstructed images under highly undersampling factors? We respond to this question through a new approach called CS+M, in which the novelty largely relies on incorporating, explicitly and simultaneously in a single model, an estimate of the scene’s motion to the algorithmic MRI reconstruction to provide higher quality images with less motion artefacts. We demonstrate the potentials of our approach based on exhaustive qualitative and quantitative analyses, in which we show that our approach outperforms traditional reconstruction schemes in terms of better preservation of fine details and organs’ shape and less blurring artefacts, resulting in less residual artefacts and reconstructions closer to the gold-standard.

    Whilst the objective of this work is to open a new line of research for further clinical investigation, future work will include experimentation with Non-Cartesian sampling and computational speed improvement based on a CUDA-based implementation. Moreover, since the CS+M model proved that motion has significant positive effects that translates to clinical potentials, future work might address how to improve the motion estimation.

    Acknowledgments

    AIAR acknowledges support from the Centre for Mathematical Imaging in Healthcare (CMIH), University of Cambridge is greatly acknowledged. We would like to thank Samar Alsaleh from the George Washington University for her discussion. CBS acknowledges support from Leverhulme Trust project on Breaking the non-convexity barrier, EPSRC grant Nr. EP/M00483X/1, the EPSRC Centre Nr. EP/N014588/1, the RISE projects CHiPS and NoMADS, the Cantab Capital Institute for the Mathematics of Information and the Alan Turing Institute.

    References

    • [1] Zaitsev M, Maclaren J and Herbst, M. Motion artifacts in MRI: a complex problem with many partial solutions. Journal of Magnetic Resonance Imaging, 2015; 42(4):887–901.
    • [2] Hollingsworth KG. Reducing acquisition time in clinical MRI by data undersampling and compressed sensing reconstruction. Physics in medicine and biology, 2015; vol. 60, no 21.
    • [3] Havsteen I, Ohlhues A, Madsen KH, Nybing JD, Christensen H and Christensen A. Are Movement Artifacts in Magnetic Resonance Imaging a Real Problem? — A Narrative Review. Frontiers in neurology, 2017; vol. 8.
    • [4] Majumdar A. Compressed sensing for magnetic resonance image reconstruction. Cambridge University Press, 2015; ch. 2, p. 50.
    • [5] Lustig M, Donoho D and Pauly JM. Sparse MRI: The application of compressed sensing for rapid MR Imaging. Magnetic Resonance in Medicine, 2007; 58:1182–1195.
    • [6] Candes EJ, Romberg JK and Tao T. Stable signal recovery from incomplete and inaccurate measurements. Communications on pure and applied mathematics, 2006; vol. 59, no 8.
    • [7] Donoho DL. Compressed sensing. IEEE Transactions on Information Theory, 2006; vol. 52, no 4.
    • [8] Lustig M, Santos JM, Donoho DL, and Pauly JM. k-t SPARSE: High frame rate dynamic MRI exploiting spatio-temporal sparsity. In Proceedings of the 13th Annual Meeting of ISMRM, 2006; vol. 2420.
    • [9] Otazo R, Kim D, Axel L and Sodickson DK. Combination of compressed sensing and parallel imaging for highly accelerated first‐pass cardiac perfusion MRI. Magnetic Resonance in Medicine, 2010; 64(3), 767-776.
    • [10] Lustig M, and Pauly JM. SPIRiT: Iterative self‐consistent parallel imaging reconstruction from arbitrary k‐space. Magnetic resonance in medicine, 2010; 64(2), 457-471.
    • [11] Liang ZP. Spatiotemporal imagingwith partially separable functions. IEEE International Symposium on Biomedical Imaging, 2007.
    • [12] Pedersen H, Kozerke S, Ringgaard S, Nehrke K and Kim WY. k-t PCA: Temporally constrained k-t BLAST reconstruction using principal component analysis. Magnetic resonance in medicine, 2009; 62(3), 706-716.
    • [13] Feng L, Srichai MB, Lim RP, Harrison A, King W., Adluru G., Dibella EVR, Sodickson DK, Otazo R, Kim D. Highly accelerated real‐time cardiac cine MRI using k-t SPARSE-SENSE. Magnetic resonance in medicine, 2013; 70(1), 64-74.
    • [14] Velikina JV and Samsonov AA. Reconstruction of dynamic image series from undersampled MRI data using data‐driven model consistency condition (MOCCO). Magnetic resonance in medicine, 2015; 74(5), 1279-1290.
    • [15] Trzasko J, Manduca A and Borisch E. Local versus global low-rank promotion in dynamic MRI series reconstruction. In Proc. Int. Symp. Magn. Reson. Med, 2011; p. 4371.
    • [16] Zhang T, Pauly JM and Levesque IR. Accelerating parameter mapping with a locally low rank constraint. Magnetic resonance in medicine, 2015; 73(2), 655-661.
    • [17] Miao X, Lingala SG, Guo Y, Jao T, Usman M, Prieto C and Nayak KS. Accelerated cardiac cine MRI using locally low rank and finite difference constraints. Magnetic resonance imaging, 2016; 34(6), 707-714.
    • [18] Saucedo A, Lefkimmiatis S, Rangwala N and Sung K. Improved Computational Efficiency of Locally Low Rank MRI Reconstruction Using Iterative Random Patch Adjustments. IEEE Transactions on Medical Imaging, 2017.
    • [19] Lingala SG, Hu Y, DiBella E and Jacob M. Accelerated dynamic MRI exploiting sparsity and low-rank structure: kt SLR. IEEE transactions on medical imaging, 2011; 30(5), 1042-1054.
    • [20] Majumdar A and Ward RK. Exploiting rank deficiency and transform domain sparsity for MR image reconstruction. Magnetic resonance imaging, 2012; 30(1), 9-18.
    • [21] Zhao B, Haldar JP, Christodoulou AG and Liang, ZP. Image reconstruction from highly undersampled (k,t)-space data with joint partial separability and sparsity constraints. IEEE transactions on medical imaging, 2012; 31(9), 1809-1820.
    • [22] Candès EJ, Li X, Ma Y and Wright J. Robust principal component analysis?. Journal of the ACM (JACM), (2011); 58(3), 11.
    • [23] Gao H, Rapacchi S, Wang D, Moriarty J, Meehan C, Sayre J, Laub G, Finn P and Hu P. Compressed sensing using prior rank, intensity and sparsity model (PRISM): applications in cardiac cine MRI. In Proceedings of the 20th Annual Meeting of ISMRM, Melbourne, Australia, 2012.
    • [24] Trémoulhéac B, Dikaios N, Atkinson D and Arridge SR. Dynamic MR Image Reconstruction–Separation From Undersampled ()-Space via Low-Rank Plus Sparse Prior. IEEE Transactions on Medical Imaging, 2014; 33(8), 1689–1701.
    • [25] Otazo R, Emmanuel C and Sodickson DK. Low-rank and sparse matrix decomposition for accelerated DCE-MRI with background and contrast separation. In Proceedings of the ISMRM Workshop on Data Sampling and Image Reconstruction, 2013.
    • [26] Otazo R, Candes E and Sodickson DK. Low-rank plus sparse matrix decomposition for accelerated dynamic MRI with separation of background and dynamic components. Magnetic Resonance in Medicine, 2015; 73(3), 1125–1136.
    • [27] Singh V, Tewfik AH and Ress DB. Under-sampled functional MRI using low-rank plus sparse matrix decomposition. IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2015; pp. 897-901.
    • [28] Burger M, Dirks H and Schönlieb C-B. A variational model for joint motion estimation and image reconstruction, SIAM Journal on Imaging Sciences, 2018; pp. 94–128.
    • [29] Aviles-Rivero AI, Williams G, Graves M and Schönlieb CB. CS+M: A Simultaneous Reconstruction and Motion Estimation Approach for Improving Undersampled MRI Reconstruction. In Proceedings of the 26th Annual Meeting ISMRM, 2018.
    • [30]

      Tomasi C and Kanade T. Shape and motion from image streams under orthography: a factorization method. International Journal of Computer Vision (IJCV), 1992; 9(2), pp. 137-154.

    • [31] Gilland DR, Mair BA, Bowsher JE and Jaszczak RJ. Simultaneous reconstruction and motion estimation for gated cardiac ECT. IEEE Transactions on Nuclear Science, 2002; 49(5), pp. 2344-2349.
    • [32] Huynh N, Lucka F., Zhang E, Betcke M, Arridge S, Beard P and Cox B. Sub-sampled Fabry-Perot photoacoustic scanner for fast 3D imaging, Proc.SPIE, 2017; vol.10064.
    • [33] Burger M, Dirks H , Frerking L, Hauptmann A, Helin T and Siltanen S. A variational reconstruction method for undersampled dynamic X-ray tomography based on physical motion models, Inverse Problems, 2017; pp. 1–24.
    • [34] Wood M, and Henkelman RM. MR image artifacts from periodic motion. Medical physics, 1985; 12.2: 143-151.
    • [35] Van de Walle R, Lemahieu I and Achten E. Magnetic resonance imaging and the reduction of motion artifacts: review of the principles. Technology and Health Care, 1997; 5(6), 419-435.
    • [36] Torr PH and Zisserman A. Feature based methods for structure and motion estimation. The International Conference on Computer Vision - Workshop on vision algorithms (ICCVW), 1999; pp. 278–294.
    • [37] Irani M and Anandan P. About direct methods. The International Conference on Computer Vision - Workshop on vision algorithms (ICCVW), 1999; pp. 267–277.
    • [38]

      Horn BK and Schunck BG. Determining optical flow. Artificial intelligence, 1981; 17(1-3), 185-203.

    • [39] Aubert G, Deriche R and Kornprobst P. Computing optical flow via variational techniques. SIAM Journal on Applied Mathematics, 1999; 60(1), 156-182.
    • [40]

      Zach C, Pock T and Bischof H. A duality based approach for realtime TV-L1 optical flow. In Joint Pattern Recognition Symposium Springer, 2007; pp. 214–223.

    • [41] Perez JS, Meinhardt-Llopis E and Facciolo G. TV-L1 optical flow estimation. Image Processing On Line (IPOL), 2013; pp. 137–150.
    • [42] Chambolle A and Pock T. A first-order primal-dual algorithm for convex problems with applications to imaging. Journal of Mathematical Imaging and Vision, 2011; vol. 40 pp. 120–145.
    • [43] Schloegl M, Holler M, Schwarzl A, Bredies K and Stollberger R. Infimal convolution of total generalized variation functionals for dynamic MRI. Magnetic resonance in medicine, 2017; 78(1), 142-155.
    • [44] Wissmann L, Santelli C, Segars WP and Kozerke S. MRXCAT: Realistic numerical phantoms for cardiovascular magnetic resonance. Journal of Cardiovascular Magnetic Resonance, 2014; 16(1), 63.
    • [45] Wang Z, Bovik AC, Sheikh HR and Simoncelli EP. Image quality assessment: from error visibility to structural similarity. IEEE Transactions on Image Processing (TIP), 2004; 13(4), pp. 600-612.
    • [46] Grosse R, Johnson MK, Adelson EH and Freeman WT. Ground truth dataset and baseline evaluations for intrinsic image algorithms. IEEE International Conference in Computer Vision (ICCV), 2009; pp. 2335-2342.
    • [47] Rudin LI, Osher S and Fatemi E. Nonlinear total variation based noise removal algorithms. Journal of Physics D, 1992; 60:259–268.
    • [48] Schlögl M, Holler M, Bredies K and Stollberger R. A variational approach for coil-sensitivity estimation for undersampled phase-sensitive dynamic MRI reconstruction. In Proceedings of the 23th Annual Meeting ISMRM, 2015.

    Supplemental Material for the Article:

    Compressed Sensing Plus Motion (CS+M): A New Perspective for Improving Undersampled MR Image Reconstruction.

    Angelica I. Aviles-Rivero, Guy Williams, Martin Graves and Carola-Bibiane Schönlieb.


    Centre for Mathematical Imaging in Healthcare, Department of Pure Mathematics and Mathematical Statistics (DPMMS), University of Cambridge, U.K.

    Wolfson Brain Imaging Centre, Department of Clinical Neurosciences, University of Cambridge, U.K.

    Department of Radiology, University of Cambridge, U.K.

    Department for Applied Mathematics and Theoretical Physics (DAMTP), University of Cambridge, U.K.


    1 Outline

    This document extends the clinical and technical discussion presented in the main paper, this, with the aim of offering further details regarding our experiments and results. The remaining document is structured as follows.

    • [noitemsep]

    • Section 2: More qualitative results for a closer inspection of the reconstructed data.

    • Section 3: We provide more details about the experimental settings and the quantitative results.

    • Section 4: We give an explicit definition and motivation of the metrics used for the quantitative analyses.

    2 Supplementary Visual Results

    To further support the quantitative and qualitative results presented in the main paper, in this section we offer visual comparison of our approach against the L+S [25] method, which had the second best approximation as reported in Table 1 in the main paper.

    Figure 8 shows the results of both methods along with the gold-standard for visual quality assessment. For this comparison, we used the cardiac datasets 3 and 4, cine and perfusion correspondingly, and selected some clinically meaningful frames. Whilst the L+S technique frequently produces a good reconstruction, a closer inspection of 8 shows that our approach performs better in terms of preserving both the anatomical shape and details. Moreover, we can observe less blurry effects from our approach (pointed out with yellow arrows in the zoom-in views). The results in Figure 11 echoes similar findings for the application of cardiac perfusion, in which our approach generates a more visually meaningful reconstructions. This is more clearly reflected with the red arrows in the zoom-in views.

    Figure 8: Visual comparisons of our and the L+S approaches using Dataset III. We display of reconstructed samples from different cardiac cycles, in which yellow arrows point at regions of interest.
    Figure 9: Reconstructed samples of cardiac perfusion (Dataset IV), in which visual comparisons between our and L+S approaches, along with the gold-standard, is shown. Zoom-in views indicate interesting regions for comparison.
    Figure 10: Comparison performance on the temporal domain of our approach against three from the body of literature. (a) and (b) display a reconstructed sample for each evaluated scheme along with the gold-standard using Datasets V and VI. (a.1) and (b.1) show the corresponding temporal signal intensity profiles of a region of interest (blue square) in the LV cavity.
    Figure 11: Reconstructed sample from the brain perfusion dataset along with few visualisation plots of the physical motion used to improve the image reconstruction.

    3 Further Details and Analysis of the Results

    For clarification purposes of the experimental results, in this section, we discuss details regarding the multi-coil reconstructions and the elapsed time of our approach.

    Figure 12: Elapsed time comparison between L+S and ours approaches.

    In this work, we used both single- and multi- coill MRI data (see Section 3.1 in the main paper for detailed description of the datasets). Whilst for the multi-coil data generated using the MRXCAT framework, the coil sensitivities were computed using the Biot-Savart law (refer to [44] for further details), for the non-simulated data, the coil sensitivities were estimated using [48].

    Computational time requirements. We ran a set of experiments to demonstrate the computational feasibility of our reconstruction scheme. In particular, we compared the elapsed time between the L+S and our approach on a CPU basis and using all datasets described in Section 3.1 in the main paper. The experiments were designed to measure only the elapsed time required by the optimisation task. We found that our approach required an average of more computational time than that of [26]. See Figure  12 for three visualisations of the elapsed time with different datasets sizes. This is a slight increment in time that comes with the benefits of having a higher image quality even during sever undersampling factors. Moreover, our implementation can be significantly improved on a CUDA-based implementation.

    4 Image-Quality Metrics.

    For clarification purposes and with the aim of avoiding some ambiguity regarding the form of how the two quality metrics were computed. We describe the explicit definition of the used metrics, which computation is as follows:

    Structural Similarity Index (SSIM). Let and be the reconstructions coming from the gold-standard and a given approximation from a reconstruction scheme, and let

    be the local mean, variances and covariance of

    and respectively. Then the SSIM is computed as:

    where and are constants to ensure numerical stability when the denominator is close to zero. In our experiments, we set standard values used for as and .

    Inverted Localised Mean Squared Error (sLMSE). It is to be noted that we are not using the most common metrics for comparing pixel-wise two reconstructions, which are Peak signal-to-noise ratio (PSNR) and MSE. Our main motivation for using the sLMSE is that PSNR and MSE are too stringent (based on a global absolute difference, where a small incorrectly reconstructed edge can dominate the error). Given and - a gold-standard and an approximation reconstruction - as input, we compute the LMSE as the MSE summed over patches, , of size and spaced in steps of 10 as:

    to be consistent with the SSIM, we compute a nomalised and inverted version so that the maximum output value is 1, meaning that, as s closer to 1 as higher quality reconstruction which reads:

4 Conclusion, Limitations and Future Work

In this paper we addressed a central question in MRI which is How to get high quality MRI reconstructed images under highly undersampling factors? We respond to this question through a new approach called CS+M, in which the novelty largely relies on incorporating, explicitly and simultaneously in a single model, an estimate of the scene’s motion to the algorithmic MRI reconstruction to provide higher quality images with less motion artefacts. We demonstrate the potentials of our approach based on exhaustive qualitative and quantitative analyses, in which we show that our approach outperforms traditional reconstruction schemes in terms of better preservation of fine details and organs’ shape and less blurring artefacts, resulting in less residual artefacts and reconstructions closer to the gold-standard.

Whilst the objective of this work is to open a new line of research for further clinical investigation, future work will include experimentation with Non-Cartesian sampling and computational speed improvement based on a CUDA-based implementation. Moreover, since the CS+M model proved that motion has significant positive effects that translates to clinical potentials, future work might address how to improve the motion estimation.

Acknowledgments

AIAR acknowledges support from the Centre for Mathematical Imaging in Healthcare (CMIH), University of Cambridge is greatly acknowledged. We would like to thank Samar Alsaleh from the George Washington University for her discussion. CBS acknowledges support from Leverhulme Trust project on Breaking the non-convexity barrier, EPSRC grant Nr. EP/M00483X/1, the EPSRC Centre Nr. EP/N014588/1, the RISE projects CHiPS and NoMADS, the Cantab Capital Institute for the Mathematics of Information and the Alan Turing Institute.

References

  • [1] Zaitsev M, Maclaren J and Herbst, M. Motion artifacts in MRI: a complex problem with many partial solutions. Journal of Magnetic Resonance Imaging, 2015; 42(4):887–901.
  • [2] Hollingsworth KG. Reducing acquisition time in clinical MRI by data undersampling and compressed sensing reconstruction. Physics in medicine and biology, 2015; vol. 60, no 21.
  • [3] Havsteen I, Ohlhues A, Madsen KH, Nybing JD, Christensen H and Christensen A. Are Movement Artifacts in Magnetic Resonance Imaging a Real Problem? — A Narrative Review. Frontiers in neurology, 2017; vol. 8.
  • [4] Majumdar A. Compressed sensing for magnetic resonance image reconstruction. Cambridge University Press, 2015; ch. 2, p. 50.
  • [5] Lustig M, Donoho D and Pauly JM. Sparse MRI: The application of compressed sensing for rapid MR Imaging. Magnetic Resonance in Medicine, 2007; 58:1182–1195.
  • [6] Candes EJ, Romberg JK and Tao T. Stable signal recovery from incomplete and inaccurate measurements. Communications on pure and applied mathematics, 2006; vol. 59, no 8.
  • [7] Donoho DL. Compressed sensing. IEEE Transactions on Information Theory, 2006; vol. 52, no 4.
  • [8] Lustig M, Santos JM, Donoho DL, and Pauly JM. k-t SPARSE: High frame rate dynamic MRI exploiting spatio-temporal sparsity. In Proceedings of the 13th Annual Meeting of ISMRM, 2006; vol. 2420.
  • [9] Otazo R, Kim D, Axel L and Sodickson DK. Combination of compressed sensing and parallel imaging for highly accelerated first‐pass cardiac perfusion MRI. Magnetic Resonance in Medicine, 2010; 64(3), 767-776.
  • [10] Lustig M, and Pauly JM. SPIRiT: Iterative self‐consistent parallel imaging reconstruction from arbitrary k‐space. Magnetic resonance in medicine, 2010; 64(2), 457-471.
  • [11] Liang ZP. Spatiotemporal imagingwith partially separable functions. IEEE International Symposium on Biomedical Imaging, 2007.
  • [12] Pedersen H, Kozerke S, Ringgaard S, Nehrke K and Kim WY. k-t PCA: Temporally constrained k-t BLAST reconstruction using principal component analysis. Magnetic resonance in medicine, 2009; 62(3), 706-716.
  • [13] Feng L, Srichai MB, Lim RP, Harrison A, King W., Adluru G., Dibella EVR, Sodickson DK, Otazo R, Kim D. Highly accelerated real‐time cardiac cine MRI using k-t SPARSE-SENSE. Magnetic resonance in medicine, 2013; 70(1), 64-74.
  • [14] Velikina JV and Samsonov AA. Reconstruction of dynamic image series from undersampled MRI data using data‐driven model consistency condition (MOCCO). Magnetic resonance in medicine, 2015; 74(5), 1279-1290.
  • [15] Trzasko J, Manduca A and Borisch E. Local versus global low-rank promotion in dynamic MRI series reconstruction. In Proc. Int. Symp. Magn. Reson. Med, 2011; p. 4371.
  • [16] Zhang T, Pauly JM and Levesque IR. Accelerating parameter mapping with a locally low rank constraint. Magnetic resonance in medicine, 2015; 73(2), 655-661.
  • [17] Miao X, Lingala SG, Guo Y, Jao T, Usman M, Prieto C and Nayak KS. Accelerated cardiac cine MRI using locally low rank and finite difference constraints. Magnetic resonance imaging, 2016; 34(6), 707-714.
  • [18] Saucedo A, Lefkimmiatis S, Rangwala N and Sung K. Improved Computational Efficiency of Locally Low Rank MRI Reconstruction Using Iterative Random Patch Adjustments. IEEE Transactions on Medical Imaging, 2017.
  • [19] Lingala SG, Hu Y, DiBella E and Jacob M. Accelerated dynamic MRI exploiting sparsity and low-rank structure: kt SLR. IEEE transactions on medical imaging, 2011; 30(5), 1042-1054.
  • [20] Majumdar A and Ward RK. Exploiting rank deficiency and transform domain sparsity for MR image reconstruction. Magnetic resonance imaging, 2012; 30(1), 9-18.
  • [21] Zhao B, Haldar JP, Christodoulou AG and Liang, ZP. Image reconstruction from highly undersampled (k,t)-space data with joint partial separability and sparsity constraints. IEEE transactions on medical imaging, 2012; 31(9), 1809-1820.
  • [22] Candès EJ, Li X, Ma Y and Wright J. Robust principal component analysis?. Journal of the ACM (JACM), (2011); 58(3), 11.
  • [23] Gao H, Rapacchi S, Wang D, Moriarty J, Meehan C, Sayre J, Laub G, Finn P and Hu P. Compressed sensing using prior rank, intensity and sparsity model (PRISM): applications in cardiac cine MRI. In Proceedings of the 20th Annual Meeting of ISMRM, Melbourne, Australia, 2012.
  • [24] Trémoulhéac B, Dikaios N, Atkinson D and Arridge SR. Dynamic MR Image Reconstruction–Separation From Undersampled ()-Space via Low-Rank Plus Sparse Prior. IEEE Transactions on Medical Imaging, 2014; 33(8), 1689–1701.
  • [25] Otazo R, Emmanuel C and Sodickson DK. Low-rank and sparse matrix decomposition for accelerated DCE-MRI with background and contrast separation. In Proceedings of the ISMRM Workshop on Data Sampling and Image Reconstruction, 2013.
  • [26] Otazo R, Candes E and Sodickson DK. Low-rank plus sparse matrix decomposition for accelerated dynamic MRI with separation of background and dynamic components. Magnetic Resonance in Medicine, 2015; 73(3), 1125–1136.
  • [27] Singh V, Tewfik AH and Ress DB. Under-sampled functional MRI using low-rank plus sparse matrix decomposition. IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2015; pp. 897-901.
  • [28] Burger M, Dirks H and Schönlieb C-B. A variational model for joint motion estimation and image reconstruction, SIAM Journal on Imaging Sciences, 2018; pp. 94–128.
  • [29] Aviles-Rivero AI, Williams G, Graves M and Schönlieb CB. CS+M: A Simultaneous Reconstruction and Motion Estimation Approach for Improving Undersampled MRI Reconstruction. In Proceedings of the 26th Annual Meeting ISMRM, 2018.
  • [30]

    Tomasi C and Kanade T. Shape and motion from image streams under orthography: a factorization method. International Journal of Computer Vision (IJCV), 1992; 9(2), pp. 137-154.

  • [31] Gilland DR, Mair BA, Bowsher JE and Jaszczak RJ. Simultaneous reconstruction and motion estimation for gated cardiac ECT. IEEE Transactions on Nuclear Science, 2002; 49(5), pp. 2344-2349.
  • [32] Huynh N, Lucka F., Zhang E, Betcke M, Arridge S, Beard P and Cox B. Sub-sampled Fabry-Perot photoacoustic scanner for fast 3D imaging, Proc.SPIE, 2017; vol.10064.
  • [33] Burger M, Dirks H , Frerking L, Hauptmann A, Helin T and Siltanen S. A variational reconstruction method for undersampled dynamic X-ray tomography based on physical motion models, Inverse Problems, 2017; pp. 1–24.
  • [34] Wood M, and Henkelman RM. MR image artifacts from periodic motion. Medical physics, 1985; 12.2: 143-151.
  • [35] Van de Walle R, Lemahieu I and Achten E. Magnetic resonance imaging and the reduction of motion artifacts: review of the principles. Technology and Health Care, 1997; 5(6), 419-435.
  • [36] Torr PH and Zisserman A. Feature based methods for structure and motion estimation. The International Conference on Computer Vision - Workshop on vision algorithms (ICCVW), 1999; pp. 278–294.
  • [37] Irani M and Anandan P. About direct methods. The International Conference on Computer Vision - Workshop on vision algorithms (ICCVW), 1999; pp. 267–277.
  • [38]

    Horn BK and Schunck BG. Determining optical flow. Artificial intelligence, 1981; 17(1-3), 185-203.

  • [39] Aubert G, Deriche R and Kornprobst P. Computing optical flow via variational techniques. SIAM Journal on Applied Mathematics, 1999; 60(1), 156-182.
  • [40]

    Zach C, Pock T and Bischof H. A duality based approach for realtime TV-L1 optical flow. In Joint Pattern Recognition Symposium Springer, 2007; pp. 214–223.

  • [41] Perez JS, Meinhardt-Llopis E and Facciolo G. TV-L1 optical flow estimation. Image Processing On Line (IPOL), 2013; pp. 137–150.
  • [42] Chambolle A and Pock T. A first-order primal-dual algorithm for convex problems with applications to imaging. Journal of Mathematical Imaging and Vision, 2011; vol. 40 pp. 120–145.
  • [43] Schloegl M, Holler M, Schwarzl A, Bredies K and Stollberger R. Infimal convolution of total generalized variation functionals for dynamic MRI. Magnetic resonance in medicine, 2017; 78(1), 142-155.
  • [44] Wissmann L, Santelli C, Segars WP and Kozerke S. MRXCAT: Realistic numerical phantoms for cardiovascular magnetic resonance. Journal of Cardiovascular Magnetic Resonance, 2014; 16(1), 63.
  • [45] Wang Z, Bovik AC, Sheikh HR and Simoncelli EP. Image quality assessment: from error visibility to structural similarity. IEEE Transactions on Image Processing (TIP), 2004; 13(4), pp. 600-612.
  • [46] Grosse R, Johnson MK, Adelson EH and Freeman WT. Ground truth dataset and baseline evaluations for intrinsic image algorithms. IEEE International Conference in Computer Vision (ICCV), 2009; pp. 2335-2342.
  • [47] Rudin LI, Osher S and Fatemi E. Nonlinear total variation based noise removal algorithms. Journal of Physics D, 1992; 60:259–268.
  • [48] Schlögl M, Holler M, Bredies K and Stollberger R. A variational approach for coil-sensitivity estimation for undersampled phase-sensitive dynamic MRI reconstruction. In Proceedings of the 23th Annual Meeting ISMRM, 2015.

Supplemental Material for the Article:

Compressed Sensing Plus Motion (CS+M): A New Perspective for Improving Undersampled MR Image Reconstruction.

Angelica I. Aviles-Rivero, Guy Williams, Martin Graves and Carola-Bibiane Schönlieb.


Centre for Mathematical Imaging in Healthcare, Department of Pure Mathematics and Mathematical Statistics (DPMMS), University of Cambridge, U.K.

Wolfson Brain Imaging Centre, Department of Clinical Neurosciences, University of Cambridge, U.K.

Department of Radiology, University of Cambridge, U.K.

Department for Applied Mathematics and Theoretical Physics (DAMTP), University of Cambridge, U.K.


1 Outline

This document extends the clinical and technical discussion presented in the main paper, this, with the aim of offering further details regarding our experiments and results. The remaining document is structured as follows.

  • [noitemsep]

  • Section 2: More qualitative results for a closer inspection of the reconstructed data.

  • Section 3: We provide more details about the experimental settings and the quantitative results.

  • Section 4: We give an explicit definition and motivation of the metrics used for the quantitative analyses.

2 Supplementary Visual Results

To further support the quantitative and qualitative results presented in the main paper, in this section we offer visual comparison of our approach against the L+S [25] method, which had the second best approximation as reported in Table 1 in the main paper.

Figure 8 shows the results of both methods along with the gold-standard for visual quality assessment. For this comparison, we used the cardiac datasets 3 and 4, cine and perfusion correspondingly, and selected some clinically meaningful frames. Whilst the L+S technique frequently produces a good reconstruction, a closer inspection of 8 shows that our approach performs better in terms of preserving both the anatomical shape and details. Moreover, we can observe less blurry effects from our approach (pointed out with yellow arrows in the zoom-in views). The results in Figure 11 echoes similar findings for the application of cardiac perfusion, in which our approach generates a more visually meaningful reconstructions. This is more clearly reflected with the red arrows in the zoom-in views.

Figure 8: Visual comparisons of our and the L+S approaches using Dataset III. We display of reconstructed samples from different cardiac cycles, in which yellow arrows point at regions of interest.
Figure 9: Reconstructed samples of cardiac perfusion (Dataset IV), in which visual comparisons between our and L+S approaches, along with the gold-standard, is shown. Zoom-in views indicate interesting regions for comparison.
Figure 10: Comparison performance on the temporal domain of our approach against three from the body of literature. (a) and (b) display a reconstructed sample for each evaluated scheme along with the gold-standard using Datasets V and VI. (a.1) and (b.1) show the corresponding temporal signal intensity profiles of a region of interest (blue square) in the LV cavity.
Figure 11: Reconstructed sample from the brain perfusion dataset along with few visualisation plots of the physical motion used to improve the image reconstruction.

3 Further Details and Analysis of the Results

For clarification purposes of the experimental results, in this section, we discuss details regarding the multi-coil reconstructions and the elapsed time of our approach.

Figure 12: Elapsed time comparison between L+S and ours approaches.

In this work, we used both single- and multi- coill MRI data (see Section 3.1 in the main paper for detailed description of the datasets). Whilst for the multi-coil data generated using the MRXCAT framework, the coil sensitivities were computed using the Biot-Savart law (refer to [44] for further details), for the non-simulated data, the coil sensitivities were estimated using [48].

Computational time requirements. We ran a set of experiments to demonstrate the computational feasibility of our reconstruction scheme. In particular, we compared the elapsed time between the L+S and our approach on a CPU basis and using all datasets described in Section 3.1 in the main paper. The experiments were designed to measure only the elapsed time required by the optimisation task. We found that our approach required an average of more computational time than that of [26]. See Figure  12 for three visualisations of the elapsed time with different datasets sizes. This is a slight increment in time that comes with the benefits of having a higher image quality even during sever undersampling factors. Moreover, our implementation can be significantly improved on a CUDA-based implementation.

4 Image-Quality Metrics.

For clarification purposes and with the aim of avoiding some ambiguity regarding the form of how the two quality metrics were computed. We describe the explicit definition of the used metrics, which computation is as follows:

Structural Similarity Index (SSIM). Let and be the reconstructions coming from the gold-standard and a given approximation from a reconstruction scheme, and let

be the local mean, variances and covariance of

and respectively. Then the SSIM is computed as:

where and are constants to ensure numerical stability when the denominator is close to zero. In our experiments, we set standard values used for as and .

Inverted Localised Mean Squared Error (sLMSE). It is to be noted that we are not using the most common metrics for comparing pixel-wise two reconstructions, which are Peak signal-to-noise ratio (PSNR) and MSE. Our main motivation for using the sLMSE is that PSNR and MSE are too stringent (based on a global absolute difference, where a small incorrectly reconstructed edge can dominate the error). Given and - a gold-standard and an approximation reconstruction - as input, we compute the LMSE as the MSE summed over patches, , of size and spaced in steps of 10 as:

to be consistent with the SSIM, we compute a nomalised and inverted version so that the maximum output value is 1, meaning that, as s closer to 1 as higher quality reconstruction which reads:

Acknowledgments

AIAR acknowledges support from the Centre for Mathematical Imaging in Healthcare (CMIH), University of Cambridge is greatly acknowledged. We would like to thank Samar Alsaleh from the George Washington University for her discussion. CBS acknowledges support from Leverhulme Trust project on Breaking the non-convexity barrier, EPSRC grant Nr. EP/M00483X/1, the EPSRC Centre Nr. EP/N014588/1, the RISE projects CHiPS and NoMADS, the Cantab Capital Institute for the Mathematics of Information and the Alan Turing Institute.

References

  • [1] Zaitsev M, Maclaren J and Herbst, M. Motion artifacts in MRI: a complex problem with many partial solutions. Journal of Magnetic Resonance Imaging, 2015; 42(4):887–901.
  • [2] Hollingsworth KG. Reducing acquisition time in clinical MRI by data undersampling and compressed sensing reconstruction. Physics in medicine and biology, 2015; vol. 60, no 21.
  • [3] Havsteen I, Ohlhues A, Madsen KH, Nybing JD, Christensen H and Christensen A. Are Movement Artifacts in Magnetic Resonance Imaging a Real Problem? — A Narrative Review. Frontiers in neurology, 2017; vol. 8.
  • [4] Majumdar A. Compressed sensing for magnetic resonance image reconstruction. Cambridge University Press, 2015; ch. 2, p. 50.
  • [5] Lustig M, Donoho D and Pauly JM. Sparse MRI: The application of compressed sensing for rapid MR Imaging. Magnetic Resonance in Medicine, 2007; 58:1182–1195.
  • [6] Candes EJ, Romberg JK and Tao T. Stable signal recovery from incomplete and inaccurate measurements. Communications on pure and applied mathematics, 2006; vol. 59, no 8.
  • [7] Donoho DL. Compressed sensing. IEEE Transactions on Information Theory, 2006; vol. 52, no 4.
  • [8] Lustig M, Santos JM, Donoho DL, and Pauly JM. k-t SPARSE: High frame rate dynamic MRI exploiting spatio-temporal sparsity. In Proceedings of the 13th Annual Meeting of ISMRM, 2006; vol. 2420.
  • [9] Otazo R, Kim D, Axel L and Sodickson DK. Combination of compressed sensing and parallel imaging for highly accelerated first‐pass cardiac perfusion MRI. Magnetic Resonance in Medicine, 2010; 64(3), 767-776.
  • [10] Lustig M, and Pauly JM. SPIRiT: Iterative self‐consistent parallel imaging reconstruction from arbitrary k‐space. Magnetic resonance in medicine, 2010; 64(2), 457-471.
  • [11] Liang ZP. Spatiotemporal imagingwith partially separable functions. IEEE International Symposium on Biomedical Imaging, 2007.
  • [12] Pedersen H, Kozerke S, Ringgaard S, Nehrke K and Kim WY. k-t PCA: Temporally constrained k-t BLAST reconstruction using principal component analysis. Magnetic resonance in medicine, 2009; 62(3), 706-716.
  • [13] Feng L, Srichai MB, Lim RP, Harrison A, King W., Adluru G., Dibella EVR, Sodickson DK, Otazo R, Kim D. Highly accelerated real‐time cardiac cine MRI using k-t SPARSE-SENSE. Magnetic resonance in medicine, 2013; 70(1), 64-74.
  • [14] Velikina JV and Samsonov AA. Reconstruction of dynamic image series from undersampled MRI data using data‐driven model consistency condition (MOCCO). Magnetic resonance in medicine, 2015; 74(5), 1279-1290.
  • [15] Trzasko J, Manduca A and Borisch E. Local versus global low-rank promotion in dynamic MRI series reconstruction. In Proc. Int. Symp. Magn. Reson. Med, 2011; p. 4371.
  • [16] Zhang T, Pauly JM and Levesque IR. Accelerating parameter mapping with a locally low rank constraint. Magnetic resonance in medicine, 2015; 73(2), 655-661.
  • [17] Miao X, Lingala SG, Guo Y, Jao T, Usman M, Prieto C and Nayak KS. Accelerated cardiac cine MRI using locally low rank and finite difference constraints. Magnetic resonance imaging, 2016; 34(6), 707-714.
  • [18] Saucedo A, Lefkimmiatis S, Rangwala N and Sung K. Improved Computational Efficiency of Locally Low Rank MRI Reconstruction Using Iterative Random Patch Adjustments. IEEE Transactions on Medical Imaging, 2017.
  • [19] Lingala SG, Hu Y, DiBella E and Jacob M. Accelerated dynamic MRI exploiting sparsity and low-rank structure: kt SLR. IEEE transactions on medical imaging, 2011; 30(5), 1042-1054.
  • [20] Majumdar A and Ward RK. Exploiting rank deficiency and transform domain sparsity for MR image reconstruction. Magnetic resonance imaging, 2012; 30(1), 9-18.
  • [21] Zhao B, Haldar JP, Christodoulou AG and Liang, ZP. Image reconstruction from highly undersampled (k,t)-space data with joint partial separability and sparsity constraints. IEEE transactions on medical imaging, 2012; 31(9), 1809-1820.
  • [22] Candès EJ, Li X, Ma Y and Wright J. Robust principal component analysis?. Journal of the ACM (JACM), (2011); 58(3), 11.
  • [23] Gao H, Rapacchi S, Wang D, Moriarty J, Meehan C, Sayre J, Laub G, Finn P and Hu P. Compressed sensing using prior rank, intensity and sparsity model (PRISM): applications in cardiac cine MRI. In Proceedings of the 20th Annual Meeting of ISMRM, Melbourne, Australia, 2012.
  • [24] Trémoulhéac B, Dikaios N, Atkinson D and Arridge SR. Dynamic MR Image Reconstruction–Separation From Undersampled ()-Space via Low-Rank Plus Sparse Prior. IEEE Transactions on Medical Imaging, 2014; 33(8), 1689–1701.
  • [25] Otazo R, Emmanuel C and Sodickson DK. Low-rank and sparse matrix decomposition for accelerated DCE-MRI with background and contrast separation. In Proceedings of the ISMRM Workshop on Data Sampling and Image Reconstruction, 2013.
  • [26] Otazo R, Candes E and Sodickson DK. Low-rank plus sparse matrix decomposition for accelerated dynamic MRI with separation of background and dynamic components. Magnetic Resonance in Medicine, 2015; 73(3), 1125–1136.
  • [27] Singh V, Tewfik AH and Ress DB. Under-sampled functional MRI using low-rank plus sparse matrix decomposition. IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2015; pp. 897-901.
  • [28] Burger M, Dirks H and Schönlieb C-B. A variational model for joint motion estimation and image reconstruction, SIAM Journal on Imaging Sciences, 2018; pp. 94–128.
  • [29] Aviles-Rivero AI, Williams G, Graves M and Schönlieb CB. CS+M: A Simultaneous Reconstruction and Motion Estimation Approach for Improving Undersampled MRI Reconstruction. In Proceedings of the 26th Annual Meeting ISMRM, 2018.
  • [30]

    Tomasi C and Kanade T. Shape and motion from image streams under orthography: a factorization method. International Journal of Computer Vision (IJCV), 1992; 9(2), pp. 137-154.

  • [31] Gilland DR, Mair BA, Bowsher JE and Jaszczak RJ. Simultaneous reconstruction and motion estimation for gated cardiac ECT. IEEE Transactions on Nuclear Science, 2002; 49(5), pp. 2344-2349.
  • [32] Huynh N, Lucka F., Zhang E, Betcke M, Arridge S, Beard P and Cox B. Sub-sampled Fabry-Perot photoacoustic scanner for fast 3D imaging, Proc.SPIE, 2017; vol.10064.
  • [33] Burger M, Dirks H , Frerking L, Hauptmann A, Helin T and Siltanen S. A variational reconstruction method for undersampled dynamic X-ray tomography based on physical motion models, Inverse Problems, 2017; pp. 1–24.
  • [34] Wood M, and Henkelman RM. MR image artifacts from periodic motion. Medical physics, 1985; 12.2: 143-151.
  • [35] Van de Walle R, Lemahieu I and Achten E. Magnetic resonance imaging and the reduction of motion artifacts: review of the principles. Technology and Health Care, 1997; 5(6), 419-435.
  • [36] Torr PH and Zisserman A. Feature based methods for structure and motion estimation. The International Conference on Computer Vision - Workshop on vision algorithms (ICCVW), 1999; pp. 278–294.
  • [37] Irani M and Anandan P. About direct methods. The International Conference on Computer Vision - Workshop on vision algorithms (ICCVW), 1999; pp. 267–277.
  • [38]

    Horn BK and Schunck BG. Determining optical flow. Artificial intelligence, 1981; 17(1-3), 185-203.

  • [39] Aubert G, Deriche R and Kornprobst P. Computing optical flow via variational techniques. SIAM Journal on Applied Mathematics, 1999; 60(1), 156-182.
  • [40]

    Zach C, Pock T and Bischof H. A duality based approach for realtime TV-L1 optical flow. In Joint Pattern Recognition Symposium Springer, 2007; pp. 214–223.

  • [41] Perez JS, Meinhardt-Llopis E and Facciolo G. TV-L1 optical flow estimation. Image Processing On Line (IPOL), 2013; pp. 137–150.
  • [42] Chambolle A and Pock T. A first-order primal-dual algorithm for convex problems with applications to imaging. Journal of Mathematical Imaging and Vision, 2011; vol. 40 pp. 120–145.
  • [43] Schloegl M, Holler M, Schwarzl A, Bredies K and Stollberger R. Infimal convolution of total generalized variation functionals for dynamic MRI. Magnetic resonance in medicine, 2017; 78(1), 142-155.
  • [44] Wissmann L, Santelli C, Segars WP and Kozerke S. MRXCAT: Realistic numerical phantoms for cardiovascular magnetic resonance. Journal of Cardiovascular Magnetic Resonance, 2014; 16(1), 63.
  • [45] Wang Z, Bovik AC, Sheikh HR and Simoncelli EP. Image quality assessment: from error visibility to structural similarity. IEEE Transactions on Image Processing (TIP), 2004; 13(4), pp. 600-612.
  • [46] Grosse R, Johnson MK, Adelson EH and Freeman WT. Ground truth dataset and baseline evaluations for intrinsic image algorithms. IEEE International Conference in Computer Vision (ICCV), 2009; pp. 2335-2342.
  • [47] Rudin LI, Osher S and Fatemi E. Nonlinear total variation based noise removal algorithms. Journal of Physics D, 1992; 60:259–268.
  • [48] Schlögl M, Holler M, Bredies K and Stollberger R. A variational approach for coil-sensitivity estimation for undersampled phase-sensitive dynamic MRI reconstruction. In Proceedings of the 23th Annual Meeting ISMRM, 2015.

Supplemental Material for the Article:

Compressed Sensing Plus Motion (CS+M): A New Perspective for Improving Undersampled MR Image Reconstruction.

Angelica I. Aviles-Rivero, Guy Williams, Martin Graves and Carola-Bibiane Schönlieb.


Centre for Mathematical Imaging in Healthcare, Department of Pure Mathematics and Mathematical Statistics (DPMMS), University of Cambridge, U.K.

Wolfson Brain Imaging Centre, Department of Clinical Neurosciences, University of Cambridge, U.K.

Department of Radiology, University of Cambridge, U.K.

Department for Applied Mathematics and Theoretical Physics (DAMTP), University of Cambridge, U.K.


1 Outline

This document extends the clinical and technical discussion presented in the main paper, this, with the aim of offering further details regarding our experiments and results. The remaining document is structured as follows.

  • [noitemsep]

  • Section 2: More qualitative results for a closer inspection of the reconstructed data.

  • Section 3: We provide more details about the experimental settings and the quantitative results.

  • Section 4: We give an explicit definition and motivation of the metrics used for the quantitative analyses.

2 Supplementary Visual Results

To further support the quantitative and qualitative results presented in the main paper, in this section we offer visual comparison of our approach against the L+S [25] method, which had the second best approximation as reported in Table 1 in the main paper.

Figure 8 shows the results of both methods along with the gold-standard for visual quality assessment. For this comparison, we used the cardiac datasets 3 and 4, cine and perfusion correspondingly, and selected some clinically meaningful frames. Whilst the L+S technique frequently produces a good reconstruction, a closer inspection of 8 shows that our approach performs better in terms of preserving both the anatomical shape and details. Moreover, we can observe less blurry effects from our approach (pointed out with yellow arrows in the zoom-in views). The results in Figure 11 echoes similar findings for the application of cardiac perfusion, in which our approach generates a more visually meaningful reconstructions. This is more clearly reflected with the red arrows in the zoom-in views.

Figure 8: Visual comparisons of our and the L+S approaches using Dataset III. We display of reconstructed samples from different cardiac cycles, in which yellow arrows point at regions of interest.
Figure 9: Reconstructed samples of cardiac perfusion (Dataset IV), in which visual comparisons between our and L+S approaches, along with the gold-standard, is shown. Zoom-in views indicate interesting regions for comparison.
Figure 10: Comparison performance on the temporal domain of our approach against three from the body of literature. (a) and (b) display a reconstructed sample for each evaluated scheme along with the gold-standard using Datasets V and VI. (a.1) and (b.1) show the corresponding temporal signal intensity profiles of a region of interest (blue square) in the LV cavity.
Figure 11: Reconstructed sample from the brain perfusion dataset along with few visualisation plots of the physical motion used to improve the image reconstruction.

3 Further Details and Analysis of the Results

For clarification purposes of the experimental results, in this section, we discuss details regarding the multi-coil reconstructions and the elapsed time of our approach.

Figure 12: Elapsed time comparison between L+S and ours approaches.

In this work, we used both single- and multi- coill MRI data (see Section 3.1 in the main paper for detailed description of the datasets). Whilst for the multi-coil data generated using the MRXCAT framework, the coil sensitivities were computed using the Biot-Savart law (refer to [44] for further details), for the non-simulated data, the coil sensitivities were estimated using [48].

Computational time requirements. We ran a set of experiments to demonstrate the computational feasibility of our reconstruction scheme. In particular, we compared the elapsed time between the L+S and our approach on a CPU basis and using all datasets described in Section 3.1 in the main paper. The experiments were designed to measure only the elapsed time required by the optimisation task. We found that our approach required an average of more computational time than that of [26]. See Figure  12 for three visualisations of the elapsed time with different datasets sizes. This is a slight increment in time that comes with the benefits of having a higher image quality even during sever undersampling factors. Moreover, our implementation can be significantly improved on a CUDA-based implementation.

4 Image-Quality Metrics.

For clarification purposes and with the aim of avoiding some ambiguity regarding the form of how the two quality metrics were computed. We describe the explicit definition of the used metrics, which computation is as follows:

Structural Similarity Index (SSIM). Let and be the reconstructions coming from the gold-standard and a given approximation from a reconstruction scheme, and let

be the local mean, variances and covariance of

and respectively. Then the SSIM is computed as:

where and are constants to ensure numerical stability when the denominator is close to zero. In our experiments, we set standard values used for as and .

Inverted Localised Mean Squared Error (sLMSE). It is to be noted that we are not using the most common metrics for comparing pixel-wise two reconstructions, which are Peak signal-to-noise ratio (PSNR) and MSE. Our main motivation for using the sLMSE is that PSNR and MSE are too stringent (based on a global absolute difference, where a small incorrectly reconstructed edge can dominate the error). Given and - a gold-standard and an approximation reconstruction - as input, we compute the LMSE as the MSE summed over patches, , of size and spaced in steps of 10 as:

to be consistent with the SSIM, we compute a nomalised and inverted version so that the maximum output value is 1, meaning that, as s closer to 1 as higher quality reconstruction which reads:

References

  • [1] Zaitsev M, Maclaren J and Herbst, M. Motion artifacts in MRI: a complex problem with many partial solutions. Journal of Magnetic Resonance Imaging, 2015; 42(4):887–901.
  • [2] Hollingsworth KG. Reducing acquisition time in clinical MRI by data undersampling and compressed sensing reconstruction. Physics in medicine and biology, 2015; vol. 60, no 21.
  • [3] Havsteen I, Ohlhues A, Madsen KH, Nybing JD, Christensen H and Christensen A. Are Movement Artifacts in Magnetic Resonance Imaging a Real Problem? — A Narrative Review. Frontiers in neurology, 2017; vol. 8.
  • [4] Majumdar A. Compressed sensing for magnetic resonance image reconstruction. Cambridge University Press, 2015; ch. 2, p. 50.
  • [5] Lustig M, Donoho D and Pauly JM. Sparse MRI: The application of compressed sensing for rapid MR Imaging. Magnetic Resonance in Medicine, 2007; 58:1182–1195.
  • [6] Candes EJ, Romberg JK and Tao T. Stable signal recovery from incomplete and inaccurate measurements. Communications on pure and applied mathematics, 2006; vol. 59, no 8.
  • [7] Donoho DL. Compressed sensing. IEEE Transactions on Information Theory, 2006; vol. 52, no 4.
  • [8] Lustig M, Santos JM, Donoho DL, and Pauly JM. k-t SPARSE: High frame rate dynamic MRI exploiting spatio-temporal sparsity. In Proceedings of the 13th Annual Meeting of ISMRM, 2006; vol. 2420.
  • [9] Otazo R, Kim D, Axel L and Sodickson DK. Combination of compressed sensing and parallel imaging for highly accelerated first‐pass cardiac perfusion MRI. Magnetic Resonance in Medicine, 2010; 64(3), 767-776.
  • [10] Lustig M, and Pauly JM. SPIRiT: Iterative self‐consistent parallel imaging reconstruction from arbitrary k‐space. Magnetic resonance in medicine, 2010; 64(2), 457-471.
  • [11] Liang ZP. Spatiotemporal imagingwith partially separable functions. IEEE International Symposium on Biomedical Imaging, 2007.
  • [12] Pedersen H, Kozerke S, Ringgaard S, Nehrke K and Kim WY. k-t PCA: Temporally constrained k-t BLAST reconstruction using principal component analysis. Magnetic resonance in medicine, 2009; 62(3), 706-716.
  • [13] Feng L, Srichai MB, Lim RP, Harrison A, King W., Adluru G., Dibella EVR, Sodickson DK, Otazo R, Kim D. Highly accelerated real‐time cardiac cine MRI using k-t SPARSE-SENSE. Magnetic resonance in medicine, 2013; 70(1), 64-74.
  • [14] Velikina JV and Samsonov AA. Reconstruction of dynamic image series from undersampled MRI data using data‐driven model consistency condition (MOCCO). Magnetic resonance in medicine, 2015; 74(5), 1279-1290.
  • [15] Trzasko J, Manduca A and Borisch E. Local versus global low-rank promotion in dynamic MRI series reconstruction. In Proc. Int. Symp. Magn. Reson. Med, 2011; p. 4371.
  • [16] Zhang T, Pauly JM and Levesque IR. Accelerating parameter mapping with a locally low rank constraint. Magnetic resonance in medicine, 2015; 73(2), 655-661.
  • [17] Miao X, Lingala SG, Guo Y, Jao T, Usman M, Prieto C and Nayak KS. Accelerated cardiac cine MRI using locally low rank and finite difference constraints. Magnetic resonance imaging, 2016; 34(6), 707-714.
  • [18] Saucedo A, Lefkimmiatis S, Rangwala N and Sung K. Improved Computational Efficiency of Locally Low Rank MRI Reconstruction Using Iterative Random Patch Adjustments. IEEE Transactions on Medical Imaging, 2017.
  • [19] Lingala SG, Hu Y, DiBella E and Jacob M. Accelerated dynamic MRI exploiting sparsity and low-rank structure: kt SLR. IEEE transactions on medical imaging, 2011; 30(5), 1042-1054.
  • [20] Majumdar A and Ward RK. Exploiting rank deficiency and transform domain sparsity for MR image reconstruction. Magnetic resonance imaging, 2012; 30(1), 9-18.
  • [21] Zhao B, Haldar JP, Christodoulou AG and Liang, ZP. Image reconstruction from highly undersampled (k,t)-space data with joint partial separability and sparsity constraints. IEEE transactions on medical imaging, 2012; 31(9), 1809-1820.
  • [22] Candès EJ, Li X, Ma Y and Wright J. Robust principal component analysis?. Journal of the ACM (JACM), (2011); 58(3), 11.
  • [23] Gao H, Rapacchi S, Wang D, Moriarty J, Meehan C, Sayre J, Laub G, Finn P and Hu P. Compressed sensing using prior rank, intensity and sparsity model (PRISM): applications in cardiac cine MRI. In Proceedings of the 20th Annual Meeting of ISMRM, Melbourne, Australia, 2012.
  • [24] Trémoulhéac B, Dikaios N, Atkinson D and Arridge SR. Dynamic MR Image Reconstruction–Separation From Undersampled ()-Space via Low-Rank Plus Sparse Prior. IEEE Transactions on Medical Imaging, 2014; 33(8), 1689–1701.
  • [25] Otazo R, Emmanuel C and Sodickson DK. Low-rank and sparse matrix decomposition for accelerated DCE-MRI with background and contrast separation. In Proceedings of the ISMRM Workshop on Data Sampling and Image Reconstruction, 2013.
  • [26] Otazo R, Candes E and Sodickson DK. Low-rank plus sparse matrix decomposition for accelerated dynamic MRI with separation of background and dynamic components. Magnetic Resonance in Medicine, 2015; 73(3), 1125–1136.
  • [27] Singh V, Tewfik AH and Ress DB. Under-sampled functional MRI using low-rank plus sparse matrix decomposition. IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2015; pp. 897-901.
  • [28] Burger M, Dirks H and Schönlieb C-B. A variational model for joint motion estimation and image reconstruction, SIAM Journal on Imaging Sciences, 2018; pp. 94–128.
  • [29] Aviles-Rivero AI, Williams G, Graves M and Schönlieb CB. CS+M: A Simultaneous Reconstruction and Motion Estimation Approach for Improving Undersampled MRI Reconstruction. In Proceedings of the 26th Annual Meeting ISMRM, 2018.
  • [30]

    Tomasi C and Kanade T. Shape and motion from image streams under orthography: a factorization method. International Journal of Computer Vision (IJCV), 1992; 9(2), pp. 137-154.

  • [31] Gilland DR, Mair BA, Bowsher JE and Jaszczak RJ. Simultaneous reconstruction and motion estimation for gated cardiac ECT. IEEE Transactions on Nuclear Science, 2002; 49(5), pp. 2344-2349.
  • [32] Huynh N, Lucka F., Zhang E, Betcke M, Arridge S, Beard P and Cox B. Sub-sampled Fabry-Perot photoacoustic scanner for fast 3D imaging, Proc.SPIE, 2017; vol.10064.
  • [33] Burger M, Dirks H , Frerking L, Hauptmann A, Helin T and Siltanen S. A variational reconstruction method for undersampled dynamic X-ray tomography based on physical motion models, Inverse Problems, 2017; pp. 1–24.
  • [34] Wood M, and Henkelman RM. MR image artifacts from periodic motion. Medical physics, 1985; 12.2: 143-151.
  • [35] Van de Walle R, Lemahieu I and Achten E. Magnetic resonance imaging and the reduction of motion artifacts: review of the principles. Technology and Health Care, 1997; 5(6), 419-435.
  • [36] Torr PH and Zisserman A. Feature based methods for structure and motion estimation. The International Conference on Computer Vision - Workshop on vision algorithms (ICCVW), 1999; pp. 278–294.
  • [37] Irani M and Anandan P. About direct methods. The International Conference on Computer Vision - Workshop on vision algorithms (ICCVW), 1999; pp. 267–277.
  • [38]

    Horn BK and Schunck BG. Determining optical flow. Artificial intelligence, 1981; 17(1-3), 185-203.

  • [39] Aubert G, Deriche R and Kornprobst P. Computing optical flow via variational techniques. SIAM Journal on Applied Mathematics, 1999; 60(1), 156-182.
  • [40]

    Zach C, Pock T and Bischof H. A duality based approach for realtime TV-L1 optical flow. In Joint Pattern Recognition Symposium Springer, 2007; pp. 214–223.

  • [41] Perez JS, Meinhardt-Llopis E and Facciolo G. TV-L1 optical flow estimation. Image Processing On Line (IPOL), 2013; pp. 137–150.
  • [42] Chambolle A and Pock T. A first-order primal-dual algorithm for convex problems with applications to imaging. Journal of Mathematical Imaging and Vision, 2011; vol. 40 pp. 120–145.
  • [43] Schloegl M, Holler M, Schwarzl A, Bredies K and Stollberger R. Infimal convolution of total generalized variation functionals for dynamic MRI. Magnetic resonance in medicine, 2017; 78(1), 142-155.
  • [44] Wissmann L, Santelli C, Segars WP and Kozerke S. MRXCAT: Realistic numerical phantoms for cardiovascular magnetic resonance. Journal of Cardiovascular Magnetic Resonance, 2014; 16(1), 63.
  • [45] Wang Z, Bovik AC, Sheikh HR and Simoncelli EP. Image quality assessment: from error visibility to structural similarity. IEEE Transactions on Image Processing (TIP), 2004; 13(4), pp. 600-612.
  • [46] Grosse R, Johnson MK, Adelson EH and Freeman WT. Ground truth dataset and baseline evaluations for intrinsic image algorithms. IEEE International Conference in Computer Vision (ICCV), 2009; pp. 2335-2342.
  • [47] Rudin LI, Osher S and Fatemi E. Nonlinear total variation based noise removal algorithms. Journal of Physics D, 1992; 60:259–268.
  • [48] Schlögl M, Holler M, Bredies K and Stollberger R. A variational approach for coil-sensitivity estimation for undersampled phase-sensitive dynamic MRI reconstruction. In Proceedings of the 23th Annual Meeting ISMRM, 2015.

1 Outline

This document extends the clinical and technical discussion presented in the main paper, this, with the aim of offering further details regarding our experiments and results. The remaining document is structured as follows.

  • [noitemsep]

  • Section 2: More qualitative results for a closer inspection of the reconstructed data.

  • Section 3: We provide more details about the experimental settings and the quantitative results.

  • Section 4: We give an explicit definition and motivation of the metrics used for the quantitative analyses.

2 Supplementary Visual Results

To further support the quantitative and qualitative results presented in the main paper, in this section we offer visual comparison of our approach against the L+S [25] method, which had the second best approximation as reported in Table 1 in the main paper.

Figure 8 shows the results of both methods along with the gold-standard for visual quality assessment. For this comparison, we used the cardiac datasets 3 and 4, cine and perfusion correspondingly, and selected some clinically meaningful frames. Whilst the L+S technique frequently produces a good reconstruction, a closer inspection of 8 shows that our approach performs better in terms of preserving both the anatomical shape and details. Moreover, we can observe less blurry effects from our approach (pointed out with yellow arrows in the zoom-in views). The results in Figure 11 echoes similar findings for the application of cardiac perfusion, in which our approach generates a more visually meaningful reconstructions. This is more clearly reflected with the red arrows in the zoom-in views.

Figure 8: Visual comparisons of our and the L+S approaches using Dataset III. We display of reconstructed samples from different cardiac cycles, in which yellow arrows point at regions of interest.
Figure 9: Reconstructed samples of cardiac perfusion (Dataset IV), in which visual comparisons between our and L+S approaches, along with the gold-standard, is shown. Zoom-in views indicate interesting regions for comparison.
Figure 10: Comparison performance on the temporal domain of our approach against three from the body of literature. (a) and (b) display a reconstructed sample for each evaluated scheme along with the gold-standard using Datasets V and VI. (a.1) and (b.1) show the corresponding temporal signal intensity profiles of a region of interest (blue square) in the LV cavity.
Figure 11: Reconstructed sample from the brain perfusion dataset along with few visualisation plots of the physical motion used to improve the image reconstruction.

3 Further Details and Analysis of the Results

For clarification purposes of the experimental results, in this section, we discuss details regarding the multi-coil reconstructions and the elapsed time of our approach.

Figure 12: Elapsed time comparison between L+S and ours approaches.

In this work, we used both single- and multi- coill MRI data (see Section 3.1 in the main paper for detailed description of the datasets). Whilst for the multi-coil data generated using the MRXCAT framework, the coil sensitivities were computed using the Biot-Savart law (refer to [44] for further details), for the non-simulated data, the coil sensitivities were estimated using [48].

Computational time requirements. We ran a set of experiments to demonstrate the computational feasibility of our reconstruction scheme. In particular, we compared the elapsed time between the L+S and our approach on a CPU basis and using all datasets described in Section 3.1 in the main paper. The experiments were designed to measure only the elapsed time required by the optimisation task. We found that our approach required an average of more computational time than that of [26]. See Figure  12 for three visualisations of the elapsed time with different datasets sizes. This is a slight increment in time that comes with the benefits of having a higher image quality even during sever undersampling factors. Moreover, our implementation can be significantly improved on a CUDA-based implementation.

4 Image-Quality Metrics.

For clarification purposes and with the aim of avoiding some ambiguity regarding the form of how the two quality metrics were computed. We describe the explicit definition of the used metrics, which computation is as follows:

Structural Similarity Index (SSIM). Let and be the reconstructions coming from the gold-standard and a given approximation from a reconstruction scheme, and let

be the local mean, variances and covariance of

and respectively. Then the SSIM is computed as:

where and are constants to ensure numerical stability when the denominator is close to zero. In our experiments, we set standard values used for as and .

Inverted Localised Mean Squared Error (sLMSE). It is to be noted that we are not using the most common metrics for comparing pixel-wise two reconstructions, which are Peak signal-to-noise ratio (PSNR) and MSE. Our main motivation for using the sLMSE is that PSNR and MSE are too stringent (based on a global absolute difference, where a small incorrectly reconstructed edge can dominate the error). Given and - a gold-standard and an approximation reconstruction - as input, we compute the LMSE as the MSE summed over patches, , of size and spaced in steps of 10 as:

to be consistent with the SSIM, we compute a nomalised and inverted version so that the maximum output value is 1, meaning that, as s closer to 1 as higher quality reconstruction which reads:

2 Supplementary Visual Results

To further support the quantitative and qualitative results presented in the main paper, in this section we offer visual comparison of our approach against the L+S [25] method, which had the second best approximation as reported in Table 1 in the main paper.

Figure 8 shows the results of both methods along with the gold-standard for visual quality assessment. For this comparison, we used the cardiac datasets 3 and 4, cine and perfusion correspondingly, and selected some clinically meaningful frames. Whilst the L+S technique frequently produces a good reconstruction, a closer inspection of 8 shows that our approach performs better in terms of preserving both the anatomical shape and details. Moreover, we can observe less blurry effects from our approach (pointed out with yellow arrows in the zoom-in views). The results in Figure 11 echoes similar findings for the application of cardiac perfusion, in which our approach generates a more visually meaningful reconstructions. This is more clearly reflected with the red arrows in the zoom-in views.

Figure 8: Visual comparisons of our and the L+S approaches using Dataset III. We display of reconstructed samples from different cardiac cycles, in which yellow arrows point at regions of interest.
Figure 9: Reconstructed samples of cardiac perfusion (Dataset IV), in which visual comparisons between our and L+S approaches, along with the gold-standard, is shown. Zoom-in views indicate interesting regions for comparison.
Figure 10: Comparison performance on the temporal domain of our approach against three from the body of literature. (a) and (b) display a reconstructed sample for each evaluated scheme along with the gold-standard using Datasets V and VI. (a.1) and (b.1) show the corresponding temporal signal intensity profiles of a region of interest (blue square) in the LV cavity.
Figure 11: Reconstructed sample from the brain perfusion dataset along with few visualisation plots of the physical motion used to improve the image reconstruction.

3 Further Details and Analysis of the Results

For clarification purposes of the experimental results, in this section, we discuss details regarding the multi-coil reconstructions and the elapsed time of our approach.

Figure 12: Elapsed time comparison between L+S and ours approaches.

In this work, we used both single- and multi- coill MRI data (see Section 3.1 in the main paper for detailed description of the datasets). Whilst for the multi-coil data generated using the MRXCAT framework, the coil sensitivities were computed using the Biot-Savart law (refer to [44] for further details), for the non-simulated data, the coil sensitivities were estimated using [48].

Computational time requirements. We ran a set of experiments to demonstrate the computational feasibility of our reconstruction scheme. In particular, we compared the elapsed time between the L+S and our approach on a CPU basis and using all datasets described in Section 3.1 in the main paper. The experiments were designed to measure only the elapsed time required by the optimisation task. We found that our approach required an average of more computational time than that of [26]. See Figure  12 for three visualisations of the elapsed time with different datasets sizes. This is a slight increment in time that comes with the benefits of having a higher image quality even during sever undersampling factors. Moreover, our implementation can be significantly improved on a CUDA-based implementation.

4 Image-Quality Metrics.

For clarification purposes and with the aim of avoiding some ambiguity regarding the form of how the two quality metrics were computed. We describe the explicit definition of the used metrics, which computation is as follows:

Structural Similarity Index (SSIM). Let and be the reconstructions coming from the gold-standard and a given approximation from a reconstruction scheme, and let

be the local mean, variances and covariance of

and respectively. Then the SSIM is computed as:

where and are constants to ensure numerical stability when the denominator is close to zero. In our experiments, we set standard values used for as and .

Inverted Localised Mean Squared Error (sLMSE). It is to be noted that we are not using the most common metrics for comparing pixel-wise two reconstructions, which are Peak signal-to-noise ratio (PSNR) and MSE. Our main motivation for using the sLMSE is that PSNR and MSE are too stringent (based on a global absolute difference, where a small incorrectly reconstructed edge can dominate the error). Given and - a gold-standard and an approximation reconstruction - as input, we compute the LMSE as the MSE summed over patches, , of size and spaced in steps of 10 as:

to be consistent with the SSIM, we compute a nomalised and inverted version so that the maximum output value is 1, meaning that, as s closer to 1 as higher quality reconstruction which reads:

3 Further Details and Analysis of the Results

For clarification purposes of the experimental results, in this section, we discuss details regarding the multi-coil reconstructions and the elapsed time of our approach.

Figure 12: Elapsed time comparison between L+S and ours approaches.

In this work, we used both single- and multi- coill MRI data (see Section 3.1 in the main paper for detailed description of the datasets). Whilst for the multi-coil data generated using the MRXCAT framework, the coil sensitivities were computed using the Biot-Savart law (refer to [44] for further details), for the non-simulated data, the coil sensitivities were estimated using [48].

Computational time requirements. We ran a set of experiments to demonstrate the computational feasibility of our reconstruction scheme. In particular, we compared the elapsed time between the L+S and our approach on a CPU basis and using all datasets described in Section 3.1 in the main paper. The experiments were designed to measure only the elapsed time required by the optimisation task. We found that our approach required an average of more computational time than that of [26]. See Figure  12 for three visualisations of the elapsed time with different datasets sizes. This is a slight increment in time that comes with the benefits of having a higher image quality even during sever undersampling factors. Moreover, our implementation can be significantly improved on a CUDA-based implementation.

4 Image-Quality Metrics.

For clarification purposes and with the aim of avoiding some ambiguity regarding the form of how the two quality metrics were computed. We describe the explicit definition of the used metrics, which computation is as follows:

Structural Similarity Index (SSIM). Let and be the reconstructions coming from the gold-standard and a given approximation from a reconstruction scheme, and let

be the local mean, variances and covariance of

and respectively. Then the SSIM is computed as:

where and are constants to ensure numerical stability when the denominator is close to zero. In our experiments, we set standard values used for as and .

Inverted Localised Mean Squared Error (sLMSE). It is to be noted that we are not using the most common metrics for comparing pixel-wise two reconstructions, which are Peak signal-to-noise ratio (PSNR) and MSE. Our main motivation for using the sLMSE is that PSNR and MSE are too stringent (based on a global absolute difference, where a small incorrectly reconstructed edge can dominate the error). Given and - a gold-standard and an approximation reconstruction - as input, we compute the LMSE as the MSE summed over patches, , of size and spaced in steps of 10 as:

to be consistent with the SSIM, we compute a nomalised and inverted version so that the maximum output value is 1, meaning that, as s closer to 1 as higher quality reconstruction which reads:

4 Image-Quality Metrics.

For clarification purposes and with the aim of avoiding some ambiguity regarding the form of how the two quality metrics were computed. We describe the explicit definition of the used metrics, which computation is as follows:

Structural Similarity Index (SSIM). Let and be the reconstructions coming from the gold-standard and a given approximation from a reconstruction scheme, and let

be the local mean, variances and covariance of

and respectively. Then the SSIM is computed as:

where and are constants to ensure numerical stability when the denominator is close to zero. In our experiments, we set standard values used for as and .

Inverted Localised Mean Squared Error (sLMSE). It is to be noted that we are not using the most common metrics for comparing pixel-wise two reconstructions, which are Peak signal-to-noise ratio (PSNR) and MSE. Our main motivation for using the sLMSE is that PSNR and MSE are too stringent (based on a global absolute difference, where a small incorrectly reconstructed edge can dominate the error). Given and - a gold-standard and an approximation reconstruction - as input, we compute the LMSE as the MSE summed over patches, , of size and spaced in steps of 10 as:

to be consistent with the SSIM, we compute a nomalised and inverted version so that the maximum output value is 1, meaning that, as s closer to 1 as higher quality reconstruction which reads: