up
PDF (letter size)
PDF (legal size)

A comparative study and simulation of HYPR based algorithms

Nasser M. Abbasi

California State University, Fullerton. Summer 2008   Compiled on November 11, 2018 at 10:55pm

Highly constrained backprojection reconstruction

1 Introduction

In the context of reconstruction of MRI images from K-space data sets, the following are two desirable properties which are difficult to achieve simultaneously: High spatial resolution and High temporal resolution.

The first requires longer acquisition time while the second requires less time be consumed acquiring each image. Hence the inherent conflict in achieving both simultaneously. One possible remedy is to undersample image acquisition (along a radial or other trajectories such as Cartesian) which results in speed up of data acquisition, hence improving the temporal resolution. Next, an appropriate image reconstruction method is applied to the acquired data which attempts to compensate for some of the effects of the image undersampling.

Due to undersampling, streaking artifacts will be present in the final image. These streaking artifacts become more visible the larger the undersampling. Radial undersampling is a preferred method of acquisition compared to using other trajectories such as Cartesian: ”the aliasing artifacts from radial under-sampling usually appear as streaks, which are visually less distracting than the wrap-around artifacts obtained with Cartesian under-sampling.”[6]

Mathematically, the problem of image reconstruction from undersampled K-space data is an inverse problem: ”Mathematically, the reconstruction problem from sparse K-space samples is an ill-posed inverse problem with infinitely many solutions[6]

Highly constrained backprojection reconstruction (HYPR) was introduced recently for the reconstruction of radially undersampled MRI images. HYPR is able to reconstruct these images with less visible artifacts while maintaining good SNR[4].

The method starts with the construction of one composite image made up from filtered backprojections of a large number of projections (each of these projections is the inverse Fourier transform of a corresponding radial lines from the K-space data set. This follows from the central slice theorem).

Since the composite image is made up of individual images collected over longer time period, it posses good spatial characteristics. In addition, its SNR is larger since a larger amount of images data is contained in it. However, the composite image temporal characteristics are poor since it combines images that were generated from varying time instances into a single image.

A weight image is then constructed from the ratio of a small number of unfiltered backprojections obtained from the original projections and from the composite image. Since the weight image is constructed from images that span a smaller time window than the case is with the composite image, the weight image posses good temporal characteristics. However, its spatial characteristic is poor due to the undersampling effect in the original data.

HYPR now generates a new image by multiplying the weight image with the composite image. This results in a HYPR image which combines the best characteristics found in the composite image and in the weight image, resulting in an image with good SNR, good temporal and spatial characteristics and with limited artifacts. The above process is repeated for the next HYPR image reconstruction until all the K-space data set is processed. The relationship between the weight image and the composite image and the resulting HYPR image is summarized in the following diagram.

pict

2 Discussion

The following HYPR based algorithms were analyzed in terms of their mathematical formulation. In addition, their properties were studied by simulation[7] under different conditions. The algorithms are: Original HYPR (O-HYPR)[4], Wright-Huang variation of HYPR (W-HYPR)[2], Iterative HYPR (I-HYPR)[1] using original HYPR as its kernel, Iterative HYPR using Wright-Huang HYPR (IW-HYPR)[7] as its kernel, and HYPR-LR[5].

For each algorithm, its mathematical formulation is given, its attributes and the situations in which the algorithm is known to work best and where the algorithm can have difficulty in terms of the quality of reconstruction are both outlined.

2.1 Original HYPR

  1. Mathematical formulation:  HYPR image for time frame \(k\) is \(J_{k}\) and is given by \(J_{k}=\frac{1}{N_{p}}C\{\displaystyle \sum \limits _{i=1}^{N_{p}}} \frac{P_{i}}{P_{c_{i}}}\) where \(N_{p}\) is the number of projections per time frame, \(C\) is the composite image, \(P_{i}\) is the unfiltered backprojection image from a projection obtained from the original data and \(P_{c_{i}}\) is the unfiltered backprojection image from a projection obtained from the composite image \(C\). The subscript \(i\) indicates that the angle \(\theta _{i}\) was used for the backprojection operation.
  2. Attributes: Generates HYPR images with good SNR and good temporal and spatial characteristics.
  3. Situations where algorithm is best applied to: Images with high sparsity and limited object movements.
  4. Situations where algorithm can face difficultly: When objects with very different temporal behavior are very close to each other: ”it has been shown that HYPR reconstruction results in inaccurate representation of signal changes in more spatially and temporally complex data[1] In addition, there is also the issue of crosstalk, which results when two objects with varying temporal characteristics are close to each (case verified by simulation): ”as the sparsity and spatiotemporal correlation deteriorate, there can be crosstalk of signals from different portions of the imaging volume[5]

