Estimating correlations among demographic parameters in population models

Abstract Estimating correlations among demographic parameters is critical to understanding population dynamics and life‐history evolution, where correlations among parameters can inform our understanding of life‐history trade‐offs, result in effective applied conservation actions, and shed light on evolutionary ecology. The most common approaches rely on the multivariate normal distribution, and its conjugate inverse Wishart prior distribution. However, the inverse Wishart prior for the covariance matrix of multivariate normal distributions has a strong influence on posterior distributions. As an alternative to the inverse Wishart distribution, we individually parameterize the covariance matrix of a multivariate normal distribution to accurately estimate variances (σ 2) of, and process correlations (ρ) between, demographic parameters. We evaluate this approach using simulated capture–mark–recapture data. We then use this method to examine process correlations between adult and juvenile survival of black brent geese marked on the Yukon–Kuskokwim River Delta, Alaska (1988–2014). Our parameterization consistently outperformed the conjugate inverse Wishart prior for simulated data, where the means of posterior distributions estimated using an inverse Wishart prior were substantially different from the values used to simulate the data. Brent adult and juvenile annual apparent survival rates were strongly positively correlated (ρ = 0.563, 95% CRI 0.181–0.823), suggesting that habitat conditions have significant effects on both adult and juvenile survival. We provide robust simulation tools, and our methods can readily be expanded for use in other capture–recapture or capture‐recovery frameworks. Further, our work reveals limits on the utility of these approaches when study duration or sample sizes are small.

