  Multivariate spatially-oriented data sets are prevalent in the environmental
and physical sciences. Scientists seek to jointly model multiple variables,
each indexed by a spatial location, to capture any underlying spatial
association for each variable and associations among the different dependent
variables. Multivariate latent spatial process models have proved effective in
driving statistical inference and rendering better predictive inference at
arbitrary locations for the spatial process. High-dimensional multivariate
spatial data, which is the theme of this article, refers to data sets where the
number of spatial locations and the number of spatially dependent variables is
very large. The field has witnessed substantial developments in scalable models
for univariate spatial processes, but such methods for multivariate spatial
processes, especially when the number of outcomes are moderately large, are
limited in comparison. Here, we extend scalable modeling strategies for a
single process to multivariate processes. We pursue Bayesian inference which is
attractive for full uncertainty quantification of the latent spatial process.
Our approach exploits distribution theory for the Matrix-Normal distribution,
which we use to construct scalable versions of a hierarchical linear model of
coregionalization (LMC) and spatial factor models that deliver inference over a
high-dimensional parameter space including the latent spatial process. We
illustrate the computational and inferential benefits of our algorithms over
competing methods using simulation studies and an analysis of a massive
vegetation index data set.
