Erin M. Schliep

Joint Species Distribution Modeling: Dimension Reduction Using Dirichlet Processes

Species distribution models are used to evaluate the variables that affect the distribution and abundance of species and to predict biodiversity. Historically, such models have been fitted to each species independently. While independent models can provide useful information regarding distribution and abundance, they ignore the fact that, after accounting for environmental covariates, residual interspecies dependence persists. With stacking of individual models, misleading behaviors, may arise. In particular, individual models often imply too many species per location. Recently developed joint species distribution models have application to presence–absence, continuous or discrete abundance, abundance with large numbers of zeros, and discrete, ordinal, and compositional data. Here, we deal with the challenge of joint modeling for a large number of species. To appreciate the challenge in the simplest way, with just presence/absence (binary) response and say, S species, we have an S-way contingency table with 2 cell probabilities. Even if S is as small as 100 this is an enormous table, infeasible to work with without some structure to reduce dimension. We develop a computationally feasible approach to accommodate a large number of species (say order 10) that allows us to: 1) assess the dependence structure across species; 2) identify clusters of species that have similar dependence patterns; and 3) jointly predict species distributions. To do so, we build hierarchical models capturing dependence between species at the first or “data” stage rather than at a second or “mean” stage. We employ the Dirichlet process for clustering in a novel way to reduce dimension in the joint covariance structure. This last step makes computation tractable. We use Forest Inventory Analysis (FIA) data in the eastern region of the United States to demonstrate our method. It consists of presence–absence measurements for 112 tree species, observed east of the Mississippi. As a proof of concept for our dimension reduction approach, we also include simulations using continuous and binary data. ∗Postdoctoral Associate, SAMSI/Duke University, Research Triangle Park, NC 27709, taylor-rodriguez@samsi.info †Postdoctoral Researcher, SAMSI/North Carolina State University, Research Triangle Park, NC 27709, kkaufeld@samsi.info ‡Assistant Professor, Department of Statistics, University of Missouri, Columbia, MO 65211, schliepe@missouri.edu §Professor, Nicholas School of the Environment, Department of Statistical Science, Duke University, Durham, NC 27708, jimclark@duke.edu ¶Professor, Department of Statistical Science, Duke University, Durham, NC 27708, alan@stat.duke.edu c

Modeling change in forest biomass across the eastern US

Predictions of above-ground biomass and the change in above-ground biomass require attachment of uncertainty due the range of reported predictions for forests. Because above-ground biomass is seldom measured, there have been no opportunities to obtain such uncertainty estimates. Standard methods involve applying an allometric equation to each individual tree on sample plots and summing the individual values. There is uncertainty in the allometry which leads to uncertainty in biomass at the tree level. Due to interdependence between competing trees, the uncertainty at the plot level that results from aggregating individual tree biomass in this way is expected to overestimate variability. That is, the variance at the plot level should be less than the sum of the individual variances. We offer a modeling strategy to learn about change in biomass at the plot level and model cumulative uncertainty to accommodate this dependence among neighboring trees. The plot-level variance is modeled using a parametric density-dependent asymptotic function. Plot-by-time covariate information is introduced to explain the change in biomass. These features are incorporated into a hierarchical model and inference is obtain within a Bayesian framework. We analyze data for the eastern United States from the Forest Inventory and Analysis (FIA) Program of the US Forest Service. This region contains roughly 25,000 FIA monitored plots from which there are measurements of approximately 1 million trees spanning more than 200 tree species. Due to the high species richness in the FIA data, we combine species into plant functional types. We present predictions of biomass and change in biomass for two plant functional types.

Alternating Gaussian process modulated renewal processes for modeling threshold exceedances and durations

It is often of interest to model the incidence and duration of threshold exceedance events for an environmental variable over a set of monitoring locations. Such data arrive over continuous time and can be considered as observations of a two-state process yielding, sequentially, a length of time in the below threshold state followed by a length of time in the above threshold state, then returning to the below threshold state, etc. We have a two-state continuous time Markov process, often referred to as an alternating renewal process. The process is observed over a truncated time window and, within this window, duration in each state is modeled using a distinct cumulative intensity specification. Initially, we model each intensity over the window using a parametric regression specification. We extend the regression specification adding temporal random effects to enrich the model using a realization of a log Gaussian process over time. With only one type of renewal, this specification is referred to as a Gaussian process modulated renewal process. Here, we introduce Gaussian process modulation to the intensity for each state. Model fitting is done within a Bayesian framework. We clarify that fitting with a customary log Gaussian process specification over a lengthy time window is computationally infeasible. The nearest neighbor Gaussian process, which supplies sparse covariance structure, is adopted to enable tractable computation. We propose methods for both generating data under our models and for conducting model comparison. The model is applied to hourly ozone data for four monitoring sites at different locations across the United States for the ozone season of 2014. For each site, we obtain estimated profiles of up-crossing and down-crossing intensity functions through time. In addition, we obtain inference regarding the number of exceedances, the distribution of the duration of exceedance events, and the proportion of time in the above and below threshold state for any time interval.

