 Research
 Open Access
Optimal compressed sensing reconstructions of fMRI using 2D deterministic and stochastic sampling geometries
 Oliver Jeromin^{1, 2},
 Marios S Pattichis^{1}Email author and
 Vince D Calhoun^{1, 3}
https://doi.org/10.1186/1475925X1125
© Jeromin et al; licensee BioMed Central Ltd. 2012
 Received: 6 June 2011
 Accepted: 2 May 2012
 Published: 20 May 2012
Abstract
Background
Compressive sensing can provide a promising framework for accelerating fMRI image acquisition by allowing reconstructions from a limited number of frequencydomain samples. Unfortunately, the majority of compressive sensing studies are based on stochastic sampling geometries that cannot guarantee fast acquisitions that are needed for fMRI. The purpose of this study is to provide a comprehensive optimization framework that can be used to determine the optimal 2D stochastic or deterministic sampling geometry, as well as to provide optimal reconstruction parameter values for guaranteeing image quality in the reconstructed images.
Methods
We investigate the use of frequencyspace (kspace) sampling based on: (i) 2D deterministic geometries of dyadic phase encoding (DPE) and spiral low pass (SLP) geometries, and (ii) 2D stochastic geometries based on random phase encoding (RPE) and random samples on a PDF (RSP). Overall, we consider over 36 frequencysampling geometries at different sampling rates. For each geometry, we compute optimal reconstructions of single BOLD fMRI ON & OFF images, as well as BOLD fMRI activity maps based on the difference between the ON and OFF images. We also provide an optimization framework for determining the optimal parameters and sampling geometry prior to scanning.
Results
For each geometry, we show that reconstruction parameter optimization converged after just a few iterations. Parameter optimization led to significant image quality improvements. For activity detection, retaining only 20.3% of the samples using SLP gave a mean PSNR value of 57.58 dB. We also validated this result with the use of the Structural Similarity Index Matrix (SSIM) image quality metric. SSIM gave an excellent mean value of 0.9747 (max = 1). This indicates that excellent reconstruction results can be achieved. Median parameter values also gave excellent reconstruction results for the ON/OFF images using the SLP sampling geometry (mean SSIM > =0.93). Here, median parameter values were obtained using meanSSIM optimization. This approach was also validated using leaveoneout.
Conclusions
We have found that compressive sensing parameter optimization can dramatically improve fMRI image reconstruction quality. Furthermore, 2D MRI scanning based on the SLP geometries consistently gave the best image reconstruction results. The implication of this result is that less complex sampling geometries will suffice over random sampling. We have also found that we can obtain stable parameter regions that can be used to achieve specific levels of image reconstruction quality when combined with specific kspace sampling geometries. Furthermore, median parameter values can be used to obtain excellent reconstruction results.
Keywords
 Compressive sensing
 MRI
 fMRI
 Numerical optimization
Background
Using compressive sensing (CS), we can recover certain noisefree signals and images exactly from limited numbers of kspace samples. In theory, to reconstruct images from a limited number of samples, we require that the signal exhibits sparsity and incoherence (lack of correlation) between the sensing basis and the transform basis used for reconstruction. For onedimensional signals with N samples that are composed of T "spikes", perfect reconstuction can be obtained from O(T ⋅ log(N)) frequencydomain samples [1].
In what follows, we demonstrate a reconstruction paradigm that explores several kspace sampling geometries. We consider both stochastic and deterministic Fourier plane sampling geometries. A limitation of our study is that it is based on retrospectively downsampled data from an original fullysampled EPI sequence. Thus, for practical application of our findings for invivo experiments, the EPI sampling pattern will have to be modified to reflect the proposed downsampling patterns. Here, we note that the use of stochastic geometries that result in nonconsecutive EPI phase encoding is prone to artifacts. This should not lead to problems with the SLP geometry that is found to perform best in our experiments. On the other hand, this observation can lead to reconstruction artifacts in stochastic geometries if they are implemented by moving around kspace randomly. We also note that actual MRI imaging time depends on how many phase encoding lines are acquired. Random geometries are included to comply with the conventional (sparsity and incoherence) compressive sensing requirements.
To measure the quality of the reconstructions, we use an optimization process that solves for the compressive sensing objective function parameters that maximize the peak signal to noise ratio (PSNR) or the structural similarity index matrix (SSIM) obtained by each geometry. Here, we note that the structural similarity index matrix (SSIM) [3] provides for a recent (and arguably far more reliable) image reconstruction quality metric.
kspace sampling geometry class examples
% Retained samples  62.5%  48.5%  40.6%  32.8%  28.1%  26.6%  25.0%  21.9%  20.3% 

Geometry class  
Dyadic Phase Encoded (DPE) 









Random Phase Encoded (RPE) 









Random Samples on a PDF (RSP) 









Spiral Low Pass (SLP) 









