Resting-state HRF Deconvolution

Under the assumption that the transformation from neural activation to BOLD response can be modeled as a linear $y(t)$ can be seen as the result of the convolution of neural states $s(t)$ with a HRF $h(t)$:
$$y(t)=s(t)\otimes h(t) + \epsilon(t)$$
Once calculated $h(t)$, we can obtain an approximation $\hat s(t)$ of the neural signal from the observed data using a Wiener filter
$$\widetilde s(t)=d(t)\otimes y(t)$$
Let $H(x)$, $Y(x)$, $E(x)$, and $D(x)$ be the Fourier transforms of $h(t)$, $y(t)$, $\epsilon(t)$, and $d(t)$, respectively.
The estimation $\widetilde s(t)$ of the neural states $s(t)$ is then given by
$$\widetilde s(t)=FT^{-1}\left(D(\omega)Y(\omega)\right)=FT^{-1}\left(\frac{H^*(\omega)Y(\omega)}{|H(\omega)|^2+|E(\omega)|^2}\right)
$$

BOLD fMRI parameters setting

1
2
3
4

data: BOLD-fMRI signal, nobs x nvar (nvar: number of variables; e.g. 200x90, 200x 50000, ....)
TR = 2; % 2seconds
thr = 1; % (mean+) thr*standard deviation threshold to detect event.
event_lag_max: maximum time from neural event to BOLD event in bins, not time (e.g. if we assume 10 seconds, and TR=2s, set the value to 10/2=5)

HRF Deconvolution

1

[data_deconv,event,HRF, adjust_global, PARA] = wgr_deconv_canonhrf_par(data,thr,event_lag_max,TR)

1
2
3
4
5

data_deconv: Deconvolved BOLD-fMRI signal
PARA: HRF parameters
PARA(1,:): response height (response magnitude of neuronal activity)
PARA(2,:): Time to peak (latency of neuronal activity)
PARA(3,:): Width / FWHM (duration of neuronal activity)

The code [download] uses the parallel for loop ‘parfor’. In case of older matlab versions, parfor can be changed to standard ‘for’.

This code uses canonical HRF with its delay and dispersion derivatives (basis function), as described in the paper [1]. The basis function can be modified to use FIR model.

It’s important that the input data are voxel containing some meaningful BOLD signal, so you can use a mask, even restricted to some specific area.

Here is a demo for HRF deconvolution (From ‘clean’ BOLD-fMRI NIFTI files to deconvolved NIFTI files), it assumed that you organize your data in following way:

The outputed files are organized in the following way:

Citation

Guo-Rong Wu, Wei Liao, Sebastiano Stramaglia, Ju-Rong Ding, Huafu Chen, Daniele Marinazzo. “A blind deconvolution approach to recover effective connectivity brain networks from resting state fMRI data.” Medical Image Analysis, 2013, 17:365-374. PDF
Guo-Rong Wu, Daniele Marinazzo. “Retrieving the Hemodynamic Response Function in resting state fMRI: methodology and applications.” PeerJ PrePrints, 2015. PDF