  Several problems in neuroimaging and beyond require inference on the
parameters of multi-task sparse hierarchical regression models. Examples
include M/EEG inverse problems, neural encoding models for task-based fMRI
analyses, and climate science. In these domains, both the model parameters to
be inferred and the measurement noise may exhibit a complex spatio-temporal
structure. Existing work either neglects the temporal structure or leads to
computationally demanding inference schemes. Overcoming these limitations, we
devise a novel flexible hierarchical Bayesian framework within which the
spatio-temporal dynamics of model parameters and noise are modeled to have
Kronecker product covariance structure. Inference in our framework is based on
majorization-minimization optimization and has guaranteed convergence
properties. Our highly efficient algorithms exploit the intrinsic Riemannian
geometry of temporal autocovariance matrices. For stationary dynamics described
by Toeplitz matrices, the theory of circulant embeddings is employed. We prove
convex bounding properties and derive update rules of the resulting algorithms.
On both synthetic and real neural data from M/EEG, we demonstrate that our
methods lead to improved performance.