Alternatively, an optimal regularization parameter can be determined by the Lcurve method [5]. Here, the “L” refers to the shape of the curve of plotting the solution norm versus the corresponding residual norm. The optimal regularization parameter is then determined by the point that gives the corner of the “L”shaped plot. At the corner point of “L”, we have the point that minimizes the solution norm while also minimizing the residual norm. In this way, the Lcurve approach avoids producing a regularization parameter that overly depends on the noise in the image [5].
We note that our parameter optimization approach is direct over the mean structural similarity index (SSIM) [3] and the PSNR. Here, we note that the PSNR is related to the residual error and it appears that the Lcurve method could be applied to find robust estimates of the regularization parameters. However, the emphasis of this paper is on finding the optimal parameters with respect to the meanSSIM [3]. We refer to [6] for a description of the problems associated with the use of the meansquared error for digital images. On the other hand, the meanSSIM provides for a robust metric that better correlates with human visual perception of the images that are being examined.
Our approach finds the optimal regularization parameters that minimize the meanSSIM over a training set of images. Here, we note that each image will produce a different set of optimal parameters depending on the noise variance and how the sampling geometry captures the structure of the activity region. In the section titled “A Fourier Model for Brain Activity,” we provide a frequencydomain model for brain activity. In this section, we carefully show how brain activity favors the 2D spiral low pass geometry (SLP) in the sense that it captures most of the image energy within the effective bandwidth, while it also allows for TVnorm optimization to reconstruct the activity region.
Estimates of the optimal regularization parameters are determined by taking median values of the optimal values estimated over the training set. The effectiveness of this approach is then evaluated over an independent testing set. As we shall describe in the Results section, the approach works very well giving meanSSIM values above 0.93 for the SLP geometry. Furthermore, note that the median is considered to be a robust, nonparametric statistic that avoids outliers [7]. The median provides for robust estimates of the regularization parameters that maximize the meanSSIM of the reconstructed image. Here, we note the conflicting goals. Clearly, the best reconstructions will require the largest number of kspace samples. Thus, we are interested in determining the minimum acceptable quality that also yields acceptable reconstruction with the minimum number of required frequencydomain samples. To this end, we determine effective PSNR ranges and associate them with different image reconstruction qualities. For the average SSIM, we consider fixed values (SSIM > 0.75 for good results, SSIM > 0.90 for excellent results).
We limit our study to optimizing for the best MRI image reconstructions, and for the best quality activity maps based on the difference image between the “ON” and “OFF” MRI images. For detecting activity regions, we focus on the parameters that lead to highest possible quality in the difference images (activation versus noactivation images). Here, we will avoid issues associated with postprocessing the difference images to better delineate between active and inactive regions of the brain. Clearly though, when the difference image is of the highest possible quality, artifact removal will simply improve upon our results. Furthermore, for very large average SSIM values (SSIM > 0.93), reconstruction quality should be considered to be excellent. For the use of postprocessing methods, we refer to the body of work described in [8–15], which utilize a wide variety of signal detection, correlation analysis, and statistical tests to infer brain activity in fMRI difference images.
Related Work
CS methods have been developed that resulted in improved spatial resolution from reduced kspace sampling. The work presented here falls into this category of kspace undersampling reconstruction methods, and is closely related to that of Lustig, Donoho, and Pauly [16]. Lustig et al. also provides a general framework for kspace sampling and reconstruction using CS theory. While the sampling geometries presented here are designed to satisfy the CS criteria, the random and pseudorandom sampling of kspace could result in slower acquisition due to scanner programming constraints.
Most of the subsequent literature on CS applications to MRI imaging explores the merging of CS theory and other fast MRI acquisition techniques, unique kspace sampling methods, and novel reconstruction algorithms. Gamper, Boesiger, and Kozerke meet the sparsity condition of CS theory by applying the Fourier transformation along the temporal dimension, assuming that some regions of the filedofview change at a high temporal rate while other parts remain stationary or change slowly. Their methods show the effectiveness of CS reconstruction for accelerated dynamic (continuous sampling) MRI reconstruction by comparing them to kt BLAST reconstructions over the same data sets [17]. Their sampling scheme can be described as randomly skipping phaseencoding lines in each dynamic frame. Jung et al. developed kt FOCUSS to provide a general framework beyond kt SENSE and kt BLAST for modelbased dynamic MRI by applying CS theory to randomly sampled reconstruction of the prediction and residual encoding that are significantly sparse [18].
Recent studies have focused on extending the work in [16] to nonFourier bases. Haldar, Hernando, and Liang use tailored spatially selective RF pulses to better satisfy the incoherence requirement of CS theory [19]. Liu, Zou, and Ying apply CS theory to parallel reconstruction using sensitivity encoding. Their extension of SENSE to CS is based upon a reconstruction method using Bregman iteration [20]. Trzakso, Manduca, and Borisch present a CS method that minimizes the L0 seminorm using a redescending Mestimator functions instead of L1 norms typically found in CS literature [21]. The extension of the sparsity measure to multiscale form allows for rapid reconstructions compared to the nontrivial solutions described in [17–21].
Methods that exploit the spatial and/or temporal redundancy of kspace data provide a way of describing another class of kspace undersampling reconstructions. A recent study by Lindquist et al. describes methods for obtaining a 3D echovolumnar fMRI imaging sequence by sampling the small central portion of kspace at a high temporal rate [4]. The sampling trajectory is sampled successively across the fourth dimension in a 4D acquisition instead of successively over each temporal slice, and is constrained to a spiral pattern. Other MRI processing techniques that utilize information across slices collected at successive time intervals can be found in [22–26]. Alternatively, the authors in [27] provide a simple iterative algorithm, termed deconvolutioninterpolation gridding (DING) for reconstructing MRI image from arbitrarilysampled kspace. Compared to the Compressive Sensing methods discussed here, it is important to note that DING does not require regularization, and can also work with kspace trajectories that are not known a priori. Thus, DING is also suitable for cases where patient motion occurs during the scan.
Contributions
The primary contributions of this paper can be summarized as follows:
• Comprehensive comparisons of 36 stochastic and deterministic frequency sampling geometries:
While the theory of compressing sensing suggests the use of stochastic geometries, we show that superior results can be achieved with deterministic frequency sampling geometries. Furthermore, we investigate the use of different sampling rates in order to determine the minimum sampling density that can still provide acceptable reconstruction results. As a consequence of this research, we provide the optimal geometry that can lead to the fastest reconstruction time among all candidate sampling geometries. This represents a contribution over existing methods that are based on a limited number of mostly stochastic sampling geometries.
• Optimal reconstruction parameter ranges for guaranteeing a minimum quality in the reconstruction images:
Current compressive sensing methods do not provide sufficient details on how to set reconstruction parameters for different applications. As a result, in current methods, there is no guidance as to which parameter values will provide reconstruction of sufficient quality. Here, we first provide a direct search optimization framework that is used to provide parameter ranges that can guarantee reconstruction of acceptable image quality. In particular, we provide parameter ranges for reconstruction BOLD fMRI ON & OFF images, as well as difference images that can be used for activity detection. Here, as evidence for the need of the proposed approach, we show that a completely different parameter range is needed for activity detection as opposed to parameters for reconstructing individual ON and OFF images. Furthermore, we found that the fast, spiral low pass sampling geometry can reach PSNR levels above 40 dB after a few iterations. In addition, the spiral low pass sampling geometry gives superior results with all reconstructions yielding mean SSIM values above 0.93.
Methods
Optimal TV Norm and Wavelet Transform Penalty
In what follows, we provide a summary of the proposed optimization method. As we describe below, the basic approach involves a nested optimization approach.
where A denotes the accuracy of the reconstructed image. Here, we will consider two measures of accuracy: (i) the average SSIM index, and (ii) PSNR. To solve (1), we will need to solve the optimization problem of reconstructing f from a limited number of kspace samples.
where ε denotes a small tolerance value, and ${\Vert .\Vert}_{{l}_{2}}$ denotes the l _{ 2 } norm. In (2), F _{ k } takes the 2D discrete Fourier transform of the input image and forces the kspace samples in y to match the ones available in the kth sampling geometry.
where x denotes the column vector of the wavelet coefficients and Ψ is the wavelet transform operator that contains the wavelet basis functions along its columns.
In (4), TV(Ψx) denotes the total variation of the reconstructed signal. In (4)(5), we have a constrained optimization problem in three parameters. This is converted to an unconstrained optimization problem through the use of a Lagrange multiplier as indicated in [16, 28]. In what follows, we will consider different values for α, β and use the software described in [16] to solve (4)(5) for the optimal wavelet coefficients. We provide details for parameter optimization in a separate section.
Data set and fMRI activity detection
In blood oxygenationlevel dependent (BOLD) fMRI, neural activity is detected by changes in the T2^{*} relaxation time due to changes in blood oxygenation levels in response to local activation. All images were acquired on a 3 T Siemens TIM Trio system with a 12channel radio frequency (RF) coil. The fMRI experiment used a standard Siemens gradientecho EPI sequence modified so that it stored real and imaginary data separately. We used a fieldofview (FOV) = 240 mm, slice thickness = 3.5 mm, slice gap = 1 mm, number of slices = 32, matrix size = 64 × 64, TE = 29 ms, and TR = 2 s. The fMRI experiment used a block design with periods of 30 s OFF and 30 s ON. The subjects who participated in this study tapped the fingers of their right hand during the ON period. There were five and a half cycles, starting with OFF and ending with the OFF period.
The BOLD fMRI data was then preprocessed to account for motion artifacts and spatially normalized into the standard Montreal Neurological Institute space. This spatial normalization was then subsampled to 3 × 3 × 4 mm, resulting in 53 × 63 × 46 voxels. An individual slice was then selected that ensured measureable regions of activity based on the task being performed by the test subjects.
Temporal smoothing also tends to better localize detected activity across all temporal slices within a single run. Instead of utilizing the temporal information, our reconstructions are of two individual ON and OFF BOLD images. The less dense collection of temporal samples provides a worstcase scenario for detecting neural activity. We detected neural activity in fMRI by calculating the difference image of the reconstructed ON and OFF images. Individual slices were reconstructed with optimal parameters computed for different sampling geometries. At this juncture, the neural activity detection problem becomes a segmentation problem. The preprocessed sum of squares imagery was used to generate the kspace data in this study by applying the 2D Fast Fourier Transform to the multicoil data.
As such, we will not consider optimization over a variety of different activity detection algorithms. We do note however that the activity detection algorithms that employ lowpass filtering will tend to favor zerofilling over interpolation by compressive sensing or any other method. The reason for this is that lowpass filtering attenuates highfrequency components that are estimated by the interpolation/reconstruction method. Ultimately, any segmentation for BOLD fMRI images will be governed by the accuracy of the reconstructed images themselves. Thus, for considering activity detection, we present the reconstructed difference images for qualitative comparisons. Thus, a more practical implementation of our approach would be studies involving structural MRI data, in which the benefit of a reduced acquisition time provides relief of the anxiety and discomfort patients may experience within an MRI scanner.
A Fourier model for brain activity
In this section, we introduce an analytical model for modeling the activity (ON) image. As we shall derive in this section, the lowerfrequency components can provide for an accurate model for the FourierEnergy concentration as well as a lowTV norm image reconstruction.
Thus, our objective in reconstructing the MRI image is to effectively reconstruct the stimulated region that differentiates f _{ ON } (.) from the rest. From Figure 2, we provide models for each function using a rectangular approximation for each function. The use of a rectangular approximation allows us to develop a model that is separable in the xy coordinate system. The assumption allows us to provide analytical models and evaluate the properties of the resulting Fourier expansions for (7). Also, note that we can always express images as a sum of separable functions (e.g., via the use of the singular value decomposition or the 2D Fourier Transform).
where we have A≫a, to signify that the oscillatory components are much smaller than the average image intensity. Similarly, B≫b for the stimulated region. We can verify that both of these assumptions are valid by visual inspection of the corresponding regions in Figure 2.
where W_{ smallest } refers to the smallest stimulated region that we wish to recover in the Fourier domain. The concept of width can easily be extended to nonrectangular regions.
In (15), it is important to note that we have a function that has a very small TVnorm. In fact, it is an ideal case since we expect that the TVoptimized reconstruction will be well represented by the righthandside of (15). The observation here is that we can use the lowfrequency components to capture most of the Fouriertransform energy of the ideal model as well as provide for a good initialization for optimization for low TV.
The kSpace downsampling geometries
We consider both deterministic and stochastic sampling geometries. For deterministic geometries, we consider the spiral lowpass and dyadic downsampling along the phase encoding direction. Stochastic methods are based on random sampling.
First, we consider geometries that restrict down sampling in the phaseencoded dimension only. Here, we wish to compare kspace sampling geometries in just the phaseencoded dimension to kspace undersampling in both dimensions around the center.
As noted in the literature review, the central region of kspace is essential in obtaining reconstruction performance that is acceptable. Almost all nonCS reconstruction of kspace undersampling included the center of kspace in the data that was sampled. A dyadic sampling geometry class was developed that considers sparse sampling along the phaseencoded dimension of kspace. All geometries of this class are shown in the first row of Table 1. This geometry includes samples from a collection of contiguous frequency encoded vectors centered over the center of kspace. The width of the central region can be varied to generate a number of unique masks that are members of this mask class. Nine unique geometries are generated based on a central region that sample 1/2, 1/3, 1/4, 1/6, 1/8, 1/10, 1/12, 1/16, and 1/32 of the phase encoded samples. Beyond the bounds of the central region of the geometry, additional frequency encoded vectors are sampled every 2^{nd}, 4^{th}, 8^{th}, etc. phase encoded sample until the entire support of the fully sampled kspace has been included. The dyadic characteristic of the gap size between subsequent highfrequency samples was intentionally designed to sample more densely near the center of kspace, and in hopes that such sampling in kspace might coincide with the wavelet transform our reconstruction method utilized. This class will be referred to as the dyadic phase encoded (DPE) geometry class in the remainder of this paper.
We also consider two additional geometry classes which attempt to increase the incoherence between the sensing and transform bases by introducing an element of randomness in how the samples were selected. The first class, which is referred to as the random phase encoded (RPE) geometry class, samples frequencyencoded vectors along random phase encoded samples. This sampling method was selected as a comparison to the dynamic CS MRI in [20]. In two images, single fMRI study, such sampling would not expect to provide acceptable reconstruction due to the possibility of excluding a portion of the central region of kspace. Each of the random phase encoded geometries contains the same number of samples as one of the dyadic phase encoded geometries. Due to the onetoone geometry correspondence across classes, there are nine total geometries in the RPE class. All of the random phase encoded geometries are shown in the second row of Table 1.
The second geometry class that incorporates randomness into the sampling scheme can be described as random sampling along a 1D probability distribution function across the phase encoded samples at each frequencyencoded sample. This geometry class will be referred to the random sampling PDF (RSP) class.
where U and V are the number of kspace samples in the frequency and phase encoded dimension. This function mathematically defines the probability density function on which the random samples are being selected. Our motivation for this geometry was to attempt a compromise between the inclusion of the center of kspace and also imparting incoherence into the CS problem through a pseudorandom sampling geometry. Representative geometries of the RSP class are depicted in the third row in Table 1.
Lastly, we include a geometry class that restricts kspace sampling in both the phase encoded and frequencyencoded dimensions. Here, we generate the nine geometries by a sampling along Cartesian spiral emanating from the center of kspace. Such sampling geometries are more typical of the classical kspace undersampling reconstruction techniques found in literature. This class is referred to as the spiral low pass (SLP) geometry class (see Table 1).
Optimization of CS penalty parameters
Parameter optimization is performed using two separate approaches. In the first approach, we want to estimate the optimal parameter ranges for α, β that can provide a certain level of image quality. This first approach is based on maximizing the PSNR in the reconstruction images. Our second approach is more direct. It is based on maximizing the average SSIM index. In our second approach, we want to recommend specific parameter values for each sampling geometry. The second approach is validated using leaveoneout.
Both approaches assume that we can find optimal parameters when ground truth is available. When ground truth is not available, our expectation is that our training sets will reflect the new datasets that will be collected. To initialize optimization, we first zero pad the images from the original size of 53 × 63 to be of size 64 × 64. We then use the 2D FFT to generate kspace data for our experiments. A baseline reconstruction is obtained by applying the inverseFourier transform to the kspace undersampling matrix. This initial reconstruction will exhibit the aliasing artifacts that the work presented in the introduction seeks to suppress. We then use the NelderMead search algorithm to estimate optimal reconstruction parameters of (1) [29]. Here, we note that the NelderMead algorithm does not require that we evaluate the derivatives. It belongs to the class of direct search methods [30]. This approach is used to generate a collection of optimal parameter values α, β associated with each image.
Then, in our first approach, we estimate the optimal parameter ranges for α, β that can provide a certain level of image quality. To understand how this works, note that minimum image quality requirement can be expressed as a minimum requirement for the PSNR level. We then determine the reconstructed images that exceed the minimum level. The optimal parameter region can then be estimated from the range of values of α, β that correspond to these images. To avoid outliers, the optimal regions need only be met by 75% of the total number of the considered images. We also consider the complexity of each optimal parameter region based on the inner and outer bounding boxes. The inner bounding box is the largest rectangle that can be totally contained inside the optimal parameter region. On the other hand, the outer bounding box represents the smallest rectangle that contains the optimal parameter region. A detailed example of our first approach is also provided in the results section.
Our second approach is more direct. Our goal is to provide specific recommendations for α, β that can be expected to exceed a minimum level of image quality. For this direct approach, we use the average SSIM index value to determine the required level of image quality. We consider an image with an average SSIM index value above 0.9 to represent excellent reconstruction quality. In this second approach, the image reconstruction quality is estimated using a leaveoneout method. In leaveoneout, we construct training sets consisting of all but one of the ON/OFF image pairs. We then report the reconstruction image quality on the remaining ON/OFF image pair that serves as our testing set. The method is then applied over all image pairs. Here, we use the median parameter values estimated over each training set for reconstruction of the testing set. In this approach, we report our results for each image serving as a member of the testing set. Furthermore, for simplicity, we provide the median parameter values over the entire dataset. Also, for completeness, we also provide optimal parameter values for each image pair. This allows the readers to see if the discrepancy between the median and the optimal parameter values will affect the reconstructed image quality.
Results
We provide optimization results independently for the ON, the OFF images, as well as for the difference image. Here, we consider optimization of the difference image for activity detection.
Optimal parameter regions for achieving different image quality bounds
Geometry class  Threshold  Inner bounding box  