Data augmentation and parameter expansion for independent or spatially correlated ordinal data

Data augmentation and parameter expansion can lead to improved iterative sampling algorithms for Markov chain Monte Carlo (MCMC). Data augmentation allows for simpler and more feasible simulation from a posterior distribution. Parameter expansion accelerates convergence of iterative sampling algorithms by increasing the parameter space. Data augmentation and parameter-expanded data augmentation MCMC algorithms are proposed for fitting probit models for independent ordinal response data. The algorithms are extended for fitting probit linear mixed models for spatially correlated ordinal data. The effectiveness of data augmentation and parameter-expanded data augmentation is illustrated using the probit model and ordinal response data, however, the approach can be used broadly across model and data types.

Biomass prediction using a density-dependent diameter distribution model

Prediction of aboveground biomass, particularly at large spatial scales, is necessary for estimating global-scale carbon sequestration. Since biomass can be measured only by sacrificing trees, total biomass on plots is never observed. Rather, allometric equations are used to convert individual tree diameter to individual biomass, perhaps with noise. The values for all trees on a plot are then summed to obtain a derived total biomass for the plot. Then, with derived total biomasses for a collection of plots, regression models, using appropriate environmental covariates, are employed to attempt explanation and prediction. Not surprisingly, when out-of-sample validation is examined, such a model will predict total biomass well for holdout data because it is obtained using exactly the same derived approach. Apart from the somewhat circular nature of the regression approach, it also fails to employ the actual observed plot level response data. At each plot, we observe a random number of trees, each with an associated diameter, producing a sample of diameters. A model based on this random number of tree diameters provides understanding of how environmental regressors explain abundance of individuals, which in turn explains individual diameters. We incorporate density dependence because the distribution of tree diameters over a plot of fixed size depends upon the number of trees on the plot. After fitting this model, we can obtain predictive distributions for individuallevel biomass and plot-level total biomass. We show that predictive distributions for plot-level biomass obtained from a density-dependent model for diameters will be much different from predictive distributions using the regression approach. Moreover, they can be more informative for capturing uncertainty than those obtained from modeling derived plot-level biomass directly. We develop a density-dependent diameter distribution model and illustrate with data from the national Forest Inventory and Analysis (FIA) database. We also describe how to scale predictions to larger spatial regions. Our predictions agree (in magnitude) with available wisdom on mean and variation in biomass at the hectare scale.

Assessing the joint behaviour of species traits as filtered by environment

Understanding and predicting how species traits are shaped by prevailing environmental conditions is an important yet challenging task in ecology. Functional trait-based approaches can replace potentially idiosyncratic species-specific response models in learning about community behaviour across environmental gradients. Customarily, models for traits given environment consider only trait means to predict species and functional diversity, as intra-taxon variability in traits is often thought to be negligible. A growing body of literature indicates that intra-taxon trait variability is substantial and critical in structuring plant communities and assessing ecosystem function. We propose flexible joint trait distribution models given environment and across species that incorporate intra-taxon variability as well as inter-site/plot variability. Using a Bayesian framework, our joint trait distribution models allow for mixed continuous, binary and ordinal trait variables and incorporate dependence among traits enabling both joint and conditional trait prediction at unobserved sites. The models can be used to inform about the well-known fourth-corner problem, which attempts to interpret trait-by-environment matrices. We demonstrate the utility of our methodology through joint predictive trait distributions for individual species as well as joint community-weighted trait distributions for environments while incorporating intra-taxon trait variability. Explicit details on the probabilistic interpretations of the random trait-by-environment matrices obtained arising under our model are also provided to address the fourth-corner problem. Finally, our joint trait distribution model is applied to simulated and real vegetation data collected from the Greater Cape Floristic Region of South Africa. The proposed methodology places a fully model-based foundation on explaining intra-taxon trait variation given environment. It extends the utility and interpretability of commonly applied techniques for investigating community-weighted traits and illuminates randomness in the fourth-corner problem.