The impact of temporal regularization on estimates of the BOLD hemodynamic response function: a comparative analysis

In fMRI data analysis it has been shown that for a wide range of situations the hemodynamic response function (HRF) can be reasonably characterized as the impulse response function of a linear and time invariant system. An accurate and robust extraction of the HRF is essential to infer quantitative information about the relative timing of the neuronal events in different brain regions. When no assumptions are made about the HRF shape, it is most commonly estimated using time windowed averaging or a least squares estimated general linear model based on either Fourier or delta basis functions. Recently, regularization methods have been employed to increase the estimation efficiency of the HRF; typically these methods produce more accurate HRF estimates than the least squares approach [Goutte, C., Nielsen, F.A., Hansen, L.K., 2000. Modeling the Haemodynamic Response in fMRI Using Smooth FIR Filters. IEEE Trans. Med. Imag. 19(12), 1188-1201.]. Here, we use simulations to clarify the relative merit of temporal regularization based methods compared to the least squares methods with respect to the accuracy of estimating certain characteristics of the HRF such as time to peak (TTP), height (HR) and width (W) of the response. We implemented a Bayesian approach proposed by Marrelec et al. [Marrelec, G., Benali, H., Ciuciu, P., Pelegrini-Issac, M., Poline, J.-B., 2003. Robust Estimation of the Hemodynamic Response Function in Event-Related BOLD fMRI Using Basic Physiological Information. Hum. Brain Mapp. 19, 1-17., Marrelec, G., Benali, H., Ciuciu, P., Poline, J.B. Bayesian estimation of the hemodynamic of the hemodynamic response function in functional MRI. In: R. F, editor; 2001; Melville. p 229-247.] and its deterministic counterpart based on a combination of Tikhonov regularization [Tikhonov, A.N., Arsenin, V.Y., 1977. Solution of ill-posed problems. Washington DC: W.H. Winston.] and generalized cross-validation (GCV) [Wahba, G., 1990. Spline Models for Observational Data. Philadelphia: SIAM.] for selecting the regularization parameter. The performance of both methods is compared with least square estimates as a function of temporal resolution, color and strength of the noise, and the type of stimulus sequences used. In almost all situations, under the considered assumptions (e.g. linearity, time invariance and smooth HRF), the regularization-based techniques more accurately characterize the HRF compared to the least-squares method. Our results clarify the effects of temporal resolution, noise color, and experimental design on the accuracy of HRF estimation.