  Background: Inference from fMRI data faces the challenge that the hemodynamic
system that relates neural activity to the observed BOLD fMRI signal is
unknown.
  New Method: We propose a new Bayesian model for task fMRI data with the
following features: (i) joint estimation of brain activity and the underlying
hemodynamics, (ii) the hemodynamics is modeled nonparametrically with a
Gaussian process (GP) prior guided by physiological information and (iii) the
predicted BOLD is not necessarily generated by a linear time-invariant (LTI)
system. We place a GP prior directly on the predicted BOLD response, rather
than on the hemodynamic response function as in previous literature. This
allows us to incorporate physiological information via the GP prior mean in a
flexible way, and simultaneously gives us the nonparametric flexibility of the
GP.
  Results: Results on simulated data show that the proposed model is able to
discriminate between active and non-active voxels also when the GP prior
deviates from the true hemodynamics. Our model finds time varying dynamics when
applied to real fMRI data.
  Comparison with Existing Method(s): The proposed model is better at detecting
activity in simulated data than standard models, without inflating the false
positive rate. When applied to real fMRI data, our GP model in several cases
finds brain activity where previously proposed LTI models does not.
  Conclusions: We have proposed a new non-linear model for the hemodynamics in
task fMRI, that is able to detect active voxels, and gives the opportunity to
ask new kinds of questions related to hemodynamics.