2.2 Wright-Huang HYPR

  1. Mathematical formulation:  HYPR image for time frame \(k\) is \(J_{k}=C\ \left ({\displaystyle \sum \limits _{i=1}^{N_{p}}} P_{i}\left /{\displaystyle \sum \limits _{i=1}^{N_{p}}} P_{c_{i}}\right . \right ) \).
  2. Attributes: In addition to the attributes of the original HYPR, this algorithm is able to maintain good SNR during reconstruction when noise is present in a way that is more acceptable than with the original HYPR. This was confirmed using simulation when this algorithm was used as the iterative step of iterative HYPR where it was observed that the algorithm has suppressed noise amplification better than the original HYPR was able to do with the same data.
  3. Situations where algorithm can face difficultly: Similar situations as with the original HYPR discussed above.

2.3 Iterative HYPR using original HYPR

  1. Mathematical formulation:  This method is an iterative HYPR. In the initial step the original HYPR algorithm is run. The set of generated HYPR images are then used in the next step as the composite images for the purpose of generating the following set of HYPR images. Hence in each step beyond the initial step, the composite image is replaced as many times as there are time frames. This composite image replacement causes better temporal characteristic in the generated HYPR images as the composite image being used each time would span smaller time window than the case would be with the original HYPR algorithm.
  2. Attributes: Improves over original HYPR by improving the temporal characteristics of the reconstructed images. The error between HYPR image and user images (as measured by relative mean squared error, RMSE) is reduced with more iterations. When no noise was present, the error reduction continued with more iterations but after about 5-10 iterations, little improvement was observed to justify the iterative process continuing for much longer.
  3. Situations where algorithm is best applied to: Images with complex or non-uniform temporal characteristics and nonsparse objects.
  4. Situations where algorithm can face difficultly: When noise is present in original MRI data, I-HYPR will result in amplification of noise as additional iterative steps are taken and hence reducing the accuracy of reconstruction.

2.4 Iterative HYPR using Wright-Huang HYPR

This new algorithm first conceived and implemented during this study. Simulation of this new algorithm confirmed that this algorithm reduces noise amplification by a much larger amount than I-HYPR could during the iterative process.

  1. Mathematical formulation:  Similar to I-HYPR but uses Wright-Huang HYPR as its kernel.
  2. Attributes: Similar to I-HYPR except it is able to improve RMSE as described above more than I-HYPR.
  3. Situations where algorithm is best applied to: When more improvement of the temporal characteristics of the HYPR images are desired but noise in data is an issue.
  4. Situations where algorithm can face difficultly: As with I-HYPR, the main difficulty is in how to determine the number of iterative steps needed to terminate the algorithm. This becomes more important when noise is present.

2.5 HYPR-LR

  1. Mathematical formulation:  HYPR image for time frame \(k\) is \(J_{k}=C\ \left [ \left ( F\otimes{\displaystyle \sum \limits _{i=1}^{N_{p}}} \tilde{P}_{i}\right ) \left / \left ( F\otimes{\displaystyle \sum \limits _{i=1}^{N_{p}}} \tilde{P}_{c_{i}}\right ) \right . \right ] \) Where \(\tilde{P}_{i}\) is the filtered backprojection from the original data and \(\tilde{P}_{c_{i}}\) is the filtered backprojection from projections obtained from the composite image. \(F\) is a low pass filter which is convolved with the backprojection images.
  2. Attributes: ”Unlike HYPR, the new HYPR LR method can be applied to images acquired with arbitrary k-Space trajectories. and reconstruction time is significantly shorter than for iterative methods and the original HYPR algorithm.”[5]. The above quote from the cited paper where it refers to shorter reconstruction was not verified by simulation due to time limitation. However, HYPR-LR was tested in simulation and was found to perform well. [7]
  3. Situations where algorithm is best applied to: ”It is suitable for a broad range of medical imaging applications involving serial changes in image sequence[5]
  4. Situations where algorithm can face difficultly: The problem of crosstalk which was mentioned in the original HYPR remains present here. However, the use of sliding window for the construction of the composite image can be used to reduce this. This was not verified by simulation due to time limitation.

3 Simulation result