Estimating correlations among demographic parameters has often proven challenging, because sampling correlations can obscure process correlations (Anderson & Burnham, 1976;Link & Barker, 2005). For instance, the estimation of band recovery probability is confounded with the estimation of survival (Anderson & Burnham, 1976), affecting inference on the relationship between survival and harvest. Previous methods have employed approaches to achieve independent samples by comparing estimates of survival from the marked sample with estimates of harvest of the total population (Anderson & Burnham, 1976) or by partitioning the capture-mark-recovery data (Nichols & Hines, 1983). Recently, work on the effects of harvest on survival has focused on understanding process correlations (ρ) between survival and harvest rates (Arnold et al., 2016;Bartzen & Dufour, 2017;Sedinger, White, Espinosa, Partee, & Braun, 2010), where a strong negative correlation suggests additive relationships between survival and harvest, and minimal or no correlation may be indicative of compensation or partial compensation, although these relationships are complex (Arnold et al., 2016;Péron, 2013). Others have focused on the correlation between natural mortality and harvest mortality (Péron, 2013;Servanty et al., 2010), where no correlation may be indicative of additive harvest, and a negative correlation may be indicative of compensatory harvest. Ecologists have both highlighted (Arnold et al., 2016;Péron, 2013;Servanty et al., 2010) and debated (Arnold, Afton, Anteau, Koons, & Nicolai, 2017;Lindberg, Boomer, Schmutz, & Walker, 2017) the utility of these approaches. Capture-recapture (Riecke, Leach, Gibson, & Sedinger, 2018) and integrated population models can also be used to examine covariation among parameters such as survival and breeding propensity. This allows researchers to examine life-history trade-offs at the individual or population level , as well as correlations among other demographic parameters (Kindsvater et al., 2018;Koons, Arnold, & Schaub, 2017;Link & Barker, 2005;Schaub, Jakober, & Stauber, 2013). Conservation biologists can subsequently use these relationships to predict the potential of management actions to affect wildlife populations.
Bayesian hierarchical models allow for the separation of sampling and process correlations (Link & Barker, 2005). However, the choice of priors can dramatically impact posterior distributions and inference when using Bayesian models, where relatively uninformative priors are often favored (Link, Cam, Nichols, & Cooch, 2002) for objective Bayesian analyses. Further, conjugate priors are desirable computationally, because their full-conditional distributions are analytically tractable with a known distribution, and there is no need for Metropolis-Hastings updates or tuning. A conjugate prior for the covariance matrix (K × K) of a multivariate normal distribution is the inverse Wishart distribution, where quantitative ecologists often use this distribution with degrees of freedom equal to K + 1 and a K × K identity matrix as a scale matrix (Haff, 1979;Kéry & Schaub, 2012;Péron, 2013;Wishart, 1928). Thus, for a bi-variate normal distribution, where I is a 2 × 2 identity matrix. Critically, simulation work has revealed that the inverse Wishart prior is strongly informative for covariance matrices (Alvarez, Niemi, & Simpson, 2014). Specifically, the inverse Wishart prior can result in strong prior influence on the posterior distribution of variances and covariances, particularly when variances are small (σ 2 < 0.25), affecting the posterior distribution of correlations (Alvarez et al., 2014). The shape of the logistic link function, and life-history theory (Stearns, 1992), dictates that the standard deviations (σ), and consequently variances (σ 2 ), of demographic parameters of interest will typically be small (Gimenez et al., 2018). For instance, typical North American dabbling duck survival rates exhibit variances of approximately σ 2 = 0.1, where σ 2 is from Equation (1).
Given the structure of covariance matrices, where the covariance is equal to the product of the square roots of the applicable variances, and a process correlation between the applicable parameters (σ i σ j ρ i,j ), the use of the Wishart, scaled inverse Wishart (O'Malley & Zaslavsky, 2008), or hierarchical half-t (Huang & Wand, 2013) priors for a precision matrix often leads to the underestimation of ρ. This has important implications for the interpretation of process correlations, affecting biological inference and management decisions.
To alleviate current issues regarding the estimation of ρ and σ 2 when using the logistic link and the inverse Wishart distribution, we propose individual parameterization for the components of the To demonstrate the utility of our approach, we use simulated data, as well as a robust longitudinal dataset, capture-mark-recapture-resight data from female Pacific black brent geese (Branta bernicla nigricans; hereafter, brent) marked on the Yukon-Kuskokwim River Delta in Alaska (1988Alaska ( -2014. Brent breeding at the Tutakoke River Colony on the Yukon-Kuskokwim River Delta, Alaska, has been monitored from 1984 through the present Sedinger, Flint, & Lindberg, 1995). Over the duration of the Tutakoke River Colony brent demographic project, researchers have collected data on every avian demographic component, where an enhanced understanding of process correlations among demographic components may lead to more effective management actions, and improve our understanding of demographic buffering and life-history tradeoffs. We use multivariate normal distributions, and novel covariance matrix parameterizations for demographic parameters (Barnard, McCulloch, & Meng, 2000), to assess the relationships among adult and juvenile survival of brent to better understand life-history tradeoffs and the demographic buffering hypothesis (Gaillard, Festa-Bianchet, Yoccoz, Loison, & Toigo, 2000;Pfister, 1998;Stearns, 1992), which predicts that vital rates with greater demographic sensitivities will have lower variances. Further, we simulate thousands of biologically plausible capture-mark-recapture datasets to explore the effects of prior influence on inference in the presence of parameter uncertainty and imperfect observation when truth is known. Our work has important implications for researchers using multivariate normal distributions to estimate process correlations among demo- we used fixed variances of 0.1 for adult and juvenile survival rates, and randomly generated a process correlation between adult and juvenile survival for each simulation (v) using a uniform distribution to represent all possible relationships between these two demographic parameters, To generate capture-mark-recapture data, we first simulated each individual's latent state (z i,t ) given its intial release occasion (τ i ), initial age-class (juvenile or adult), number of release occasions (T), and time-varying survival probabilities (Φ), We then simulated capture histories for each individual (y i,t ) given its current latent state (z i,t ) and a detection probability (p = .5), We simulated 1,000 populations and capture-mark-recapture datasets for each pairwise comparison of 10, 20, and 30 occasions, and 100, 1,000, and 5,000 releases per occasion, for a total of 9,000 simulated capture-mark-recapture datasets.

| Analyzing the simulated data
For each simulated dataset, we built age-specific Bayesian CJS models (Cormack, 1964;Jolly, 1965;Seber, 1965), identical to the model used to simulate the data. We first modeled each individual's latent state (z i,t ) given its release occasion (τ i ) and age at release, given Equation (4). We then modeled each individual's capture history as a function of its latent state, given Equation (5). For computational efficiency, we structured the capture-mark-recapture data in m arrays for juveniles and adults (Data availability statement). To examine process correlations between adult (ϕ ad ) and juvenile (ϕ ad ) survival, we drew occasion specific adult and juvenile survival probabilities from a multivariate normal distribution, where = ( ad , juv ) and Φ = ( ad,t , juv,t ). We specified vague priors for the means of adult and juvenile survival, F I G U R E 2 Correlations among standard deviations (σ k ) and process correlations (ρ) drawn from an inverse Wishart prior with an identity scale matrix and degrees of freedom equal to K + 1 (left), and from a multivariate normal distribution where hyperpriors of the covariance matrix were individually specified, where we used a half-normal distribution with a variance of ten for the standard deviations, and a uniform prior with lower and upper bounds of −1 and 1 for the process correlation (right). Note the strong correlation between values of (σ k ) and values of (ρ) when using the inverse Wishart prior

| Brent capture-mark-recapture data analysis
We built models to estimate adult and juvenile survival, and the correlation between these parameters, in exactly the same way as for the simulated data with one exception. We modeled age-related and temporal heterogeneity in detection probability. Brent do not breed until after their second year, and they are less likely to be detected at the breeding colony until this time. However, we do resight second- year females at the colony and recapture these individuals during ringing efforts. Thus, we modeled detection probabilities of second- year and adult brent separately based on age-class (k), where we estimated detection probability for both age-classes (p k ) as random variables drawn from a normal distribution with a mean ( p k ) and vari-

| Computational details
We simulated data using R 3.5.1 (R Core Team, 2018), and all analyses were conducted in JAGS (Plummer, 2003) using the 'jagsUI' package (Kellner, 2016). We ran two MCMC chains of 25,000 iterations for each model for simulated data, where we discarded the first 15,000 iterations and retained every fifth saved iteration (Schaub & Fletcher, 2015). For the brent data, we ran two MCMC chains of 100,000 iterations, where we discarded the first 50,000 iterations and retained every fifth saved iteration.

| Capture-mark-recapture simulation
When we used the same model to simulate and analyze the capturemark-recapture data with an inverse Wishart prior on the covariance matrix of the multivariate normal distribution, the means of the posterior distributions of the process correlation (ρ) between adult and juvenile survival were different than the values used to simulate the data for all combinations of releases and occasions (Table 2). Further, 95% Bayesian credible intervals only covered truth for approximately 75% of simulations (Table 3). Critically, coverage declined as |ρ| increased (Figure 3). However, when we used hyperpriors on the (7) ad ∼ Normal (0,1), juv ∼ Normal (0,1).  (Table 2). Further, credible interval coverage was adequate for all combinations of releases and occasions (Table 3).

| Correlations among survival rates of female brent
When we used an inverse Wishart distribution as a prior for the co-

| Implications for black brent populations
The estimated variance of adult survival rates was significantly less than the variance of juvenile survival rates ( Figure 5). This is consistent with the demographic buffering hypothesis, where we would expect population growth rates of long-lived organisms to be most sensitive to adult survival (Gaillard et al., 2000;Pfister, 1998).

TA B L E 2
Mean absolute value of the difference between the simulated process correlation used to generate the data, and the mean of the posterior distributions of the estimated process correlations (ρ), as a function of the number of occasions (t), the number of individuals released per occasion (n), and the parameterization of the covariance matrix used for the capturerecapture models

TA B L E 3
The proportion of simulations (out of 1,000) for which the value of ρ used to simulate the capture-mark-recapture data was included in the 95% Bayesian credible interval of the estimate of ρ as a function of the number of occasions (t), the number of individuals released per occasion (n), and the parameterization of the covariance matrix used for the capture-recapture models Therefore, adult survival should be the most consistent demographic component and be buffered from environmental variation or reproductive allocation (Gaillard et al., 2000;Pfister, 1998). As predicted, and indicated by previous research (Leach et al., 2017), we observed positive correlations between adult and juvenile survival.
Brent family groups remain together from hatch through late spring (Sedinger, Nicolai, Lensink, Wentworth, & Conant, 2007), where we would expect environmental variation which impacts adult survival to affect juvenile survival as well, despite strong carry-over effects of environmental conditions during growth on first-year survival, and lifetime fitness, of juvenile brent (Riecke et al., 2018;Sedinger et al., 1995Sedinger et al., , 2007. Long-term declines in both adult (Leach et al., 2017;Riecke et al., 2018) and juvenile survival (Leach et al., 2017) have critical implications for brent populations, which are also experiencing declines in fecundity (Ward, Amundson, Stehn, & Dau, 2018) and population size .

F I G U R E 3
Means of the posterior distributions for process correlations (ρ) among adult (ϕ ad ) and juvenile (ϕ juv ) survival estimated from simulated capture-mark-recapture data using inverse Wishart (red) and hyperprior (blue) parameterizations for the covariance matrix of the multivariate normal distribution from which the parameters were estimated for different combinations of releases (n = 100, 1,000, or 5,000) and occasions (t = 10, 20, or 30). The dashed line represents 1:1 congruence between the simulated data and the estimate

| Implications for future research
We examine the behavior of multivariate normal distributions for demographic parameters in the presence of parameter uncertainty for a variety of biologically plausible scenarios. We suggest that investigators simulate data and conduct power analyses prior to drawing inference from their data. We provide simple R script to perform these analyses, where short-term (<10 years) datasets with small sample sizes (<100 releases per year) generally yield inaccurate estimates of process correlations when using capture-mark-recapture or capture-mark-recovery data. Our work has important implications for previous (Péron, 2013;Servanty et al., 2010) and future research, and our parameterizations can be expanded to covariance matrices of K × K dimensions (Alvarez et al., 2014). Thus, our preliminary analyses highlight the power and importance of robust longitudinal datasets and the importance of model parameterization.
Additionally, our simulations do not induce variation or heterogeneity in detection probability (p), where variation in these probabilities  can affect the absolute bias and precision of estimates of ρ (e.g., increased detection probabilities or recovery and reporting rates lead to increased precision and accuracy in the estimation of ρ).  . We would also like to thank two anonymous reviewers for providing positive contributions to the manuscript, as well as Martyn Plummer for developing the "dmnorm.vcov()" function in JAGS.

CO N FLI C T O F I NTE R E S T
None declared. VvALv fJQ WI qWnIfH0E.