Min α  Max α  Min β  Max β  
DPE 62.5%  45 dB  8.0808e05  2.8687e04  2.5758e04  3.9697e04 
DPE 48.5%  40 dB  1.4949e04  2.5253e04  2.3434e04  3.9697e04 
DPE 32.8%  35 dB  2.1111e04  3.0404e04  3.8990e04  4.9293e04 
DPE 21.9%  30 dB  4.2424e04  6.3030e04  2.8081e04  3.2727e04 
SLP 62.5%  45 dB  2.5253e04  4.2424e04  1.4141e04  1.8788e04 
SLP 40.6%  40 dB  1.4949e04  2.8687e04  7.1717e05  2.1111e04 
SLP 20.3%  35 dB  8.0808e05  1.4949e04  4.8485e05  9.4949e05 
SLP 20.3%  30 dB  1.2121e05  1.4949e04  2.5253e05  9.4949e05 
Geometry class  Threshold  Outer bounding box  
Min α  Max α  Min β  Max β  
DPE 62.5%  45 dB  2.2222e05  3.8990e04  1.8788e04  4.6667e04 
DPE 48.5%  40 dB  4.6465e05  3.2121e04  1.6465e04  5.1313e04 
DPE 32.8%  35 dB  1.1818e04  3.7374e04  1.4949e04  5.9596e04 
DPE 21.9%  30 dB  1.4141e04  3.9697e04  2.5253e04  7.3333e04 
SLP 62.5%  45 dB  4.6465e05  4.5859e04  7.1717e05  3.0404e04 
SLP 40.6%  40 dB  2.8687e04  4.2424e04  2.0202e06  3.0404e04 
SLP 20.3%  35 dB  2.2222e05  2.5253e04  2.1212e05  1.8788e04 
SLP 20.3%  30 dB  5.6566e05  2.5253e04  2.1212e05  2.5758e04 
Our experiments reveal that the DPE and SLP geometry class outperform the RPE and RSP classes. It was expected that the RPE class exhibit poor reconstruction, since the geometries do not contain sufficient samples from the center of kspace. And only a single geometry, the RSP geometry including 62.5% of samples resulted in a suggested operating parameter space for achieving “acceptable” reconstruction.
On the other hand, the DPE class achieved excellent reconstruction from two geometries (62.5% and 48.5%) and the SLP class achieved excellent reconstruction from the three geometries that retained the most number of kspace samples (60.5%, 48.5%, and 40.6%). When lowerquality reconstructions are acceptable, we can use sampling geometries that use fewer kspace samples. Twelve additional DPE geometries were found to satisfy the parameter space threshold and intersection procedure. All eighteen (nine each from the ON and OFF imagery) SLP geometries that were included in this experiment satisfied our requirement for acceptable reconstruction.
For activity detection, we also compared the difference images of the original ON and OFF samples, the CS reconstructed ON and OFF samples, and also the reconstruction values of zerofilling the missing kspace samples. Based on our investigations, we have found that it was possible to obtain excellent reconstruction results with very fast kspace sampling geometries. In what follows, we present our best results, obtained with the SLP geometry using only 20.3% of kspace.
Dyadic phase geometry: Mean SSIM results (OFF & ON) over each image pair and median over all pairs
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image10  Image11  Image12  Median 

