Spatial statistical modeling of arsenic accumulation in microsites of diverse soils

}, volume={411}, ISSN={["1872-6259"]}, DOI={10.1016/j.geoderma.2022.115697}, abstractNote={Determining reaction mechanisms that control the mobility of nutrients and toxic elements in soil matrices is confounded by complex assemblages of minerals, non-crystalline solids, organic matter, and biota. Our objective was to infer the chemical elements and solids that contribute to As binding in matrices of soil samples from different pedogenic environments at the micrometer spatial scale. Arsenic was reacted with and imaged in thin weathering coatings on eight quartz sand grains separated from soils of different drainage classes to vary contents of Fe and Al (hydr)oxides, organic carbon (OC), and other elements. The grains were analyzed using X-ray fluorescence microprobe (µ-XRF) imaging and microscale X-ray absorption near edge structure (μ-XANES) spectroscopy before and after treatment with 0.1 mM As(V) solution. Partial correlation analyses and regression models developed from multi-element µ-XRF signals collected across 100 × 100 µm2 areas of sand-grain coatings inferred augmenting effects of Fe, Zn, Ti, Mn, or Cu on As retention. Significant partial correlations (r′ > 0.11) between Fe and Al from time-of-flight secondary ion mass spectrometry (TOF-SIMS) analysis of most samples suggested that Fe and Al (hydr)oxides were partially co-localized at the microscale. Linear combination fitting (LCF) results for As K-edge μ-XANES spectra collected across grain coatings typically included >80% of As(V) adsorbed on goethite, along with varying proportions of standards of As(V) adsorbed on boehmite, As(V) or As(III) bound to Fe(III)-treated peat, and dimethylarsinic acid. Complementary fits for Fe K-edge μ-XANES spectra included ≥50% of the Fe(III)-treated peat standard for all samples, along with goethite. Our collective results inferred a dominance of Fe and possibly Al (hydr)oxides in controlling As immobilization, with variable contributions from Zn, Ti, Cu, or Mn, both across the coating of a single sand grain and between grains from soils developed under different pedogenic environments. Overall, these results highlight the extreme heterogeneity of soils on the microscale and have implications on soil management for mitigating the adverse environmental impacts of As.}, journal={GEODERMA}, author={Sharma, Aakriti and Guinness, Joseph and Muyskens, Amanda and Polizzotto, Matthew L. and Fuentes, Montserrat and Hesterberg, Dean}, year={2022}, month={Apr} } @article{lan_reich_guinness_bandyopadhyay_ma_moeller_2022, title={Geostatistical modeling of positive-definite matrices: An application to diffusion tensor imaging}, volume={78}, ISSN={["1541-0420"]}, DOI={10.1111/biom.13445}, abstractNote={Geostatistical modeling for continuous point-referenced data has extensively been applied to neuroimaging because it produces efficient and valid statistical inference. However, diffusion tensor imaging (DTI), a neuroimaging technique characterizing the brain's anatomical structure, produces a positive-definite (p.d.) matrix for each voxel. Currently, only a few geostatistical models for p.d. matrices have been proposed because introducing spatial dependence among p.d. matrices properly is challenging. In this paper, we use the spatial Wishart process, a spatial stochastic process (random field), where each p.d. matrix-variate random variable marginally follows a Wishart distribution, and spatial dependence between random matrices is induced by latent Gaussian processes. This process is valid on an uncountable collection of spatial locations and is almost-surely continuous, leading to a reasonable way of modeling spatial dependence. Motivated by a DTI data set of cocaine users, we propose a spatial matrix-variate regression model based on the spatial Wishart process. A problematic issue is that the spatial Wishart process has no closed-form density function. Hence, we propose an approximation method to obtain a feasible Cholesky decomposition model, which we show to be asymptotically equivalent to the spatial Wishart process model. A local likelihood approximation method is also applied to achieve fast computation. The simulation studies and real data application demonstrate that the Cholesky decomposition process model produces reliable inference and improved performance, compared to other methods.}, number={2}, journal={BIOMETRICS}, author={Lan, Zhou and Reich, Brian J. and Guinness, Joseph and Bandyopadhyay, Dipankar and Ma, Liangsuo and Moeller, F. Gerard}, year={2022}, month={Jun}, pages={548–559} } @article{roy_reich_guinness_shinohara_staicu_2021, title={Spatial Shrinkage Via the Product Independent Gaussian Process Prior}, volume={6}, ISSN={["1537-2715"]}, DOI={10.1080/10618600.2021.1923512}, abstractNote={We study the problem of sparse signal detection on a spatial domain. We propose a novel approach to model continuous signals that are sparse and piecewise-smooth as the product of independent Gaussian (PING) processes with a smooth covariance kernel. The smoothness of the PING process is ensured by the smoothness of the covariance kernels of the Gaussian components in the product, and sparsity is controlled by the number of components. The bivariate kurtosis of the PING process implies that more components in the product results in the thicker tail and sharper peak at zero. We develop an efficient computation algorithm based on spectral methods. The simulation results demonstrate superior estimation using the PING prior over Gaussian process prior for different image regressions. We apply our method to a longitudinal magnetic resonance imaging dataset to detect the regions that are affected by multiple sclerosis computation in this domain. Supplementary materials for this article are available online.}, journal={JOURNAL OF COMPUTATIONAL AND GRAPHICAL STATISTICS}, author={Roy, Arkaprava and Reich, Brian J. and Guinness, Joseph and Shinohara, Russell T. and Staicu, Ana-Maria}, year={2021}, month={Jun} } @article{brantley_guinness_chi_2020, title={BASELINE DRIFT ESTIMATION FOR AIR QUALITY DATA USING QUANTILE TREND FILTERING}, volume={14}, ISSN={["1932-6157"]}, DOI={10.1214/19-AOAS1318}, abstractNote={We address the problem of estimating smoothly varying baseline trends in time series data. This problem arises in a wide range of fields, including chemistry, macroeconomics and medicine; however, our study is motivated by the analysis of data from low cost air quality sensors. Our methods extend the quantile trend filtering framework to enable the estimation of multiple quantile trends simultaneously while ensuring that the quantiles do not cross. To handle the computational challenge posed by very long time series, we propose a parallelizable alternating direction method of multipliers (ADMM) algorithm. The ADMM algorthim enables the estimation of trends in a piecewise manner, both reducing the computation time and extending the limits of the method to larger data sizes. We also address smoothing parameter selection and propose a modified criterion based on the extended Bayesian information criterion. Through simulation studies and our motivating application to low cost air quality sensor data, we demonstrate that our model provides better quantile trend estimates than existing methods and improves signal classification of low-cost air quality sensor output.}, number={2}, journal={ANNALS OF APPLIED STATISTICS}, author={Brantley, Halley L. and Guinness, Joseph and Chi, Eric C.}, year={2020}, month={Jun}, pages={585–604} } @article{sahoo_guinness_reich_2019, title={A TEST FOR ISOTROPY ON A SPHERE USING SPHERICAL HARMONIC FUNCTIONS}, volume={29}, ISSN={["1996-8507"]}, DOI={10.5705/ss.202017.0475}, abstractNote={Analysis of geostatistical data is often based on the assumption that the spatial random field is isotropic. This assumption, if erroneous, can adversely affect model predictions and statistical inference. Nowadays many applications consider data over the entire globe and hence it is necessary to check the assumption of isotropy on a sphere. In this paper, a test for spatial isotropy on a sphere is proposed. The data are first projected onto the set of spherical harmonic functions. Under isotropy, the spherical harmonic coefficients are uncorrelated whereas they are correlated if the underlying fields are not isotropic. This motivates a test based on the sample correlation matrix of the spherical harmonic coefficients. In particular, we use the largest eigenvalue of the sample correlation matrix as the test statistic. Extensive simulations are conducted to assess the Type I errors of the test under different scenarios. We show how temporal correlation affects the test and provide a method for handling temporal correlation. We also gauge the power of the test as we move away from isotropy. The method is applied to the near-surface air temperature data which is part of the HadCM3 model output. Although we do not expect global temperature fields to be isotropic, we propose several anisotropic models with increasing complexity, each of which has an isotropic process as model component and we apply the test to the isotropic component in a sequence of such models as a method of determining how well the models capture the anisotropy in the fields.}, number={3}, journal={STATISTICA SINICA}, author={Sahoo, Indranil and Guinness, Joseph and Reich, Brian J.}, year={2019}, month={Jul}, pages={1253–1276} } @article{sharma_muyskens_guinness_polizzotto_fuentes_tappero_chen-wiegart_thieme_williams_acerbo_et al._2019, title={Multi-element effects on arsenate accumulation in a geochemical matrix determined using mu-XRF, mu-XANES and spatial statistics}, volume={26}, ISSN={["1600-5775"]}, DOI={10.1107/S1600577519012785}, abstractNote={Soils regulate the environmental impacts of trace elements, but direct measurements of reaction mechanisms in these complex, multi-component systems can be challenging. The objective of this work was to develop approaches for assessing effects of co-localized geochemical matrix elements on the accumulation and chemical speciation of arsenate applied to a soil matrix. Synchrotron X-ray fluorescence microprobe (µ-XRF) images collected across 100 µm × 100 µm and 10 µm × 10 µm regions of a naturally weathered soil sand-grain coating before and after treatment with As(V) solution showed strong positive partial correlations ( r ′ = 0.77 and 0.64, respectively) between accumulated As and soil Fe, with weaker partial correlations ( r ′ > 0.1) between As and Ca, and As and Zn in the larger image. Spatial and non-spatial regression models revealed a dominant contribution of Fe and minor contributions of Ca and Ti in predicting accumulated As, depending on the size of the sample area analyzed. Time-of-flight secondary ion mass spectrometry analysis of an area of the sand grain showed a significant correlation ( r = 0.51) between Fe and Al, so effects of Fe versus Al (hydr)oxides on accumulated As could not be separated. Fitting results from 25 As K -edge microscale X-ray absorption near-edge structure (µ-XANES) spectra collected across a separate 10 µm × 10 µm region showed ∼60% variation in proportions of Fe(III) and Al(III)-bound As(V) standards, and fits to µ-XANES spectra collected across the 100 µm × 100 µm region were more variable. Consistent with insights from studies on model systems, the results obtained here indicate a dominance of Fe and possibly Al (hydr)oxides in controlling As(V) accumulation within microsites of the soil matrix analyzed, but the analyses inferred minor augmentation from co-localized Ti, Ca and possibly Zn.}, journal={JOURNAL OF SYNCHROTRON RADIATION}, author={Sharma, Aakriti and Muyskens, Amanda and Guinness, Joseph and Polizzotto, Matthew L. and Fuentes, Montserrat and Tappero, Ryan V. and Chen-Wiegart, Yu-chen K. and Thieme, Juergen and Williams, Garth J. and Acerbo, Alvin S. and et al.}, year={2019}, month={Nov}, pages={1967–1979} } @article{guinness_hammerling_2018, title={Compression and Conditional Emulation of Climate Model Output}, volume={113}, ISSN={0162-1459 1537-274X}, url={http://dx.doi.org/10.1080/01621459.2017.1395339}, DOI={10.1080/01621459.2017.1395339}, abstractNote={Numerical climate model simulations run at high spatial and temporal resolutions generate massive quantities of data. As our computing capabilities continue to increase, storing all of the data is not sustainable, and thus it is important to develop methods for representing the full datasets by smaller compressed versions. We propose a statistical compression and decompression algorithm based on storing a set of summary statistics as well as a statistical model describing the conditional distribution of the full dataset given the summary statistics. The statistical model can be used to generate realizations representing the full dataset, along with characterizations of the uncertainties in the generated data. Thus, the methods are capable of both compression and conditional emulation of the climate models. Considerable attention is paid to accurately modeling the original dataset--one year of daily mean temperature data--particularly with regard to the inherent spatial nonstationarity in global fields, and to determining the statistics to be stored, so that the variation in the original data can be closely captured, while allowing for fast decompression and conditional emulation on modest computers.}, number={521}, journal={Journal of the American Statistical Association}, publisher={Informa UK Limited}, author={Guinness, Joseph and Hammerling, Dorit}, year={2018}, month={Jan}, pages={56–67} } @article{reich_guinness_vandekar_shinohara_staicu_2018, title={Fully Bayesian spectral methods for imaging data}, volume={74}, ISSN={["1541-0420"]}, DOI={10.1111/biom.12782}, abstractNote={Medical imaging data with thousands of spatially correlated data points are common in many fields. Methods that account for spatial correlation often require cumbersome matrix evaluations which are prohibitive for data of this size, and thus current work has either used low-rank approximations or analyzed data in blocks. We propose a method that accounts for nonstationarity, functional connectivity of distant regions of interest, and local signals, and can be applied to large multi-subject datasets using spectral methods combined with Markov Chain Monte Carlo sampling. We illustrate using simulated data that properly accounting for spatial dependence improves precision of estimates and yields valid statistical inference. We apply the new approach to study associations between cortical thickness and Alzheimer's disease, and find several regions of the cortex where patients with Alzheimer's disease are thinner on average than healthy controls.}, number={2}, journal={BIOMETRICS}, author={Reich, Brian J. and Guinness, Joseph and Vandekar, Simon N. and Shinohara, Russell T. and Staicu, Ana-Maria}, year={2018}, month={Jun}, pages={645–652} } @article{guinness_2018, title={Permutation and Grouping Methods for Sharpening Gaussian Process Approximations}, volume={60}, ISSN={["1537-2723"]}, DOI={10.1080/00401706.2018.1437476}, abstractNote={Vecchia's approximate likelihood for Gaussian process parameters depends on how the observations are ordered, which can be viewed as a deficiency because the exact likelihood is permutation-invariant. This article takes the alternative standpoint that the ordering of the observations can be tuned to sharpen the approximations. Advantageously chosen orderings can drastically improve the approximations, and in fact, completely random orderings often produce far more accurate approximations than default coordinate-based orderings do. In addition to the permutation results, automatic methods for grouping calculations of components of the approximation are introduced, having the result of simultaneously improving the quality of the approximation and reducing its computational burden. In common settings, reordering combined with grouping reduces Kullback-Leibler divergence from the target model by a factor of 80 and computation time by a factor of 2 compared to ungrouped approximations with default ordering. The claims are supported by theory and numerical results with comparisons to other approximations, including tapered covariances and stochastic partial differential equation approximations. Computational details are provided, including efficiently finding the orderings and ordered nearest neighbors, and profiling out linear mean parameters and using the approximations for prediction and conditional simulation. An application to space-time satellite data is presented.}, number={4}, journal={TECHNOMETRICS}, author={Guinness, Joseph}, year={2018}, pages={415–429} } @article{matli_fang_guinness_rabalais_craig_obenour_2018, title={Space-Time Geostatistical Assessment of Hypoxia in the Northern Gulf of Mexico}, volume={52}, ISSN={["1520-5851"]}, url={https://doi.org/10.1021/acs.est.8b03474}, DOI={10.1021/acs.est.8b03474}, abstractNote={Nearly every summer, a large hypoxic zone forms in the northern Gulf of Mexico. Research on the causes and consequences of hypoxia requires reliable estimates of hypoxic extent, which can vary at submonthly time scales due to hydro-meteorological variability. Here, we use an innovative space-time geostatistical model and data collected by multiple research organizations to estimate bottom-water dissolved oxygen (BWDO) concentrations and hypoxic area across summers from 1985 to 2016. We find that 27% of variability in BWDO is explained by deterministic trends with location, depth, and date, while correlated stochasticity accounts for 62% of observational variance within a range of 185 km and 28 days. Space-time modeling reduces uncertainty in estimated hypoxic area by 30% when compared to a spatial-only model, and results provide new insights into the temporal variability of hypoxia. For years with shelf-wide cruises in multiple months, hypoxia is most severe in July in 59% of years, 29% in August, and 12% in June. Also, midsummer cruise estimates of hypoxic area are only modestly correlated with summer-wide (June-August) average estimates ( r2 = 0.5), suggesting midsummer cruises are not necessarily reflective of seasonal hypoxic severity. Furthermore, summer-wide estimates are more strongly correlated with nutrient loading than midsummer estimates.}, number={21}, journal={ENVIRONMENTAL SCIENCE & TECHNOLOGY}, publisher={American Chemical Society (ACS)}, author={Matli, V. Rohith Reddy and Fang, Shiqi and Guinness, Joseph and Rabalais, Nancy. N. and Craig, J. Kevin and Obenour, Daniel R.}, year={2018}, month={Nov}, pages={12484–12493} } @article{castruccio_guinness_2017, title={An evolutionary spectrum approach to incorporate large-scale geographical descriptors on global processes}, volume={66}, ISSN={["1467-9876"]}, DOI={10.1111/rssc.12167}, abstractNote={We introduce a nonstationary spatio-temporal statistical model for gridded data on the sphere. The model specifies a computationally convenient covariance structure that depends on heterogeneous geography. Widely used statistical models on a spherical domain are nonstationary for different latitudes, but stationary at the same latitude (axial symmetry). This assumption has been acknowledged to be too restrictive for quantities such as surface temperature, whose statistical behavior is influenced by large scale geographical descriptors such as land and ocean. We propose an evolutionary spectrum approach that is able to account for different regimes across the Earth’s geography, and results in a more general and flexible class of models that vastly outperforms axially symmetric models and captures longitudinal patterns that would otherwise be assumed constant. The model can be estimated with in a multi-step conditional likelihood approximation that preserves the nonstationary features while allowing for easily distributed computations: we show how the fit of a data sets larger than 20 million data can be performed in less than one day on a state-of-the-art workstation. Once the parameters are estimated, it is possible to instantaneously generate surrogate runs from a common laptop. Further, the resulting estimates from the statistical model can be regarded as a synthetic description (i.e. a compression) of the space-time characteristics of an entire initial condition ensemble. Compared to traditional algorithms aiming at compressing the bit-by-bit information on each climate model run, the proposed approach achieves vastly superior compression rates.}, number={2}, journal={JOURNAL OF THE ROYAL STATISTICAL SOCIETY SERIES C-APPLIED STATISTICS}, author={Castruccio, Stefano and Guinness, Joseph}, year={2017}, month={Feb}, pages={329–344} } @article{guinness_fuentes_2017, title={Circulant Embedding of Approximate Covariances for Inference From Gaussian Data on Large Lattices}, volume={26}, ISSN={["1537-2715"]}, DOI={10.1080/10618600.2016.1164534}, abstractNote={Recently proposed computationally efficient Markov chain Monte Carlo (MCMC) and Monte Carlo expectation–maximization (EM) methods for estimating covariance parameters from lattice data rely on successive imputations of values on an embedding lattice that is at least two times larger in each dimension. These methods can be considered exact in some sense, but we demonstrate that using such a large number of imputed values leads to slowly converging Markov chains and EM algorithms. We propose instead the use of a discrete spectral approximation to allow for the implementation of these methods on smaller embedding lattices. While our methods are approximate, our examples indicate that the error introduced by this approximation is small compared to the Monte Carlo errors present in long Markov chains or many iterations of Monte Carlo EM algorithms. Our results are demonstrated in simulation studies, as well as in numerical studies that explore both increasing domain and fixed domain asymptotics. We compare the exact methods to our approximate methods on a large satellite dataset, and show that the approximate methods are also faster to compute, especially when the aliased spectral density is modeled directly. Supplementary materials for this article are available online.}, number={1}, journal={JOURNAL OF COMPUTATIONAL AND GRAPHICAL STATISTICS}, author={Guinness, Joseph and Fuentes, Montserrat}, year={2017}, pages={88–97} } @article{farjat_reich_guinness_whetten_mckeand_isik_2017, title={Optimal Seed Deployment Under Climate Change Using Spatial Models: Application to Loblolly Pine in the Southeastern US}, volume={112}, DOI={10.1080/01621459.2017.1292179}, abstractNote={Provenance tests are a common tool in forestry designed to identify superior genotypes for planting at specific locations. The trials are replicated experiments established with seed from parent trees collected from different regions and grown at several locations. In this work, a Bayesian spatial approach is developed for modeling the expected relative performance of seed sources using climate variables as predictors associated with the origin of seed source and the planting site. The proposed modeling technique accounts for the spatial dependence in the data and introduces a separable Matérn covariance structure that provides a flexible means to estimate effects associated with the origin and planting site locations. The statistical model was used to develop a quantitative tool for seed deployment aimed to identify the location of superior performing seed sources that could be suitable for a specific planting site under a given climate scenario. Cross-validation results indicate that the proposed spatial models provide superior predictive ability compared to multiple linear regression methods in unobserved locations. The general trend of performance predictions based on future climate scenarios suggests an optimal assisted migration of loblolly pine seed sources from southern and warmer regions to northern and colder areas in the southern USA. Supplementary materials for this article are available online.}, number={519}, journal={Journal of the American Statistical Association}, publisher={Informa UK Limited}, author={Farjat, Alfredo and Reich, Brian J. and Guinness, Joseph and Whetten, Ross and McKeand, Steven and Isik, Fikret}, year={2017}, month={Feb}, pages={909–920} } @article{guinness_fuentes_2016, title={Isotropic covariance functions on spheres: Some properties and modeling considerations}, volume={143}, ISSN={0047-259X}, url={http://dx.doi.org/10.1016/J.JMVA.2015.08.018}, DOI={10.1016/J.JMVA.2015.08.018}, abstractNote={Introducing flexible covariance functions is critical for interpolating spatial data since the properties of interpolated surfaces depend on the covariance function used for Kriging. An extensive literature is devoted to covariance functions on Euclidean spaces, where the Matérn covariance family is a valid and flexible parametric family capable of controlling the smoothness of corresponding stochastic processes. Many applications in environmental statistics involve data located on spheres, where less is known about properties of covariance functions, and where the Matérn is not generally a valid model with great circle distance metric. In this paper, we advance the understanding of covariance functions on spheres by defining the notion of and proving a characterization theorem for m times mean square differentiable processes on d-dimensional spheres. Stochastic processes on spheres are commonly constructed by restricting processes on Euclidean spaces to spheres of lower dimension. We prove that the resulting sphere-restricted process retains its differentiability properties, which has the important implication that the Matérn family retains its full range of smoothness when applied to spheres so long as Euclidean distance is used. The restriction operation has been questioned for using Euclidean instead of great circle distance. To address this question, we construct several new covariance functions and compare them to the Matérn with Euclidean distance on the task of interpolating smooth and non-smooth datasets. The Matérn with Euclidean distance is not outperformed by the new covariance functions or the existing covariance functions, so we recommend using the Matérn with Euclidean distance due to the ease with which it can be computed.}, journal={Journal of Multivariate Analysis}, publisher={Elsevier BV}, author={Guinness, Joseph and Fuentes, Montserrat}, year={2016}, month={Jan}, pages={143–152} } @article{guinness_fuentes_2015, title={Likelihood approximations for big nonstationary spatial temporal lattice data}, volume={25}, number={1}, journal={Statistica Sinica}, author={Guinness, J. and Fuentes, M.}, year={2015}, pages={329–349} } @article{guinness_fuentes_hesterberg_polizzotto_2014, title={Multivariate spatial modeling of conditional dependence in microscale soil elemental composition data}, volume={9}, ISSN={2211-6753}, url={http://dx.doi.org/10.1016/J.SPASTA.2014.03.009}, DOI={10.1016/J.SPASTA.2014.03.009}, abstractNote={The mobility and environmental impacts of toxic trace elements are regulated by their reactions with soils, which are complex heterogeneous mixtures of minerals and organic matter. We describe an experiment that maps the composition of elements on an individual soil sand grain using X-ray fluorescence microprobe analyses, after the grain is treated with arsenic solutions, resulting in multivariate spatial lattice maps of elemental abundance. To understand the behavior of arsenic in soils, it is important to disentangle the complex multivariate relationships among the elements in the sample. The abundance of most elements, including arsenic, correlates strongly with that of iron; but conditional on the amount of iron, some elements mitigate or potentiate the accumulation of arsenic. This problem motivates our work to define conditional correlation in spatial lattice models and give general conditions under which two components are conditionally uncorrelated given the rest. We describe how to enforce that two components are conditionally uncorrelated given a third in parametric models, which provides a basis for likelihood ratio tests for conditional correlation between arsenic and chromium given iron. We show how to apply our results to big datasets using the Whittle likelihood, and we demonstrate through simulation that tapering improves Whittle likelihood parameter estimates governing cross covariance.}, number={C}, journal={Spatial Statistics}, publisher={Elsevier BV}, author={Guinness, Joseph and Fuentes, Montserrat and Hesterberg, Dean and Polizzotto, Matthew}, year={2014}, month={Aug}, pages={93–108} } @article{guinness_stein_2013, title={INTERPOLATION OF NONSTATIONARY HIGH FREQUENCY SPATIAL-TEMPORAL TEMPERATURE DATA}, volume={7}, ISSN={["1932-6157"]}, DOI={10.1214/13-aoas633}, abstractNote={The Atmospheric Radiation Measurement program is a U.S. Department of Energy project that collects meteorological observations at several locations around the world in order to study how weather processes affect global climate change. As one of its initiatives, it operates a set of fixed but irregularly-spaced monitoring facilities in the Southern Great Plains region of the U.S. We describe methods for interpolating temperature records from these fixed facilities to locations at which no observations were made, which can be useful when values are required on a spatial grid. We interpolate by conditionally simulating from a fitted nonstationary Gaussian process model that accounts for the time-varying statistical characteristics of the temperatures, as well as the dependence on solar radiation. The model is fit by maximizing an approximate likelihood, and the conditional simulations result in well-calibrated confidence intervals for the predicted temperatures. We also describe methods for handling spatial-temporal jumps in the data to interpolate a slow-moving cold front.}, number={3}, journal={ANNALS OF APPLIED STATISTICS}, author={Guinness, Joseph and Stein, Michael L.}, year={2013}, month={Sep}, pages={1684–1708} }