  The generalized extreme value (GEV) distribution is a popular model for
analyzing and forecasting extreme weather data. To increase prediction
accuracy, spatial information is often pooled via a latent Gaussian process
(GP) on the GEV parameters. Inference for GEV-GP models is typically carried
out using Markov chain Monte Carlo (MCMC) methods, or using approximate
inference methods such as the integrated nested Laplace approximation (INLA).
However, MCMC becomes prohibitively slow as the number of spatial locations
increases, whereas INLA is only applicable in practice to a limited subset of
GEV-GP models. In this paper, we revisit the original Laplace approximation for
fitting spatial GEV models. In combination with a popular sparsity-inducing
spatial covariance approximation technique, we show through simulations that
our approach accurately estimates the Bayesian predictive distribution of
extreme weather events, is scalable to several thousand spatial locations, and
is several orders of magnitude faster than MCMC. A case study in forecasting
extreme snowfall across Canada is presented.