62.50  0.9938  0.9935  0.9890  0.9932  0.9937  0.9928  0.9880  0.9933  0.9922  0.9921  0.9949  0.9931  0.99315 
48.50  0.9839  0.9806  0.9721  0.9803  0.9842  0.9781  0.9748  0.9749  0.9773  0.9784  0.9844  0.9817  0.97935 
40.60  0.9661  0.9696  0.9512  0.9657  0.9633  0.9558  0.9575  0.9608  0.9596  0.9613  0.9668  0.9658  0.9623 
32.80  0.9347  0.9428  0.9207  0.9358  0.9452  0.9027  0.9223  0.9371  0.9331  0.9223  0.9518  0.9422  0.93525 
28.10  0.9032  0.9277  0.8863  0.9020  0.9173  0.8918  0.9163  0.9068  0.9024  0.8963  0.9268  0.9206  0.905 
26.60  0.8758  0.9097  0.8612  0.8865  0.8928  0.8537  0.8440  0.8918  0.8922  0.8678  0.9068  0.8987  0.88915 
25.00  0.8798  0.8957  0.8379  0.8891  0.8815  0.8705  0.8799  0.8935  0.8764  0.8666  0.8990  0.8972  0.8807 
21.90  0.8452  0.8624  0.8087  0.8569  0.8394  0.8409  0.8507  0.8559  0.8483  0.8380  0.8649  0.8653  0.8495 
20.30  0.7779  0.7816  0.7852  0.7822  0.7928  0.6842  0.6782  0.8073  0.8019  0.7844  0.8179  0.8110  0.7848 
Spiral Low Pass geometry: Mean SSIM results (OFF & ON) over each image pair and median over all pairs
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image 11  Image 12  Median 

62.50  0.9923  0.9938  0.9916  0.9936  0.9942  0.9929  0.9897  0.9934  0.9931  0.9941  0.9952  0.9945  0.9935 
48.50  0.9888  0.9877  0.9878  0.9891  0.9889  0.9870  0.9801  0.9891  0.9869  0.9884  0.9908  0.9850  0.9881 
40.60  0.9856  0.9829  0.9824  0.9848  0.9860  0.9824  0.9726  0.9860  0.9820  0.9841  0.9867  0.9791  0.9835 
32.80  0.9789  0.9762  0.9727  0.9784  0.9811  0.9741  0.9679  0.9815  0.9761  0.9804  0.9828  0.9698  0.9773 
28.10  0.9697  0.9659  0.9624  0.9703  0.9773  0.9649  0.9574  0.9728  0.9671  0.9731  0.9767  0.9682  0.96895 
26.60  0.9660  0.9616  0.9541  0.9674  0.9748  0.9632  0.9552  0.9681  0.9637  0.9699  0.9718  0.9652  0.9656 
25.00  0.9639  0.9567  0.9528  0.9643  0.9716  0.9600  0.9497  0.9657  0.9585  0.9687  0.9713  0.9640  0.96395 
21.90  0.9512  0.9478  0.9342  0.9543  0.9630  0.9468  0.9427  0.9523  0.9456  0.9556  0.9606  0.9531  0.95175 
20.30  0.9424  0.9356  0.9244  0.9442  0.9586  0.9285  0.9300  0.9385  0.9316  0.9467  0.9545  0.9503  0.94045 
Random sampling on PDF geometry: Mean SSIM results (OFF & ON) over each image pair and median over all pairs
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  0.9917  0.3613  0.9946  0.9915  0.9965  0.4869  0.8859  0.4286  0.8429  0.9873  0.8771  0.9968  0.9366 
48.50  0.8839  0.7967  0.3797  0.9873  0.8862  0.7779  0.9661  0.8454  0.8591  0.4017  0.8990  0.9863  0.8715 
40.60  0.9538  0.5490  0.9520  0.9526  0.9626  0.7717  0.7206  0.9405  0.9593  0.9579  0.9668  0.9618  0.9532 
32.80  0.9302  0.8089  0.4618  0.8772  0.8536  0.3084  0.4843  0.4917  0.9077  0.3930  0.6974  0.9397  0.75315 
28.10  0.7025  0.8010  0.4048  0.4846  0.9073  0.6879  0.7178  0.7995  0.4558  0.6144  0.8163  0.5519  0.6952 
26.60  0.7432  0.4272  0.4323  0.5289  0.7348  0.7207  0.3768  0.8983  0.4436  0.5258  0.8413  0.4048  0.52735 
25.00  0.5855  0.9114  0.6392  0.5203  0.4872  0.5883  0.2242  0.8795  0.3602  0.5417  0.9151  0.5172  0.5636 
21.90  0.4465  0.8086  0.4637  0.5236  0.4568  0.5584  0.5955  0.7593  0.8266  0.5380  0.7414  0.5028  0.5482 
20.30  0.7635  0.4661  0.3008  0.5635  0.5954  0.6655  0.5174  0.5689  0.4259  0.4993  0.5045  0.8082  0.54045 
Random phase encoding geometry: Mean SSIM results (OFF & ON) over each image pair and median over all pairs
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image 11  Image 12  Median 

