Restoration of fMRI Signal using Wiener Filter in a Wavelet Domain
D Asefa, D Mital, S Haque, S Srinivasan
Keywords
activation points, fmri, noise, wavelet, wiener filter
Citation
D Asefa, D Mital, S Haque, S Srinivasan. Restoration of fMRI Signal using Wiener Filter in a Wavelet Domain. The Internet Journal of Medical Informatics. 2005 Volume 3 Number 1.
Abstract
Waveletbased denoising methods have the advantage over lowpass filtering in that relevant detail information is retained, while small details due to noise are discarded. This paper reports a novel technique of removing noise, using a wiener filter in wavelet domain, from an fMRI data and display selectively event related voxels in a spatial domain.
After the fMRI signal is wavelet transformed to temporal domain, its median and mean become symmetric as described in Section and its median of absolute deviation (MAD) is calibrated with the standard deviation of a Gaussian distribution to approximately determine the standard deviation of the noise in the fMRI signal. Once the noise level is determined, the result is used to determine the power spectrum of the coherent signal and the corresponding standard variable. The whole process is used to accurately estimate the power of the coherent signal and the associated noise at a given voxel location. The estimated power spectrum is used to approximate the optimal wiener filter coefficient for removing the noise from the fMRI signal.
In this experiment, only the signal to be restored is available and all prior knowledge about the ideal signal has to be estimated from it. Though an fMRI is contaminated by both Gaussian and Rician distribution related noises, it is safe to assume noise in an fMRI signal to be additive white Gaussian noise. In this paper, we showed that the power spectrum estimated from this single copy of degraded signal is a true power spectrum of the signal, and as a result, the restoration filter or noise removing filter is optimal, though there is a lack of accurate prior information. The method successfully removes noises and exposes activated voxels in fMRI signal all the time as shown in this paper.
Introduction
Each activity a person performs is managed by a certain location of the brain. The location of a brain directly related to an activity can be visualized using images from functional magnetic resonance imaging (fMRI) instrument. These images are obtained using the changes between active and nonactive state of location of a brain. The image contrast obtained this way is very small, and fMRI instrument is so sensitive that it picks unwanted signals or noise, that induce distortions of the actual experimental signals, which can be interpreted as false brain activities. The objective of this paper is to remove noise or distortions, which are unrelated to the experimental fMRI signals. The methods, both periodic discrete and translation invariant packet wavelets are suggested to remove noise from fMRI data due to their decorrelating effects in a temporal domain.
The main sources of noise are not fully understood. A number of possible sources have been suggested, for example, slow phase variations in the MR images due to respiration movements, cardiac and other physiological processes, patient movement, and local changes in the magnetic field due to scanner instabilities. During fMRI image acquisition process, high frequency components like heart rate (0.6 –1.2 Hz) and respiration (0.1 – 0.5 Hz) are under sampled with typical repeat times (TRs ) of 3 to 7s and can, according to Nyquist's theorem, be expressed as lowfrequency (0.1 Hz) signal components or aliased higher frequency signals [_{1}].
Due to the fact that Rician distribution is used in the physics of magnetic resonance and the Gaussian distribution in functional neuroimaging, the noise, in the blood oxygenation dependent (BOLD) response of an fMRI data, is both Gaussian and Rician distributed. For high SNR fMRI data, Rician distributed noise is symmetric, thus, can be considered as Gaussian distributed. For low SNR fMRI data, there is a difference between Gaussian and Rician distributed noise, i.e, an image with low intensity and Rician distributed has probability density which is asymmetric. It was further shown that the difference between two Rician distributed images is symmetric and Gaussian distributed. Since an fMRI BOLD response is the difference between two BOLD responses, the active and passive state, thus, the noise in the resulting BOLD response is symmetrically distributed, which is characteristic of Gaussian distribution.
To remove noise from an fMRI signal, general wavelet method is presented. This method decorrelates and extracts noise from a BOLD response leaving a clean event related BOLD response of the fMRI data. To remove the noise, first, the noise is decorrelated and then the intensity of the noise is estimated statistically. Second, the estimated noise is removed from the fMRI data in a temporal domain. The wavelet methods, both periodic and packet based wavelets, a wiener filter is used with coefficients estimated from the noisy fMRI data adaptively. A wavelet transform approximately acts as a detrender, which facilitates the removal of noise from fMRI data in a temporal domain without affecting an event related BOLD response. Real and simulated fMRI data is used to confirm that the methods remove noise and preserve activities related BOLD response of an fMRI data.
After denoising, to reduce the false discovery rate, clustering is performed in spatial domain based on optimal minimal cluster size obtained from Monte Carlo simulation for a userdefined confidence level. In this work, it is shown that the suggested methods can recover a signal from a data with a noise that has intensity of any standard deviation away from the center of the data but an optimal performance is obtained when the noise level is less than seven standard deviations from the center.
Modeling of fMRI data
Overall, we need to model and remove deterministic components from the time series before proceeding with the statistical analysis. The fMRI signal
Where
where
Background
Statistical Parametric Mapping (SPM) tool and many research works use Gaussian filtering to remove nosie from an fMRI data. Gaussian filtering is lowpass filtering, that is being used traditionally to process fMRI data, but it can remove relevant detail information. Gaussian filtering requires a window or a kernel size to be derived from the fMRI data itself to avoid processing errors . The other most widely used methods are Fast Fourier Transform (FFT) and SPline, but FFT unable to identify sharp transient events that are similar to fMRI data,. SPline is FFT and wavelet based, its sharp frequency characterization makes a good fit to process an fMRI or time related signal[_{2}] but it needs optimization [_{3}].
Most of the existing literatures suggest performing both the analysis and filtering in the time domain or temporal domain, and after the analysis phase in the time domain, an inverse transform is applied to reconstruct an activation map from the coefficients that are designated as significant. While this reconstructed map is very useful for visualization purposes, it does not have a direct statistical interpretation, that is, the statistical parameters, such as, t or z values are computed in the time domain and there is no straightforward way to map the statistics to the spatial domain.
In this paper both statistical analysis and filtering are performed separately, this approach helps optimize the detection of false positive active points of the voxel time series.
Methodology – Estimation of Noise and its removal from fMRI
For a continuous variate
This means that one half of the observations are below and the other half above the median. For a sample of ordered variables
For a normal distribution, the probability density function is
As for any symmetrical distribution, mean and median coincide, thus
The DICOM data with activation points is collected form patients and transformed into AFNI [_{4}] format. The signal from AFNI is processed as specified below and ported to AFNI for visualization and wavelet processed to determine the activation points in the brain and to compare and contrast the processed signal with othe fMRI processing methods. The wavelet filter used is symmlet with size four(4). The procedures employed to separate the coherent signal form the noise are as follows:

For processing 1D fMRI data :Divide the signal into 1D representations or into 1D timeseries; if the 1D signal is not power of two, padding with zeros several hundred positions is necessary

For Processing 2D fMRI data: Divide the 3D fMRI data into slices and if the size of a slice is not power of two, padding with zeros several hundred positions is necessary

Wavelet decompose the fMRI signal, this step makes the noise distribution N(0, σ2)


Decompose the fMRI data into discrete wavelet coefficients (wi)

Decompose the fMRI data into packet wavelet coefficients (wi)


Estimation of the noise level σ from its distribution N(0, σ2): σ = MAD(wi) / 0.6745 for both periodic discrete and packet wavelets, separately

Select the representative coefficients  wj form the wavelet packet coefficients  wi of step 2b; to select the representative coefficients, Shannon entropy is used . The selection process involves the following steps :


Start from the bottom of the decomposition tree and mark every thing

Determine the entropy of each element of the tree, and then if the entropy of a parent node is less than or equal to the sum of entropy of its the two children, them unmark the children and mark the parent, otherwise leave the children marked, continue this way until you reach the top. Shanon entropy is determined using the formula E(X) = Σnpnlog2pn), where pn = (xn)2/ x2 , where xn is the nth component of x, x is the signal length and E(X) is the entropy of the fMRI signal x.


Normalize the noise in the temporal domain to make the noise distribution N(0,1) : wi = wi/σ for discrete wavelet and for packet wavelet wj = wj / σ, and σ is a result from discrete or packet wavelets used in step 4

Determine the coefficient of wiener filtering of the wavelet coefficient: Φ = wi 2 / (wi 2 + σ2) for discrete wavelet and for packet wavelet Φ = wj 2 / (wj 2 + σ2)

Get the wavelet coefficients from the noisy fMRI time series and or slice from step 3, and then apply wiener filter from step 7 to determine the clean wavelet coefficients: wi clean = wi⋅Φ

Reconstruct the clean wavelet coefficient (wi clean ) into an fMRI time series or slice

Go to step 3 and repeat the steps until the entire time series or fMRI slices are processed.
Results
The brain processes different information in different way, for example, when the subject was instructed to process items according to their meanings ( is the word hot or cold ?) or when the subject is instructed to perform bilateral finger tapping at a set of interval of time. The former task was associated with activations in a set of brain regions including left lateral prefrontal cortex (PFC) and left medial temporal cortex. The latter showed relatively greater activation in right and left PFC.
The mapping of functions to brain regions presupposes, among other things, a thorough understanding of the cognitive functions that are to be mapped onto the brain regions. This understanding, however, while developing, is still rudimentary. Someone with schizophrenia, when confronted with a psychological task, might tackle it in a very different way, in terms of the cognitive strategies used, and from a healthy person confronted with the same task. The observation that brain activity differs across the two individuals would only be interpretable insofar as one thoroughly understood the processes that each individual invoked in response to the task demands.
Recent studies have found that the spatial extent of fMRI activation in healthy older adults is approximately half that of younger subjects, for both visual and task related activities. Also, it is very important to note this fact: the difference, in activation, between healthy older and young subject, is not associated with reduced hemodynamic response (HDR) amplitude in older subjects, because young and elderly adults have similar HDR amplitudes [_{5}], and have similar distributions of HDR amplitudes across voxels. Furthermore, greater head motion in the elderly does not cause this difference, although head motion differences can affect spatial extent of activation [_{5}]. Instead, spatial extent differences are associated with higher voxelwise noise levels in elderly adults, perhaps due to changes in cardiac or respiratory effects upon the fMRI signal .
In this experiment no comparison across subject is performed and an individual subject data is analyzed independently, and the subject was instructed to perform bilateral finger tapping as brain activation mapping of this action can be seen on both sides of lateral prefrontal cortex, and this finding, Figure ,is consistent with theoretical basis of neural functions for a healthy individual, the complete result of the processing an fMRI obtained from a healthy person is shown in Figure.
Figure 4
The other most important thing to notice is, during the experiment, sometimes the well instructed subjects do all kinds of things in addition to the task they were asked to perform and this showed up in the fMRI images, as seen in Figure, and interpretation of the results should take into consideration this kind of events. Furthermore, the best strategy is to mask the data to include only the areas or regions that corresponding to the activities in the brain.
Acknowledgements
For the data used in this research work, we would like to acknowledge Prof. B. Biswal, director of the fMRI research lab. of UMDNJ and Rutgers Universities. For their encouragements, I will like to acknowledge Prof. Melkamu Zeleke chairman of Mathematics dept. at Paterson University New Jersey , Prof. Wolde Woubneh of Kean University New Jersey, Associate Prof . Abebe Kebede of A&T University –North Carolina.