HYPR simulation software allows one to execute many scenarios and test cases. Here we show the result of two studies that used a set of images (the phantom clip) supplied to us by GE Healthcare where the images exhibit large degree of temporal and spatial dynamics. The HYPR algorithms were run using this clip as input both under the presence of noise and without noise. Noise was Gaussian with zero mean and standard deviation was set at \(5\%\) of the maximum projection from all the original projections. For all test cases, 8 time frames and 8 projections per time frame was used. For the iterative algorithms, 10 iterations were used. The results below are the average RMSE value, which represent the average error in the reconstruction of HYPR images. The smaller this value, the more accurate the algorithm is considered1 .







test O-HYPR W-HYPR HYPR-LR
I-HYPR
\(1^{st}\)\(10^{th}\)
IW-HYPR
\(1^{st}\)\(10^{th}\)






No noise \(6.83\) \(6.77\) \(6.7\)
\(6.83\)\(5.41\)
\(6.77\)\(5.55\)






With noise \(10.76\) \(9.55\) \(13.7\)
\(10.76\)\(14.22\)
\(9.55\)\(11.436\)






3.1 Noise analysis

Handling of noise by each algorithm was analyzed as follows. A copy of the noise signal being added to each projection was used on its own as the input to each HYPR algorithm. In other words, each noise signal was treated as a projection on its own. When each test starts, two separate computations are started: one which process the original projection with the noise signal added to it (in quadrature), and another which process the noise signal only. At the end of the above two separate computations, 2 sets of HYPR images would result. The mean \(\mu \) and the standard deviation \(\sigma \) of each HYPR image generated from the noise computation was then computed and the \(rmse=\sqrt{\mu ^{2}+\sigma ^{2}}\)for each HYPR image was found. The following table contains the average of the above rmse values over all the HYPR images that was generated.




O-HYPR W-HYPR HYPR-LR



\(0.1004\) \(0.0004\) \(0.0938\)



3.2 Analysis of simulation results

The first table above shows that without noise, O-HYPR, W-HYPR and HYPR-LR performed equally well. When iterative HYPR was run, we observe that I-HYPR and IW-HYPR performed equally well.

When noise was added, W-HYPR was more accurate than O-HYPR. HYPR-LR did not perform as well. But we must note that HYPR-LR can be used with different low pass filters and the size of each filter can be altered as well. Hence it is possible that there exist different low pass filter which can do better than the one used in this particular test. We notice also that when iterative HYPR was run with noise present, the error became larger with more iterations. This is because noise was being amplified in the process. Notice however that IW-HYPR had less noise amplification than I-HYPR. The second table above shows how each algorithm responded to the noise signal only. We observe that W-HYPR had a much smaller rmse. This correlates well with the findings of using IW-HYPR vs. I-HYPR given in the first table above.

4 Conclusion

Five HYPR based MRI image reconstruction algorithms were analyzed and simulated. Each algorithm has different attributes that need to be examined based on the type of data and the type of acquisition before selecting which algorithm to use. Therefore, the choice of which algorithm to select needs to be examined on a case by case basis. However, there are general guidelines that we can propose in selecting an algorithm.  When noise is present, maintaining a good SNR is a requirement that leads one to select the W-HYPR. When the images are less sparse and the temporal characteristics are more dynamic, one can choose the HYPR-LR algorithm. When the images are more sparse and motion of objects is less prevalent and noise is limited, then one can select the O-HYPR.

Finally, when improvement of the temporal characteristics of the generated HYPR images are needed, I-HYPR can be used. If noise is present, IW-HYPR is the preferred method since it can suppress noise amplification more than I-HYPR.

References

[1]    Iterative projection reconstruction of time-resolved images using highly-constrained back-projection (HYPR) by Rafael L. O’Halloran, Zhifei Wen, James H. Holmes, Sean B. Fain

[2]    Time-Resolved MR Angiography With Limited Projections Yuexi Huang and Graham A. Wright

[3]    Principles of computerized Tomographic imaging by Kak and Staney

[4]    Highly Constrained Backprojection for Time-Resolved MRI by C. A. Mistretta, Wieben,z J, Velikina,W. Block,J. Perry,Y. Wu. K. ohnson and Y. Wu.

[5]    Improved Waveform Fidelity Using Local HYPR Reconstruction (HYPR LR). Kevin M. Johnson,Julia Velikina, Yijing Wu, Steve Kecskemeti, Oliver Wieben Charles A. Mistretta

[6]    Projection Reconstruction MR Imaging Using FOCUSS. Jong Chul Ye, Sungho Tak, Yeji Han,and Hyun Wook Park

[7]    HYPR reports by Nasser M. Abbasi HTML