62.50  0.4977  0.7761  0.8038  0.7710  0.9385  0.3035  0.3178  0.5755  0.8920  0.5611  0.9340  0.7678  0.7694 
48.50  0.7205  0.4283  0.4978  0.4190  0.5157  0.3792  0.4124  0.4222  0.5150  0.7992  0.6382  0.4757  0.48675 
40.60  0.8368  0.4455  0.6163  0.3631  0.5303  0.2675  0.6773  0.6491  0.6260  0.6262  0.5328  0.7319  0.62115 
32.80  0.7470  0.4417  0.5449  0.5325  0.5322  0.2600  0.2553  0.7156  0.4721  0.5817  0.5961  0.3384  0.53235 
28.10  0.3020  0.5741  0.3409  0.3296  0.3299  0.2593  0.4213  0.7007  0.3138  0.6112  0.4863  0.3430  0.34195 
26.60  0.2954  0.4523  0.5428  0.4348  0.4034  0.2812  0.2324  0.6164  0.3333  0.4492  0.5069  0.4524  0.442 
25.00  0.6157  0.4613  0.3285  0.4298  0.3138  0.2766  0.6734  0.4857  0.4460  0.6840  0.3523  0.3663  0.4379 
21.90  0.5199  0.2968  0.3911  0.3072  0.4502  0.2662  0.2366  0.4205  0.3047  0.3096  0.5264  0.4307  0.35035 
20.30  0.3137  0.2704  0.4546  0.3005  0.4328  0.3067  0.2069  0.3503  0.2660  0.4082  0.2628  0.3512  0.3102 
Optimal alpha values for dyadic phase encoding geometry (SSIM optimization)
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  5.187*1E04  4.382*1E04  4.185*1E04  6.067*1E04  2.487*1E04  5.630*1E04  4.282*1E04  4.180*1E04  3.834*1E04  2.264*1E04  2.934*1E04  1.333*1E04  4.18*1E04 
48.50  4.872*1E04  7.729*1E04  3.677*1E04  4.199*1E05  7.021*1E04  3.698*1E04  2.637*1E05  3.794*1E04  1.545*1E04  5.969*1E05  6.555*1E04  3.687*1E04  3.69*1E04 
40.60  4.577*1E04  2.174*1E04  6.738*1E04  1.573*1E04  3.090*1E04  2.712*1E04  2.207*1E04  1.597*1E04  1.580*1E04  2.485*1E04  7.190*1E04  9.229*1E05  2.35*1E04 
32.80  1.446*1E04  2.568*1E04  5.670*1E04  1.686*1E04  1.308*1E04  3.551*1E04  3.744*1E04  1.821*1E04  4.293*1E04  1.641*1E04  4.952*1E04  7.518*1E04  3.06*1E04 
28.10  4.738*1E04  2.575*1E04  1.294*1E04  1.245*1E04  9.180*1E05  4.849*1E04  1.502*1E04  3.984*1E04  4.354*1E04  5.536*1E04  6.348*1E05  4.926*1E04  3.28*1E04 
26.60  5.961*1E04  2.208*1E04  2.764*1E04  8.724*1E05  2.422*1E04  1.694*1E03  3.760*1E04  1.450*1E04  1.747*1E04  2.461*1E04  3.492*1E04  1.405*1E04  2.44*1E04 
25.00  1.576*1E04  2.046*1E04  1.709*1E05  1.746*1E04  3.734*1E04  5.181*1E04  2.947*1E04  5.724*1E04  5.362*1E04  1.410*1E04  3.343*1E04  3.248*1E04  3.10*1E04 
21.90  2.257*1E04  6.291*1E04  1.065*1E04  1.095*1E04  2.912*1E04  5.319*1E04  2.600*1E04  1.323*1E04  4.067*1E04  6.829*1E04  9.795*1E04  3.029*1E04  2.97*1E04 
20.30  1.563*1E05  2.289*1E04  1.435*1E05  1.641*1E04  −2.664*1E05  2.285*1E04  3.615*1E04  1.240*1E04  1.230*1E04  3.426*1E05  1.007*1E05  7.462*1E06  7.86*1E05 
Optimal alpha values for spiral lowpass geometry (SSIM optimization)
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  6.742*1E05  2.440*1E04  3.171*1E04  3.987*1E04  2.606*1E04  2.529*1E04  2.591*1E04  1.550*1E04  1.728*1E04  4.957*1E04  3.440*1E04  1.499*1E04  2.56E04 
48.50  1.805*1E04  1.322*1E04  2.967*1E04  1.799*1E04  2.294*1E04  3.108*1E04  1.297*1E04  2.942*1E04  2.484*1E04  4.276*1E04  4.420*1E04  1.798*1E04  2.39E04 
40.60  2.994*1E04  1.038*1E06  2.844*1E04  1.086*1E04  2.632*1E04  2.026*1E04  3.734*1E04  2.843*1E04  6.503*1E04  7.424*1E05  1.903*1E04  2.900*1E04  2.74E04 
32.80  2.261*1E05  9.766*1E07  2.617*1E04  1.784*1E05  1.651*1E04  7.080*1E05  6.673*1E05  1.561*1E04  1.329*1E05  4.016*1E05  9.338*1E05  5.145*1E05  5.91E05 
28.10  5.980*1E06  0.0E + 00  7.971*1E05  1.282*1E06  2.110*1E05  8.163*1E05  2.121*1E04  9.640*1E06  3.081*1E06  1.708*1E05  5.890*1E06  1.438*1E06  7.81E06 
26.60  3.779*1E06  −2.441*1E07  8.136*1E05  9.766*1E07  2.803*1E05  1.820*1E04  3.906*1E06  7.959*1E06  9.220*1E06  1.842*1E05  6.043*1E06  4.747*1E05  8.59E06 
25.00  3.829*1E06  0.0E + 00  1.626*1E04  3.052*1E08  9.683*1E05  1.134*1E04  6.404*1E05  4.971*1E05  5.096*1E06  1.678*1E05  5.286*1E06  4.147*1E05  2.91E05 
21.90  1.008*1E05  3.052*1E08  2.252*1E05  4.578*1E08  9.543*1E05  7.703*1E05  1.178*1E04  3.884*1E05  3.132*1E05  1.179*1E05  1.053*1E05  2.166*1E05  2.21E05 
20.30  5.893*1E06  −5.341*1E08  1.347*1E05  2.441*1E07  9.809*1E06  3.146*1E05  3.906*1E06  4.454*1E05  6.435*1E05  3.173*1E05  1.164*1E05  4.158*1E05  1.26E05 
Optimal alpha values for random sampling on PDF geometry (SSIM optimization)
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  2.704*1E04  3.669*1E04  3.817*1E04  3.978*1E04  2.543*1E04  4.777*1E04  3.090*1E04  2.522*1E04  4.710*1E04  5.479*1E04  2.637*1E04  2.531*1E04  3.38*1E04 
48.50  1.672*1E03  2.404*1E03  2.460*1E03  1.266*1E03  0E + 00  2.185*1E03  1.924*1E03  1.627*1E03  1.606*1E03  8.298*1E04  1.865*1E03  1.848*1E03  1.76*1E03 
40.60  1.098*1E03  2.381*1E03  1.467*1E03  2.400*1E03  4.761*1E06  −5.186*1E05  1.389*1E03  1.380*1E03  1.284*1E03  2.162*1E03  1.273*1E03  3.481*1E03  1.38*1E03 
32.80  1.201*1E03  0E + 00  1.190*1E03  1.107*1E03  1.484*1E03  1.959*1E03  1.928*1E03  1.174*1E03  1.759*1E03  1.375*1E03  1.003*1E03  2.092*1E03  1.29*1E03 
28.10  3.753*1E03  9.521*1E06  4.616*1E03  9.082*1E04  2.677*1E03  −6.645*1E05  1.673*1E03  1.933*1E03  2.526*1E03  1.460*1E03  2.445*1E03  3.229*1E03  2.19*1E03 
26.60  1.980*1E03  1.522*1E03  1.535*1E03  1.040*1E03  1.393*1E03  5.993*1E04  6.140*1E04  8.419*1E04  1.679*1E03  1.061*1E03  1.171*1E03  1.469*1E03  1.28*1E03 
25.00  8.867*1E04  4.137*1E04  1.997*1E03  1.074*1E05  6.785*1E04  1.424*1E03  2.551*1E03  1.506*1E03  4.349*1E04  7.230*1E04  4.732*1E04  1.346*1E03  8.05*1E04 
21.90  1.003*1E03  1.731*1E03  1.562*1E03  2.544*1E03  1.717*1E03  1.979*1E03  1.400*1E03  2.576*1E04  2.223*1E03  9.923*1E04  9.351*1E04  2.016*1E03  1.64*1E03 
20.30  1.388*1E03  2.289*1E03  2.344*1E05  1.739*1E03  3.741*1E03  2.142*1E03  2.476*1E03  1.484*1E03  −6.022*1E05  −2.979*1E05  1.257*1E03  1.934*1E03  1.61*1E03 
Optimal alpha values for random phase encoding (SSIM optimization)
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  3.432*1E04  5.724*1E04  1.449*1E04  2.051*1E04  −3.137*1E05  7.519*1E04  4.396*1E04  4.749*1E04  −5.880*1E06  2.859*1E04  7.813*1E06  6.028*1E04  3.15*1E04 
48.50  −1.633*1E06  4.126*1E04  −2.972*1E05  −6.086*1E05  −1.590*1E05  −1.863*1E04  −1.957*1E04  4.404*1E04  −5.236*1E05  −2.679*1E04  7.813*1E06  5.633*1E04  −2.28*1E05 
40.60  9.952*1E04  4.536*1E04  2.723*1E04  2.109*1E04  2.979*1E05  2.140*1E04  1.763*1E04  4.544*1E04  5.251*1E04  4.368*1E04  −6.852*1E05  8.507*1E04  3.55*1E04 
32.80  5.493*1E07  −4.938*1E04  −6.073*1E06  −1.334*1E04  −5.382*1E05  5.497*1E04  2.771*1E05  −7.382*1E05  −4.323*1E05  −4.177*1E04  −3.918*1E05  −5.866*1E04  −4.85*1E05 
28.10  −3.632*1E05  2.559*1E04  −3.033*1E05  −1.074*1E05  −3.271*1E05  4.256*1E04  7.668*1E04  1.416*1E04  −7.233*1E06  3.029*1E06  −1.099*1E05  6.348*1E06  −2.10*1E06 
26.60  −3.906*1E06  1.526*1E07  −3.009*1E05  −3.311*1E05  −4.181*1E05  2.896*1E05  2.068*1E03  −4.675*1E05  −5.447*1E06  −3.906*1E06  5.127*1E06  −7.563*1E05  −4.68*1E06 
25.00  −1.844*1E04  −4.150*1E06  −1.160*1E04  −4.858*1E05  8.374*1E05  1.191*1E03  6.500*1E04  1.036*1E04  9.277*1E06  2.930*1E06  −7.629*1E08  −1.639*1E05  1.43*1E06 
21.90  9.766*1E07  7.047*1E04  1.367*1E05  1.121*1E03  3.906*1E06  2.021*1E03  6.571*1E03  9.662*1E04  1.068*1E06  9.479*1E05  1.271*1E03  7.481*1E04  7.26*1E04 
20.30  1.318*1E04  2.374*1E05  2.148*1E05  6.586*1E05  2.138*1E04  5.205*1E06  1.572*1E06  1.978*1E05  6.895*1E06  4.661*1E04  4.562*1E05  4.753*1E04  3.47*1E05 
Optimal beta values for dyadic phase encoding geometry (SSIM optimization)
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  1.084*1E03  2.698*1E04  4.536*1E04  −7.898*1E05  2.957*1E04  3.546*1E04  1.299*1E04  5.352*1E04  3.703*1E04  6.006*1E04  5.883*1E04  2.869*1E04  3.62*1E04 
48.50  5.988*1E04  2.676*1E04  3.268*1E04  3.989*1E04  5.205*1E04  5.235*1E04  3.940*1E04  3.740*1E04  5.166*1E04  1.891*1E04  2.092*1E04  3.794*1E04  3.87*1E04 
40.60  9.831*1E04  2.217*1E04  1.020*1E03  1.729*1E04  5.130*1E04  2.909*1E04  −5.859*1E05  4.841*1E04  5.129*1E04  4.080*1E04  7.499*1E04  5.183*1E04  4.99*1E04 
32.80  4.906*1E04  4.087*1E04  4.010*1E04  7.449*1E04  7.143*1E04  4.673*1E04  1.331*1E03  5.408*1E04  1.058*1E03  7.070*1E04  6.733*1E04  4.987*1E04  6.07*1E04 
28.10  8.401*1E04  2.506*1E04  4.290*1E04  4.456*1E04  6.025*1E04  2.600*1E04  5.224*1E04  5.955*1E04  2.392*1E04  1.120*1E03  6.011*1E04  1.297*1E03  5.59*1E04 
26.60  2.075*1E03  3.907*1E04  3.184*1E04  2.024*1E03  5.989*1E04  1.206*1E03  7.061*1E04  2.813*1E04  4.880*1E04  6.309*1E04  6.032*1E04  7.237*1E04  6.17*1E04 
25.00  5.646*1E04  3.037*1E04  3.977*1E04  6.276*1E04  6.339*1E04  2.351*1E04  1.793*1E04  9.519*1E04  4.954*1E04  4.382*1E04  9.376*1E04  7.264*1E04  5.30*1E04 
21.90  1.004*1E03  1.034*1E03  1.033*1E03  2.980*1E04  9.097*1E04  9.206*1E04  1.079*1E03  7.043*1E04  −9.717*1E05  1.310*1E03  1.435*1E03  7.850*1E04  9.62*1E04 
20.30  4.453*1E04  6.589*1E04  3.739*1E04  5.820*1E04  1.353*1E03  4.033*1E04  8.078*1E04  7.510*1E04  5.616*1E04  1.066*1E03  2.823*1E04  2.136*1E04  5.72*1E04 
Optimal beta values for spiral lowpass geometry (SSIM optimization)
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  2.929*1E04  4.407*1E04  5.365*1E04  3.320*1E04  5.587*1E04  4.524*1E04  4.398*1E04  4.448*1E04  1.932*1E04  8.479*1E04  4.432*1E04  1.973*1E04  4.42*1E04 
48.50  4.101*1E04  3.171*1E04  2.837*1E04  3.176*1E04  5.204*1E04  4.530*1E04  4.610*1E04  5.018*1E04  4.352*1E04  9.829*1E04  6.130*1E04  2.600*1E04  4.44*1E04 
40.60  5.247*1E04  2.430*1E04  2.495*1E04  4.220*1E04  5.691*1E04  2.830*1E04  7.014*1E04  3.369*1E04  7.813*1E04  1.039*1E04  3.052*1E04  1.489*1E04  3.21*1E04 
32.80  3.198*1E05  3.418*1E06  2.612*1E04  2.844*1E04  2.844*1E04  3.682*1E04  3.340*1E04  2.250*1E04  2.093*1E05  5.343*1E05  1.320*1E04  1.717*1E04  1.98*1E04 
28.10  9.243*1E06  7.813*1E06  4.480*1E05  7.660*1E06  3.377*1E05  3.461*1E04  3.067*1E04  1.434*1E05  1.066*1E05  1.429*1E05  7.380*1E06  5.000*1E06  1.25*1E05 
26.60  9.680*1E06  9.888*1E06  5.209*1E05  7.324*1E06  4.948*1E05  6.414*1E04  9.766*1E06  1.607*1E05  6.538*1E06  1.148*1E05  5.924*1E06  1.136*1E04  1.07*1E05 
25.00  7.717*1E06  7.813*1E06  9.662*1E05  9.750*1E06  1.960*1E04  2.863*1E04  7.765*1E05  8.038*1E05  5.798*1E06  1.576*1E05  6.060*1E06  9.034*1E05  4.67*1E05 
21.90  1.051*1E05  9.750*1E06  1.438*1E05  5.714*1E06  2.244*1E04  1.726*1E04  1.944*1E04  7.759*1E05  7.979*1E05  1.394*1E05  1.187*1E05  4.749*1E05  3.09*1E05 
20.30  8.590*1E06  7.473*1E06  7.542*1E06  3.784*1E06  3.100*1E05  8.866*1E05  1.953*1E06  1.166*1E04  1.164*1E04  5.748*1E05  1.123*1E05  9.008*1E05  2.11*1E05 
Optimal beta values for random sampling on PDF geometry (SSIM optimization)
Sample rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  3.879*1E04  4.778*1E04  3.463*1E04  2.347*1E04  2.575*1E04  5.832*1E04  3.699*1E04  2.552*1E04  6.620*1E04  6.740*1E04  2.959*1E04  3.955*1E04  3.79*1E04 
48.50  2.957*1E03  4.742*1E03  4.121*1E03  1.513*1E03  5.078*1E04  3.129*1E03  2.873*1E03  3.145*1E03  2.784*1E03  1.926*1E03  2.829*1E03  3.197*1E03  2.92*1E03 
40.60  2.834*1E03  3.961*1E03  3.234*1E03  4.632*1E03  4.970*1E04  4.112*1E04  2.602*1E03  3.376*1E03  1.889*1E03  3.994*1E03  2.581*1E03  4.480*1E03  3.03*1E03 
32.80  2.285*1E03  2.578*1E04  2.069*1E03  2.080*1E03  2.182*1E03  2.415*1E03  2.275*1E03  2.245*1E03  3.317*1E03  3.334*1E03  1.757*1E03  2.649*1E03  2.26*1E03 
28.10  4.355*1E03  4.886*1E04  6.020*1E03  2.377*1E03  3.631*1E03  7.985*1E04  2.793*1E03  2.806*1E03  4.510*1E03  3.320*1E03  3.175*1E03  4.770*1E03  3.25*1E03 
26.60  2.809*1E03  2.626*1E03  3.631*1E03  1.898*1E03  1.739*1E03  5.134*1E04  1.089*1E03  1.712*1E03  3.226*1E03  1.510*1E03  1.889*1E03  2.007*1E03  1.89*1E03 
25.00  2.474*1E03  9.774*1E04  4.994*1E03  6.392*1E04  1.381*1E03  3.023*1E03  4.750*1E03  2.253*1E03  7.196*1E04  2.347*1E03  1.936*1E03  3.353*1E03  2.30*1E03 
21.90  2.254*1E03  2.776*1E03  2.575*1E03  5.158*1E03  2.319*1E03  2.012*1E03  1.015*1E03  3.837*1E04  3.279*1E03  1.706*1E03  2.336*1E03  2.642*1E03  2.33*1E03 
20.30  2.022*1E03  4.004*1E03  3.867*1E04  3.612*1E03  4.730*1E03  2.663*1E03  1.801*1E03  2.231*1E03  5.073*1E04  3.918*1E04  2.242*1E03  3.168*1E03  2.24*1E03 
Optimal beta values for random phase encoding (SSIM optimization)
Sample Rate  Image1  Image2  Image3  Image4  Image5  Image6  Image7  Image8  Image9  Image 10  Image11  Image 12  Median 

62.50  7.155*1E04  7.113*1E04  5.765*1E04  8.666*1E04  5.501*1E04  −1.006*1E03  −6.614*1E04  −2.959*1E04  3.854*1E05  5.893*1E04  0.000  6.427*1E04  5.63*1E04 
48.50  1.199*1E05  1.786*1E03  2.434*1E04  4.558*1E04  1.400*1E04  1.468*1E03  1.243*1E03  1.103*1E03  2.851*1E04  1.590*1E03  0.000  1.707*1E03  7.79*1E04 
40.60  2.279*1E03  8.726*1E04  9.424*1E04  5.586*1E04  4.216*1E04  −3.205*1E04  −4.021*1E04  1.538*1E03  1.124*1E03  7.060*1E04  1.046*1E03  1.522*1E03  9.08*1E04 
32.80  2.561*1E04  4.782*1E03  5.028*1E05  8.742*1E04  5.847*1E04  −7.650*1E04  −1.547*1E04  4.421*1E04  2.294*1E04  3.368*1E03  2.354*1E04  4.703*1E03  3.49*1E04 
28.10  5.722*1E04  2.451*1E04  1.896*1E04  9.717*1E05  1.935*1E04  8.850*1E05  9.451*1E04  −7.617*1E05  5.061*1E05  6.785*1E05  8.459*1E05  4.883*1E07  9.28*1E05 
26.60  9.766*1E06  7.360*1E04  1.344*1E04  2.263*1E04  3.434*1E04  5.932*1E04  1.190*1E03  5.092*1E04  2.915*1E05  9.766*1E06  −4.639*1E06  3.732*1E04  2.85*1E04 
25.00  1.500*1E03  3.746*1E04  1.015*1E03  8.685*1E04  6.149*1E04  1.737*1E03  1.722*1E03  1.105*1E03  3.997*1E04  5.024*1E04  2.420*1E04  8.242*1E04  8.46*1E04 
21.90  7.324*1E06  1.136*1E03  −3.906*1E06  6.226*1E04  1.953*1E06  7.024*1E04  −3.088*1E03  1.023*1E03  7.736*1E06  −1.375*1E04  2.603*1E03  1.717*1E03  3.15*1E04 
20.30  −4.746*1E05  8.667*1E06  3.418*1E05  4.228*1E05  8.818*1E05  −9.306*1E06  −9.048*1E06  5.127*1E06  6.141*1E06  −4.362*1E04  2.938*1E05  −3.554*1E04  5.63*1E06 
For each sampling geometry, we provide the median α, βvalues in the last column of Tables 7, 8, 9, 10, 11, 12, 13, 14. The recommendation to use the median α, βvalues for reconstruction was validated using leaveoneout. Here, note that the mean SSIM values reported in Tables 3, 4, 5, 6 are based on performing the image reconstruction using the median values from the training set of the remaining 11 image pairs. There is thus a mismatch between the median α, β values and the optimal values that can only be estimated when the entire image is available. Yet, despite this mismatch, we find that the second approach worked very well. For example, all of the SLP reconstructions for 20.3% sampling gave average SSIM index values above 0.93.
We next discuss the relationship between PSNR and mean SSIM values. Numerically, an SSIM value of 1 corresponds to an infinite value for the PSNR. Furthermore, an SSIM value of 0 corresponds to a negative infinity value for PSNR. Ideally, an increase in the SSIM index will also correspond to an increase in PSNR. We have verified this fact for the meanSSIM values of Table 4 and their corresponding PSNR values. On the other hand, this property does not hold for individual images. For example, image #2 for sampling at 20.30% has a PSNR value that is above the median, while its SSIM value is below the median. In such cases, the SSIM trend is considered to be superior to the PSNR. Overall though, all of the SLP reconstructions gave mean SSIM values above 0.93. Thus, the SLP geometry reconstructions are considered to be excellent.
Discussion
We have found that parameter optimization achieved significant image quality improvements using a relatively small number of iterations. Typically, five iterations were required to improve PSNR by over 10 dB from the initial values of α=0, β=0 and achieve reconstructions that were within 5% of the optimal quality value.
We take a closer look at the α and β ranges for select geometries from the DPE and SLP classes. We selected both the DPE and SLP geometries with the least number of samples that still provided reasonable performance regions in the optimal quality regions for 30 dB, 35 dB, 40 dB, and 45 dB. Examining the range of the parameters provides a more detailed description of the parameter space, and which constraint in (3) is essential for acceptable reconstructions. The minimum and maximum parameter values that define the vertices of the inner and outer bounding boxes are listed in Table 2.
Three instances in the above table stand out because the minimum alpha and/or beta parameters are negative. These values are a byproduct of the interpolation and intersection method used to calculate the parameter ranges. When negative parameter values were inserted into the optimization algorithm, the reconstruction algorithm breaks down. In these cases, we replace the negative value with zero to achieve an acceptable result.
It is imperative that our discussion turn to the effect of a zero parameter value for either constraint in the reconstruction problem. Typically, the effect of α=0 is a more prominent presence of highfrequency errors in the spatial domain, as well as a higher amount of artifacts from the sampling geometry. Conversely, the effect of β=0 is a loss of highfrequency spatial components in the spatial domain reconstruction. In this case, the TVnorm tends to drive the solution to result whose finite difference in each dimension is minimized. While removing one of the penalty terms in the CS objective function may result in a usable reconstruction, our observations indicated that this case will not be optimal.
When both parameters are zero, the proposed optimization framework breaks down. It is possible to achieve acceptable image quality when one of the two parameters is zero, but not both. While completely removing one of the penalty terms in (3) may result in acceptable reconstructions in terms of PSNR, an increased PSNR can be achieved by using the (nonzero) optimal values.
For activity detection, we have found that optimization of the difference image led to excellent results. In particular, the SLP geometry gave excellent average SSIM values of 0.9747 (max = 1). The use of the leaveoneout validation with mean SSIM optimization also verified the fact that the SLP family of geometries gives the best result. Furthermore, recommended parameter values for the SLP geometries are given by the median values of Tables 8 and 12.
Conclusions
In conclusion, we have found that CS parameter optimization can dramatically improve fMRI image reconstruction quality. Furthermore, deterministic SLP scanning geometries (similar to the 3D geometries of [4]) consistently gave the best image reconstruction results. The implication of this result is that less complex sampling geometries will suffice over random sampling, provided the TVNorm and Transform penalty parameters are selected from the range of values we have calculated in our proposed methodology. We have also found that we can obtain stable parameter regions that can be used to achieve specific levels of image reconstruction quality when combined with specific kspace sampling geometries. More importantly, median parameter values can be used for obtaining excellent reconstructions. This observation has been validated using leaveoneout and optimization based on the mean SSIM index.
Declarations
Acknowledgements
None.
Authors’ Affiliations
References
 Candès E, Romberg J, Tao T: Robust uncertainty principles: exact signal recovery from incomplete frequency information. IEEE Trans. Information Theory 2006, 52: 489–509.View ArticleGoogle Scholar
 Chan T, Esedoglu S, Park F, Yip A: Recent developments in total variation image restoration. Springer Verlag: Mathematical Models of Computer Vision; 2005.Google Scholar
 Wang Z, Bovik AC, Sheikh HR, Simoncelli EP: Image quality assessment from error visibility to structural similarity. IEEE Transactions on Image Processing 2004, 13: 600–612. 10.1109/TIP.2003.819861View ArticleGoogle Scholar
 Lindquist MA, Zhang C, Glover G, Schepp L: Rapid threedimensional functional magnetic resonance imaging of the initial negative BOLD response. J of Magnetic Resonance 2008, 191: 100–111. 10.1016/j.jmr.2007.12.016View ArticleGoogle Scholar
 Hansen PC: The Lcurve and its use in the numerical treatment of inverse problems. 2012. Available Retrieved http://www.unimed.sintef.no/project/eVITAmeeting/2005/Lcurve.pdf Google Scholar
 Girod B: What's wrong with mean squared error? In Visual Factors of Electronic Image Communications. Edited by: Watson AB. MIT Press; 1993:207–220.Google Scholar
 Huber PJ: The 1972 Wald Lecture: Robust Statistics: A Review. Ann Math Stat 1972, 43(4):1041–1067. 10.1214/aoms/1177692459View ArticleGoogle Scholar
 Friman O, Cedefamn J, Lundberg P, Borga M, Knutsson H: Detection of neural activity in functional MRI using canonical correlation analysis. Magn Reson Med 2001, 45: 323–330. 10.1002/15222594(200102)45:2<323::AIDMRM1041>3.0.CO;2#View ArticleGoogle Scholar
 Nan F, Nowak R: Generalized likelihood ratio detection for fMRI using complex data. IEEE Transactions on Med Imaging 1999, 18(4):320–329. Available Retrieved http://dsp.rice.edu/sites/dsp.rice.edu/files/publications/journalarticle/1999/generalizetmi1999.pdf 10.1109/42.768841View ArticleGoogle Scholar
 Berstein MA, Thomasson DM, Perman WH: Improved detectability in low signaltonoise ration magnetic resonance images by means of a phasecorrected real reconstruction. Med Phys 1989, 16: 813–817. 10.1118/1.596304View ArticleGoogle Scholar
 Frank LR, Buxton RB, Wong EC: Probabilistic analysis and functional magnetic resonance imaging data. Med Phys 1998, 39: 132–148.Google Scholar
 Lai S, Glover GH: Detection of BOLD fMRI signals using complex data. Proc. Int. Society Magnetic Resonance Medicine Meeting 1997. 1671.Google Scholar
 Sekihara K, Koizumi H: Detecting cortical activities from fMRI timecourse data using MUSIC algorithm with forward and backward covariance averaging. Magn Reson Med 1996, 35(6):807–813. 10.1002/mrm.1910350604View ArticleGoogle Scholar
 Siewert B, Bly BM, Schlaug G, Darby DG, Thangaraj V, Warach S, Edelman R: Comparison of the bold and epistar technique for functional brain imaging using signal detection theory. Magn Reson Med 1996, 36: 249–255. 10.1002/mrm.1910360212View ArticleGoogle Scholar
 Strother SC, Kanno I, Rottenberg DA: Principal component analysis, variance portioning, and functional connectivity. J Cerebral Blood Flow Metabol 1995, 15: 353–360. 10.1038/jcbfm.1995.44View ArticleGoogle Scholar
 Lusting M, Donoho D, Pauly JM: Sparse MRI: the application of compressed sensing for rapid MR Imaging. Magnetic Resonance in Medicine 2007, 58(6):1182–1195. Available http://groups.bme.gatech.edu/groups/cfmg/group/JC_Fall2010/Chris.pdf Retrieved. Associated code available at http://www.stanford.edu/~mlustig/SparseMRI.html, retrieved April 9th, 2012. 10.1002/mrm.21391View ArticleGoogle Scholar
 Gamper U, Boesiger P, Kozerke S: Compressed sensing in dynamic MRI. Magn Reson Med 2008, 59: 365–373. 10.1002/mrm.21477View ArticleGoogle Scholar
 Jung H, Ye JC, Kim EY: Improved kt BLAST and kt SENSE using FOCUSS. Phys Med Biol 2007, 52: 3201–3226. 10.1088/00319155/52/11/018View ArticleGoogle Scholar
 Haldar JP, Hernando D, Liang Z: CompressedSensing MRI with Random Encoding. IEEE Transactions on Medical Imaging 30:893–903, 2011. 10.1109/TMI.2010.2085084View ArticleGoogle Scholar
 Liu B, Zou YM, Ying L: Sparsesense: an application of compressed sensing in parallel MRI. Proc. of the 2008 Int Conf on Technology and Information in Biomedicine 2008, 127–130.Google Scholar
 Trzasko J, Manduca A, Borisch E: Sparse MRI reconstruction via multiscale L0continuation. Proc. of the IEEE/SP 14th Workshop on Statistical Signal Processing. 2007, 176–180.Google Scholar
 Greiser A, von Kienlin M: Efficient kspace sampling by densityweighted phaseencoding. Magn Reson Med 2003, 50: 1266–1275. 10.1002/mrm.10647View ArticleGoogle Scholar
 Peters DC, Korosec FR, Grist TM, Block WM, Holden JE, Vigen KK, Mistretta CA: Undersampled projection reconstruction applied to MR angiography. Magn Reson Med 2000, 43: 91–101. 10.1002/(SICI)15222594(200001)43:1<91::AIDMRM11>3.0.CO;24View ArticleGoogle Scholar
 McGibney G, Smith MR, Nichols ST, Crawley A: A Quantitative evaluation of several partial Fourier reconstruction algorithms used in MRI. Magn Reson Med 1993, 30: 51–59. 10.1002/mrm.1910300109View ArticleGoogle Scholar
 Pruessman KP, Weiger M, Scheidegger MB, Boesiger P: SENSE: sensitivity encoding for fast MRI. Magn Reson Med 1999, 42: 952–962. 10.1002/(SICI)15222594(199911)42:5<952::AIDMRM16>3.0.CO;2SView ArticleGoogle Scholar
 Sodickson DK, Manning WJ: Simultaneous acquisition of spatial harmonics (SMASH): fast imaging with radiofrequency coil arrays. Magn Reson Med 1997, 38: 591–603. 10.1002/mrm.1910380414View ArticleGoogle Scholar
 Gabr RE, Aksit P, Bottomley PA, Youssef AM, Kadah YM: DeconvolutionInterpolation Gridding (DING): accurate Reconstruction for Arbitrary kSpace Trajectories. Magn Reson Med 2006, 56(6):1182–1191. 10.1002/mrm.21095View ArticleGoogle Scholar
 Donoho D: Compressed sensing. IEEE Trans Inform Theory 2006, 52(4):1289–1306.MathSciNetView ArticleGoogle Scholar
 Nelder JA, Mead R: A simplex method for function minimization. J Comput 1965, 7: 308–313. Available Retrieved http://comjnl.oxfordjournals.org/content/7/4/308.full.pdf+html 10.1093/comjnl/7.4.308View ArticleGoogle Scholar
 Powell MJD: Direct search algorithms for optimization calculations. In Acta Numerica. Edited by: Iserles A. Cambridge, UK: Cambridge University Press; 1998:287–336.Google Scholar
Copyright
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.