Supplemental materials - Combined statistical analyses

Repository:
https://github.com/egage/EVMP

Introduction

This document provides supporting information on statistical analyses of Elk Vegetation Management Plan (EVMP) data collected through the 2018 sampling season in the NPS NRR Report “Monitoring of Vegetation Response to Elk Population and Habitat Management in Rocky Mountain National Park - Analysis of Elk Vegetation Management Plan monitoring Data: 2008–2018.” Separate analyses were run for plots established in aspen, willow and upland communities and for the elk winter range (core and noncore areas) and the Kawuneeche Valley. For information on sampling and the broader analysis and interpretation of the results, refer to the body of the report, past analyses (Linda C. Zeigenfuss and Johnson 2015), and the original EVMP monitoring plan (Linda C. Zeigenfuss, Johnson, and Wiebe 2011).

Methods

The code and results presented here start from and build upon derived data produced in a separate code document focused on ingesting, compiling, and cleaning raw data provided by RMNP staff. In general, packages included in the “Tidyverse” ecosystem of R packages were used for data transformation and visualization (Wickham et al. 2019). Specialized packages particular to specific tasks (e.g., the “bayesplot” plot for visualization of Bayesian posterior distributions) were also used (Gabry et al. 2019).

Bayesian repeated measures analyses were fit separately for combined core and noncore winter range, core winter range, and noncore winter range plots with plots treated as random factors using the ‘stan_glmer’ function in the rstanrm package (Goodrich et al. 2020)(Brilleman et al. 2018). In addition, models were fit to particular subgroups of data (e.g., short saplings vs tall saplings).

Bayesian estimation was performed via MCMC adding independent weakly informed priors specific to the data type being modeled. Count data like aspen stem counts were modeled as poisson processes, while proportion data (e.g., cover) were modeled using a beta distribution (Ferrari and Cribari-Neto 2004), which constrains values from 0 to 1.

Bayesian models view a model parameter \(\theta\) as a [random variable]. In contrast, frequentist models treat model parameters as unknown constant. Rather than estimating an unknown constant, Bayesian modeling focuses on an unknown distribution of parameter values.

Using Bayes’ Law, model parameters are estimated:

\[ \begin{equation} \label{eq:bayeslaw} P(\theta | y, X) = \frac{P(y|\theta,X) \cdot P(\theta)}{P(y)} \propto P(y|\theta,X) \cdot P(\theta). \end{equation} \]

Where:

• \(X\): The observed data
• \(y\): The outcomes
• \(P(\theta|y,X)\): Our view of the possible value of our model parameters after seeing the data
• \(P(y|\theta,X)\): The likelihood of seeing the current outcome conditioned on the data in hand and a specific parameter value
• \(P(\theta)\): Our view of the possible values of our model parameters before observing any actual data, i.e. the prior
• \(P(y)\): The unconditional probability of the outcome, i.e. the probability after considering all possible parameter.

The objective is to estimate the posterior density of parameter \(P(\theta|y,X)\), typically described in terms of a point estimate of the expected value of posterior density (e.g., the median of the distribution) and a measure of the variance. To solve for \(P(\theta|y,X)\) analytically, the data likelihood \(P(y|\theta,X)\), prior \(P(\theta)\), and marginal likelihood of outcome \(P(y)\) are needed, but because there is usually no closed-form analytic solution, Markov chain Monte Carlo (MCMC) methods are used to numerically solve the posterior density by directly generating random draws of parameters.
General steps in a Bayesian analysis include:

Specify a joint distribution for the outcome(s) and unknowns, which is proportional to a posterior distribution of the unknowns conditional on the observed data.
Use Markov Chain Monte Carlo (MCMC) to draw from posterior distribution.
*Evaluation of model fit and revise the model as appropriate.

Aspen data

For an observation \(y\) that is assumed to follow a Poisson distribution (e.g., aspen stem count data), the likelihood for one observation can be written as:

\[\tfrac{1}{y!} \lambda^y e^{-\lambda},\]

where \(\lambda = E(y | \mathbf{x}) = g^{-1}(\eta)\) and \(\eta = \alpha + \mathbf{x}^\top \boldsymbol{\beta}\) is a linear predictor. For the Poisson distribution it is also true that \(\lambda = Var(y | \mathbf{x})\), i.e. the mean and variance are both \(\lambda\). The rate parameter \(\lambda\) must be positive, so with a Poisson GLM, the link function \(g\) maps between the positive real numbers \(\mathbb{R}^+\) (thesupport of \(\lambda\)) and the set of all real numbers \(\mathbb{R}\). When applied to a linear predictor \(\eta\) with values in \(\mathbb{R}\), the inverse link function \(g^{-1}(\eta)\) returns a positive real number.

The standard link function for a Poisson GLM is the log link \(g(x) = \ln{(x)}\). With the log link, the inverse link function is the exponential function and the likelihood for a single observation is:

\[\frac{g^{-1}(\eta)^y}{y!} e^{-g^{-1}(\eta)} = \frac{e^{\eta y}}{y!} e^{-e^\eta}.\]

With independent prior distributions, the joint posterior distribution for \(\alpha\) and \(\boldsymbol{\beta}\) in the Poisson model is proportional to the product of the priors and the \(N\) likelihood contributions:

\[f\left(\alpha,\boldsymbol{\beta} | \mathbf{y},\mathbf{X}\right) \propto f\left(\alpha\right) \times \prod_{k=1}^K f\left(\beta_k\right) \times \prod_{i=1}^N { \frac{g^{-1}(\eta_i)^{y_i}}{y_i!} e^{-g^{-1}(\eta_i)}}.\]

This is the posterior distribution drawn from when using MCMC.
The Probability of Direction (PD) was used to represent certainty associated with the most probable direction (positive or negative) of the effect (e.g., time class, fencing). The PD is correlated with the frequentist p-value, with a two-sided p-value of respectively .1, .05, and .01 approximated by a PD of 95%, 97.5%, and 99.5% (Makowski, Ben-Shachar, and Lüdecke 2019). The “region of practical equivalence” (ROPE) was used to evaluate the probability of a parameter being outside a range that can be considered as “practically no effect,” i.e., a region enclosing values that are equivalent to the null value for practical purposes (Kruschke 2018). The proportion of the 95% HDI credible interval falling within the ROPE was used as an index for an analog to frequentist “null-hypothesis” testing.

For analyses of aspen count data, entailing separate models for combined winter range, core and noncore winter range subsets of plots, a random effects model was used for inference:

_i = (+1x{1,i} +2x{2,i} + 3x{3,i} + 4x{4,i} + 5x{5,i} + _i)

_i  Normal (0 , ^2)

y_i  Poisson(_i)

where yi is the count of aspen stems in the ith plot; x1,i is a indicator variable that equals one if the observation in plot i was from 2013 and zero otherwise; x2,i is an indicator variable that equals one if the observation in plot i was from 2018 and zero otherwise; x3 is an indicator variable that equals one if the fenced and zero otherwise. Separate models were developed with inclusion of burning as a factor. Otherwise, analyses included both burned and unburned plots.

Interpretation of model coefficients: The intercept (exp(α)) represents the mean count of aspen in unfenced plots during the baseline year; exp(β1) is the multiplicative change in mean counts that occurred during the time interval between the baseline year and 2013; exp(β2) is the multiplicative change in mean counts that occurred during the time interval between the baseline year and 2018; and exp(β2) is the multiplicative change in mean counts that was caused by fencing.

Inferences were made on the difference between means.

So, for example, the mean of counts in year 2013 in the fenced plots could be computed as µˆfence, 2013 = exp(α + β1 + β3) and in unfenced plots as µˆopen, 2013 = exp(α + β1), and a contrast showing the difference between the means is then µˆfence, 2013 − µˆopen, 2013. Differences between year classes were computed similarly.

Willow plot and Upland plot data

Height

Bayesian repeated measures analyses of shrub height were fit separately for combined core and non-core winter range, core winter range, and non-core winter range plots. Separate models were also fit to different groups of species. Shrub height models were fit using weekly informative priors and a gamma distribution with the “stan_glmer” function in the “rstanarm” package (Goodrich et al. 2020)(Brilleman et al. 2018). Gamma regression is commonly used when the response variable is continuous and positive.

For the simplest case a GLM for a continuous outcome is simply a linear model and the likelihood for one observation is a conditionally normal probability density function.
\[\frac{1}{\sigma \sqrt{2 \pi}} e^{-\frac{1}{2} \left(\frac{y - \mu}{\sigma}\right)^2},\] where \(\mu = \alpha + \mathbf{x}^\top \boldsymbol{\beta}\) is a linear predictor and \(\sigma\) is the standard deviation of the error in predicting the outcome, \(y\). A linear predictor \(\eta = \alpha + \mathbf{x}^\top \boldsymbol{\beta}\) can more generally be related to the conditional mean of the outcome via a link function \(g\). This maps the range of values on which the outcome is defined. For the gamma distribution used in the modeling of willow height, the log link function was used. The likelihood is the product of the likelihood contributions of individual observations.

Bayesian analysis requires specifying prior distributions \(f(\alpha)\) and \(f(\boldsymbol{\beta})\) for the intercept and vector of regression coefficients. Prior distributions were represented by normal distributions with a mean zero and a small standard deviation (scale)). The joint posterior distribution for \(\alpha\) and \(\boldsymbol{\beta}\) is proportional to the product of the priors and the \(N\) likelihood contributions:

\[f\left(\boldsymbol{\beta} | \mathbf{y},\mathbf{X}\right) \propto f\left(\alpha\right) \times \prod_{k=1}^K f\left(\beta_k\right) \times \prod_{i=1}^N {f(y_i|\eta_i)},\]

where \(\mathbf{X}\) is the matrix of predictors and \(\eta\) the linear predictor, i.e., the posterior distribution that stan_glmer draws from using MCMC.

Cover

Willow and upland shrub cover data were modeled using the beta distribution with the likelihood: \[ f(y_i | a, b) = \frac{y_i^{(a-1)}(1-y_i)^{(b-1)}}{B(a,b)} \] where \(B(\cdot)\) is the beta function. The shape parameters for the distribution are \(a\) and \(b\) and enter into the model according to the following transformations, \[ a = \mu\cdot\phi \\ b = (1-\mu)\cdot\phi \]

If \(g_1(\cdot)\) is some link function, the specification of the shape parameters, \(\mu = g_1^{-1}(\mathbf{X}\boldsymbol{\beta})\), where \(\boldsymbol{X}\) is a \(N\times K\) dimensional matrix of predictors, and \(\boldsymbol{\beta}\), is a \(K\) dimensional vector of parameters associated with each predictor. \(\phi\) is a scalar parameter. With the shape parameter values included, the likelihood takes the form: \[ f(y_i | \mu, \phi) = \frac{y_i^{(\mu\phi-1)}(1-y_i)^{((1-\mu)\phi-1)}}{B(\mu\phi,(1-\mu)\phi)} \]

Bayesian analysis requires specifying prior distributions \(f(\boldsymbol{\beta})\) and \(f(\phi)\) for the vector of regression coefficients and \(\phi\). When modeling \(\phi\) with a linear predictor a full Bayesian analysis requires specifying the prior distributions \(f(\boldsymbol{\beta})\) and \(f(\boldsymbol{\gamma})\).

With a single set of explanatory variables, the posterior distribution of \(\boldsymbol{\beta}\) and \(\phi\) is proportional to the product of the likelihood contributions, the \(K\) priors on the \(\beta_k\) parameters, and \(\phi\), \[ f(\boldsymbol{\beta},\phi|\mathbf{y},\mathbf{X}) \propto \prod_{i=1}^N f(y_i | a, b) \times \prod_{k=1}^K f(\beta_k) \times f(\phi) \]

With two sets of explanatory variables, the posterior distribution of \(\boldsymbol{\beta}\) and \(\boldsymbol{\gamma}\) is proportional to the product of the likelihood contribution, the \(K\) priors on the \(\beta_k\) parameters, and the \(J\) priors on the \(\gamma_j\) parameters,

\[ f(\boldsymbol{\beta},\boldsymbol{\gamma}|\mathbf{y},\mathbf{X}) \propto \prod_{i=1}^N f(y_i | a, b) \times \prod_{k=1}^K f(\beta_k) \times \prod_{j=1}^J f(\gamma_j) \]

Results

Aspen

Comparison of Sapling Classes

Aspen Saplings
Time class	site_type	sap_size	mean	sd	min	q1	med	q3	max	mad	iqr	n
BL	AC	small saplings	25.4	30.3	0	7.0	17.5	30.0	139	16.3	21.5	40
2013	AC	small saplings	45.0	82.8	0	11.0	22.0	45.0	523	22.2	34.0	45
2018	AC	small saplings	25.0	29.3	0	5.0	11.5	34.0	111	12.6	27.2	46
BL	AK	small saplings	10.1	6.6	3	4.5	9.0	14.5	22	6.7	9.0	8
2013	AK	small saplings	12.2	9.8	3	4.0	11.5	16.0	32	8.9	12.0	8
2018	AK	small saplings	14.8	11.3	0	6.0	15.0	20.5	35	8.9	11.8	8
BL	ANC	small saplings	11.9	12.4	0	3.5	7.5	16.0	45	6.7	11.8	20
2013	ANC	small saplings	14.4	12.3	1	4.0	14.0	20.0	53	14.8	16.0	21
2018	ANC	small saplings	21.0	14.8	2	6.0	18.0	33.0	51	19.3	27.0	21
BL	AC	tall saplings	0.1	0.3	0	0.0	0.0	0.0	1	0.0	0.0	40
2013	AC	tall saplings	3.6	10.1	0	0.0	0.0	1.0	49	0.0	1.0	45
2018	AC	tall saplings	5.3	12.6	0	0.0	0.0	3.0	73	0.0	3.0	46
BL	AK	tall saplings	0.0	0.0	0	0.0	0.0	0.0	0	0.0	0.0	8
2013	AK	tall saplings	0.0	0.0	0	0.0	0.0	0.0	0	0.0	0.0	8
2018	AK	tall saplings	0.0	0.0	0	0.0	0.0	0.0	0	0.0	0.0	8
BL	ANC	tall saplings	0.8	2.9	0	0.0	0.0	0.0	13	0.0	0.0	20
2013	ANC	tall saplings	1.0	2.3	0	0.0	0.0	1.0	10	0.0	1.0	21
2018	ANC	tall saplings	1.9	2.8	0	0.0	0.0	3.0	10	0.0	3.0	21

Short Saplings (<=1.5m Tall)) - Combined Winter Range (AC+ANC)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.acanc.sm' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.27,5.25,6.16,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	Median	Mean	MAP	CI	CI_low	CI_high	pd	ps	Rhat	ESS
(Intercept)	1.42	1.42	1.43	0.9	1.17	1.66	1	1	1.01	752.18
time_class2013	0.51	0.51	0.50	0.9	0.44	0.59	1	1	1.00	6833.66
time_class2018	0.22	0.22	0.23	0.9	0.15	0.30	1	1	1.00	7402.64
fencedFenced	1.01	1.01	1.00	0.9	0.48	1.52	1	1	1.00	1121.80
time_class2013:fencedFenced	-1.04	-1.04	-1.05	0.9	-1.17	-0.89	1	1	1.00	8774.85
time_class2018:fencedFenced	-1.94	-1.94	-1.94	0.9	-2.13	-1.76	1	1	1.00	9878.32

Effect Description

Effect Existence and Significance Testing
Parameter	Median	CI	CI_low	CI_high	Direction	Significance	Large
(Intercept)	1.42	0.9	1.17	1.66	1	1	1.00
time_class2013	0.51	0.9	0.44	0.59	1	1	1.00
time_class2018	0.22	0.9	0.15	0.30	1	1	0.06
fencedFenced	1.01	0.9	0.48	1.52	1	1	0.99
time_class2013:fencedFenced	-1.04	0.9	-1.17	-0.89	1	1	1.00
time_class2018:fencedFenced	-1.94	0.9	-2.13	-1.76	1	1	1.00

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

## # A tibble: 6 x 5
##   contrast    fenced   estimate lower.HPD upper.HPD
##   <fct>       <fct>       <dbl>     <dbl>     <dbl>
## 1 BL - 2013   Unfenced   -0.507    -0.597    -0.420
## 2 BL - 2018   Unfenced   -0.224    -0.320    -0.137
## 3 2013 - 2018 Unfenced    0.283     0.216     0.352
## 4 BL - 2013   Fenced      0.532     0.393     0.673
## 5 BL - 2018   Fenced      1.72      1.52      1.92 
## 6 2013 - 2018 Fenced      1.19      0.998     1.40

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	fenced	estimate	lower.HPD	upper.HPD
BL - 2013	Unfenced	-0.51	-0.60	-0.42
BL - 2018	Unfenced	-0.22	-0.32	-0.14
2013 - 2018	Unfenced	0.28	0.22	0.35
BL - 2013	Fenced	0.53	0.39	0.67
BL - 2018	Fenced	1.72	1.52	1.92
2013 - 2018	Fenced	1.19	1.00	1.40

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
BL / 2013	Unfenced	0.60	0.55	0.66
BL / 2018	Unfenced	0.80	0.73	0.87
2013 / 2018	Unfenced	1.33	1.24	1.42
BL / 2013	Fenced	1.70	1.48	1.95
BL / 2018	Fenced	5.58	4.53	6.78
2013 / 2018	Fenced	3.28	2.65	3.97

EMMIP

Short Saplings (<=1.5m Tall) - Core Winter Range (AC)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.ac.sm' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.26,5.23,5.42,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	1.605	0.006	0.211	1.335	1.607	1.869	1207	1.002
time_class2013	0.569	0.001	0.054	0.499	0.569	0.638	6635	1.000
time_class2018	0.091	0.001	0.057	0.016	0.092	0.164	6958	1.000
fencedFenced	0.829	0.010	0.383	0.359	0.819	1.321	1460	1.002
time_class2013:fencedFenced	-1.102	0.001	0.089	-1.215	-1.101	-0.988	7011	1.000
time_class2018:fencedFenced	-1.814	0.001	0.118	-1.966	-1.814	-1.665	7745	1.000
mean_PPD	10.670	0.002	0.233	10.374	10.672	10.964	14547	1.000
log-posterior	-2887.734	0.156	6.998	-2897.009	-2887.421	-2879.003	2006	1.001

## [1] "(Intercept)" "2013"        "2018"        "Fenced"      "2013:Fenced"
## [6] "2018:Fenced"

Effect Description

Effect Existence and Significance Testing
Parameter	Median	CI	CI_low	CI_high	Direction	Significance	Large
(Intercept)	1.61	0.9	1.28	1.96	1.00	1.00	1.00
time_class2013	0.57	0.9	0.48	0.66	1.00	1.00	1.00
time_class2018	0.09	0.9	0.00	0.19	0.94	0.76	0.00
fencedFenced	0.82	0.9	0.22	1.46	0.99	0.98	0.92
time_class2013:fencedFenced	-1.10	0.9	-1.24	-0.95	1.00	1.00	1.00
time_class2018:fencedFenced	-1.81	0.9	-2.00	-1.62	1.00	1.00	1.00

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	fenced	estimate	lower.HPD	upper.HPD
BL - 2013	Unfenced	-0.57	-0.67	-0.46
BL - 2018	Unfenced	-0.09	-0.20	0.02
2013 - 2018	Unfenced	0.48	0.40	0.56
BL - 2013	Fenced	0.53	0.39	0.67
BL - 2018	Fenced	1.72	1.52	1.92
2013 - 2018	Fenced	1.19	0.99	1.39

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
BL / 2013	Unfenced	0.57	0.51	0.63
BL / 2018	Unfenced	0.91	0.82	1.02
2013 / 2018	Unfenced	1.61	1.49	1.75
BL / 2013	Fenced	1.70	1.48	1.95
BL / 2018	Fenced	5.59	4.54	6.82
2013 / 2018	Fenced	3.28	2.67	3.97

EMMIP

Short Saplings (<=1.5m Tall) - Noncore Winter Range (ANC)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.anc.sm' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.27,5.27])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	1.123	0.004	0.176	0.900	1.123	1.349	2128	1.001
time_class2013	0.204	0.001	0.089	0.090	0.204	0.319	10674	1.000
time_class2018	0.586	0.001	0.081	0.482	0.585	0.690	10369	1.000
mean_PPD	5.283	0.002	0.239	4.978	5.280	5.591	13733	1.000
log-posterior	-661.287	0.094	4.685	-667.383	-660.977	-655.549	2509	1.001

Effect Description

Effect Existence and Significance Testing
Parameter	Median	CI	CI_low	CI_high	Direction	Significance	Large
(Intercept)	1.12	0.9	0.84	1.41	1.00	1.00	1.00
time_class2013	0.20	0.9	0.06	0.35	0.99	0.96	0.14
time_class2018	0.59	0.9	0.45	0.72	1.00	1.00	1.00

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	estimate	lower.HPD	upper.HPD
BL - 2013	-0.20	-0.38	-0.03
BL - 2018	-0.59	-0.74	-0.43
2013 - 2018	-0.38	-0.52	-0.23

Pairwise contrasts
Results are given on the response scale.
contrast	ratio	lower.HPD	upper.HPD
BL / 2013	0.82	0.68	0.96
BL / 2018	0.56	0.47	0.65
2013 / 2018	0.68	0.58	0.78

EMMIP

Short Saplings (<=1.5m Tall) - Kawuneeche Valley (AK)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.ak.sm' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.27,5.27])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	1.066	0.004	0.251	0.758	1.070	1.371	3290	1.001
time_class2013	0.191	0.002	0.150	-0.003	0.191	0.383	9405	1.000
time_class2018	0.379	0.001	0.143	0.195	0.380	0.562	9091	1.000
mean_PPD	4.124	0.003	0.339	3.694	4.111	4.569	14475	1.000
log-posterior	-352.719	0.051	3.106	-356.824	-352.366	-349.045	3672	1.000

Panel plot of above MCMC area plots

Effect Description

Effect Existence and Significance Testing
Parameter	Median	CI	CI_low	CI_high	Direction	Significance	Large
(Intercept)	1.07	0.9	0.66	1.46	1.0	1.00	1.00
time_class2013	0.19	0.9	-0.05	0.45	0.9	0.82	0.23
time_class2018	0.38	0.9	0.14	0.61	1.0	0.99	0.71

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	estimate	lower.HPD	upper.HPD
BL - 2013	-0.19	-0.48	0.11
BL - 2018	-0.38	-0.66	-0.11
2013 - 2018	-0.19	-0.45	0.09

Pairwise contrasts
Results are given on the response scale.
contrast	ratio	lower.HPD	upper.HPD
BL / 2013	0.83	0.60	1.09
BL / 2018	0.68	0.51	0.89
2013 / 2018	0.83	0.62	1.06

EMMIP

Effects of Burning on Short Saplings (<1.5m Tall) - Combined Winter Range (AC+ANC)

Summary statistics
time_class	range_type	burned	fenced	variable	min	q1	med	q3	max
BL	core winter range	Unburned	Unfenced	stem_den_ac	0.0	0.0	161.9	1295.0	8741.2
BL	core winter range	Unburned	Fenced	stem_den_ac	0.0	0.0	809.4	1942.5	21529.3
BL	core winter range	Burned	Unfenced	stem_den_ac	0.0	0.0	161.9	161.9	161.9
BL	core winter range	Burned	Fenced	stem_den_ac	0.0	0.0	566.6	4046.9	4370.6
BL	non-core winter range	Unburned	Unfenced	stem_den_ac	0.0	0.0	323.7	809.4	6798.7
BL	non-core winter range	Burned	Unfenced	stem_den_ac	161.9	161.9	323.7	809.4	809.4
2013	core winter range	Unburned	Unfenced	stem_den_ac	0.0	0.0	485.6	2104.4	20234.3
2013	core winter range	Unburned	Fenced	stem_den_ac	0.0	161.9	485.6	1295.0	12302.5
2013	core winter range	Burned	Unfenced	stem_den_ac	0.0	1375.9	4532.5	18615.6	43058.6
2013	core winter range	Burned	Fenced	stem_den_ac	0.0	647.5	3075.6	3885.0	5341.9
2013	non-core winter range	Unburned	Unfenced	stem_den_ac	0.0	161.9	566.6	1133.1	5341.9
2013	non-core winter range	Burned	Unfenced	stem_den_ac	0.0	0.0	0.0	161.9	161.9
2018	core winter range	Unburned	Unfenced	stem_den_ac	0.0	0.0	485.6	1942.5	12140.6
2018	core winter range	Unburned	Fenced	stem_den_ac	0.0	0.0	161.9	647.5	2751.9
2018	core winter range	Burned	Unfenced	stem_den_ac	0.0	404.7	3642.2	6717.8	9712.5
2018	core winter range	Burned	Fenced	stem_den_ac	0.0	323.7	647.5	1133.1	1456.9
2018	non-core winter range	Unburned	Unfenced	stem_den_ac	0.0	323.7	1052.2	1699.7	4694.4
2018	non-core winter range	Burned	Unfenced	stem_den_ac	161.9	161.9	161.9	323.7	323.7

Effect Existence and Significance Testing
Parameter	Median	CI	CI_low	CI_high	Direction	Significance	Large
(Intercept)	1.32	0.9	1.06	1.57	1.00	1.00	1.00
fencedFenced	1.10	0.9	0.50	1.65	1.00	1.00	0.99
time_class2013	0.32	0.9	0.24	0.41	1.00	1.00	0.69
time_class2018	0.42	0.9	0.35	0.50	1.00	1.00	1.00
burnedBurned	1.07	0.9	0.05	2.01	0.96	0.96	0.90
fencedFenced:time_class2013	-1.16	0.9	-1.32	-1.00	1.00	1.00	1.00
fencedFenced:time_class2018	-2.37	0.9	-2.58	-2.15	1.00	1.00	1.00
fencedFenced:burnedBurned	-1.28	0.9	-2.78	0.18	0.92	0.91	0.86
time_class2013:burnedBurned	0.43	0.9	-0.17	1.07	0.88	0.85	0.63
time_class2018:burnedBurned	-0.75	0.9	-1.37	-0.13	0.97	0.96	0.87
fencedFenced:time_class2013:burnedBurned	0.96	0.9	0.24	1.61	0.98	0.98	0.94
fencedFenced:time_class2018:burnedBurned	1.89	0.9	1.13	2.62	1.00	1.00	1.00

Probability of Direction

Pairwise contrasts
Results are given on the response scale
contrast	burned	fenced	estimate	lower.HPD	upper.HPD
2013 - BL	Unburned	Unfenced	1.3834719	1.2578086	1.5236966
2018 - BL	Unburned	Unfenced	1.5249948	1.3921403	1.6772559
2018 - 2013	Unburned	Unfenced	1.1022873	1.0105488	1.1940177
2013 - BL	Burned	Unfenced	2.1237722	1.0644596	4.5978834
2018 - BL	Burned	Unfenced	0.7217646	0.3594195	1.5912529
2018 - 2013	Burned	Unfenced	0.3399078	0.2973678	0.3877800
2013 - BL	Unburned	Fenced	0.4354735	0.3674030	0.5121292
2018 - BL	Unburned	Fenced	0.1428292	0.1133788	0.1809495
2018 - 2013	Unburned	Fenced	0.3283338	0.2591522	0.4185422
2013 - BL	Burned	Fenced	1.7397530	1.2831883	2.4205241
2018 - BL	Burned	Fenced	0.4471539	0.2910151	0.6614768
2018 - 2013	Burned	Fenced	0.2578845	0.1804910	0.3727204

EMMIP

Tall Saplings (1.5-2.5m Tall) - Combined Winter Range (AC+ANC)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.acanc.lg' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.26,5.24,6.16,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	-3.799	0.008	0.558	-4.516	-3.770	-3.109	5184	1.001
time_class2013	-0.237	0.003	0.317	-0.637	-0.241	0.172	13217	1.000
time_class2018	1.582	0.002	0.260	1.251	1.577	1.918	13698	1.000
fencedFenced	0.225	0.019	1.116	-1.172	0.234	1.622	3587	1.001
time_class2013:fencedFenced	4.767	0.009	0.810	3.790	4.700	5.821	8138	1.000
time_class2018:fencedFenced	2.746	0.009	0.790	1.796	2.682	3.762	8169	1.000
mean_PPD	1.258	0.001	0.080	1.155	1.256	1.360	13938	1.000
log-posterior	-429.114	0.142	7.690	-439.142	-428.756	-419.505	2950	1.001

Effect Description

Effect Existence and Significance Testing
Parameter	Median	CI	CI_low	CI_high	Direction	Significance	Large
(Intercept)	-3.77	0.9	-4.70	-2.88	1.00	1.00	1.00
time_class2013	-0.24	0.9	-0.73	0.30	0.78	0.73	0.43
time_class2018	1.58	0.9	1.15	2.01	1.00	1.00	1.00
fencedFenced	0.23	0.9	-1.61	2.06	0.59	0.57	0.48
time_class2013:fencedFenced	4.70	0.9	3.41	6.01	1.00	1.00	1.00
time_class2018:fencedFenced	2.68	0.9	1.47	4.00	1.00	1.00	1.00

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	fenced	estimate	lower.HPD	upper.HPD
BL - 2013	Unfenced	0.24	-0.36	0.88
BL - 2018	Unfenced	-1.58	-2.08	-1.06
2013 - 2018	Unfenced	-1.81	-2.26	-1.40
BL - 2013	Fenced	-4.46	-6.05	-3.18
BL - 2018	Fenced	-4.25	-5.80	-2.93
2013 - 2018	Fenced	0.20	-0.04	0.44

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
BL / 2013	Unfenced	1.27	0.60	2.19
BL / 2018	Unfenced	0.21	0.11	0.32
2013 / 2018	Unfenced	0.16	0.10	0.24
BL / 2013	Fenced	0.01	0.00	0.03
BL / 2018	Fenced	0.01	0.00	0.04
2013 / 2018	Fenced	1.22	0.95	1.53

EMMIP

pl.emmip.sm <- plotF.emmip_FxTC_response(stmod_stally2.acanc.sm) +
  labs(title="Short saplings") +
  theme(legend.position = "none")

pl.emmip.lg <- plotF.emmip_FxTC_response(stmod_stally2.acanc.lg) +
  labs(title="Tall saplings") +
  theme(legend.position = c(.15, .85))

pl.emmip.sm + pl.emmip.lg +
  plot_annotation(title = "Combined winter range")

ggsave("./output/figures_exported/emmip_TCxF_acanc_lg_sm_sap.png", width = 7.25, height = 4.25) # save plot

Replacement Figure 13 (#1)

# revision2022
#### Replacement Figure 13 (#1):
## bw
pl.emmip.sm.bw <- plotF.emmip_FxTC_response_bw(stmod_stally2.acanc.sm) +
  # labs(title="Short saplings") +
  labs(y = "Stem count") +
  theme(legend.position = "none") +
  ylim(0,20)


pl.emmip.lg.bw <- plotF.emmip_FxTC_response_bw(stmod_stally2.acanc.lg) +
  # labs(title="Tall saplings") +
  # labs(y = "Stem count") +
  theme(legend.position = c(.2, .85))

pl.emmip.sm.bw + pl.emmip.lg.bw +
  # plot_annotation(title = "Combined winter range") +
  plot_annotation(tag_levels = "A")

ggsave("./output/figures_202202/emmip/Fig13_1_emmip_TCxF_acanc_lg_sm_sap_bw.png", width = 7.25, height = 3.75) # save plot

ggsave("./output/figures_202202/emmip/Fig13_1_emmip_TCxF_acanc_lg_sm_sap_bw.pdf", width = 7.25, height = 3.75) # save plot

Tall Saplings (1.5-2.5m Tall) - Core Winter Range (AC)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.ac.lg' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.25,5.23,5.42,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	-6.161	0.016	1.126	-7.629	-6.085	-4.769	4717	1.000
time_class2013	0.352	0.015	1.084	-0.972	0.279	1.762	5356	1.000
time_class2018	3.758	0.013	0.973	2.620	3.656	5.048	5250	1.000
fencedFenced	2.541	0.024	1.450	0.717	2.515	4.434	3628	1.000
time_class2013:fencedFenced	4.227	0.018	1.301	2.561	4.246	5.864	5133	1.000
time_class2018:fencedFenced	0.618	0.017	1.210	-0.960	0.655	2.111	5275	1.000
mean_PPD	1.557	0.001	0.110	1.416	1.553	1.698	14013	1.000
log-posterior	-291.934	0.118	6.258	-300.133	-291.666	-284.054	2813	1.002

## [1] "(Intercept)" "2013"        "2018"        "Fenced"      "2013:Fenced"
## [6] "2018:Fenced"

Effect Description

Effect Existence and Significance Testing
Parameter	Median	CI	CI_low	CI_high	Direction	Significance	Large
(Intercept)	-6.08	0.9	-7.98	-4.33	1.00	1.00	1.00
time_class2013	0.28	0.9	-1.43	2.07	0.61	0.59	0.49
time_class2018	3.66	0.9	2.30	5.38	1.00	1.00	1.00
fencedFenced	2.51	0.9	0.24	4.98	0.96	0.96	0.94
time_class2013:fencedFenced	4.25	0.9	2.02	6.31	1.00	1.00	1.00
time_class2018:fencedFenced	0.65	0.9	-1.31	2.66	0.71	0.70	0.62

Probability of Direction

##                    Parameter n   percent
##                  (Intercept) 1 0.1666667
##                 fencedFenced 1 0.1666667
##               time_class2013 1 0.1666667
##  time_class2013:fencedFenced 1 0.1666667
##               time_class2018 1 0.1666667
##  time_class2018:fencedFenced 1 0.1666667

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	fenced	estimate	lower.HPD	upper.HPD
BL - 2013	Unfenced	-0.28	-2.56	1.66
BL - 2018	Unfenced	-3.66	-5.64	-1.90
2013 - 2018	Unfenced	-3.37	-4.41	-2.40
BL - 2013	Fenced	-4.49	-6.14	-3.19
BL - 2018	Fenced	-4.29	-5.92	-2.97
2013 - 2018	Fenced	0.20	-0.03	0.43

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
BL / 2013	Unfenced	0.76	0.00	3.67
BL / 2018	Unfenced	0.03	0.00	0.09
2013 / 2018	Unfenced	0.03	0.01	0.07
BL / 2013	Fenced	0.01	0.00	0.03
BL / 2018	Fenced	0.01	0.00	0.04
2013 / 2018	Fenced	1.22	0.95	1.52

EMMIP

#### Replacement Figure 13 (#2)
# (Replacement Fig 13 revsion 2022)

pl.emmip.ac.sm <- plotF.emmip_FxTC_response_bw(stmod_stally2.ac.sm) +
  # labs(title="Short saplings") +
  theme(legend.position = "none") +
  scale_color_manual(values = c("black","grey 75")) + labs(x = "Year") +
  labs(y="Stem count")

pl.emmip.ac.lg <- plotF.emmip_FxTC_response_bw(stmod_stally2.ac.lg) +
  # labs(title="Tall saplings") +
  theme(legend.position = c(.20, .85)) +
  scale_color_manual(values = c("black","grey 75")) + labs(x = "Year") +
  labs(y="Stem count")

pl.emmip.ac.sm + pl.emmip.ac.lg +
  # plot_annotation(title = "Core winter range") +
  plot_annotation(tag_levels = "A")

ggsave("./output/figures_202202/emmip/Fig13_2_emmip_TCxF_ac_lg_sm_sap.png", width = 7.25, height = 3.75) # save plot

ggsave("./output/figures_202202/emmip/Fig13_2_emmip_TCxF_ac_lg_sm_sap.pdf", width = 7.25, height = 3.75) # save plot

pl.emmip.ac.sm.bw <- plotF.emmip_FxTC_response_bw(stmod_stally2.ac.sm) +
  labs(title="Short saplings") +
  theme(legend.position = "none")

pl.emmip.ac.lg.bw <- plotF.emmip_FxTC_response_bw(stmod_stally2.ac.lg) +
  labs(title="Tall saplings") +
  theme(legend.position = c(.15, .85))

pl.emmip.ac.sm.bw + pl.emmip.ac.lg.bw +
  plot_annotation(title = "Core winter range")

ggsave("./output/figures_exported/emmip_TCxF_ac_lg_sm_sap_bw.png", width = 7.25, height = 4.25) # save plot

Tall Saplings (1.5-2.5m Tall) - Noncore Winter Range (ANC)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.anc.lg' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.26,5.26])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	-2.601	0.011	0.706	-3.519	-2.546	-1.750	3920	1.000
time_class2013	0.206	0.003	0.330	-0.214	0.204	0.633	11797	1.000
time_class2018	0.859	0.003	0.290	0.496	0.851	1.234	12114	1.000
mean_PPD	0.632	0.001	0.101	0.508	0.629	0.766	14084	1.000
log-posterior	-126.719	0.074	4.347	-132.493	-126.444	-121.384	3456	1.001

Effect Description

Effect Existence and Significance Testing
Parameter	Median	CI	CI_low	CI_high	Direction	Significance	Large
(Intercept)	-2.55	0.9	-3.70	-1.45	1.00	1.00	1.00
time_class2013	0.20	0.9	-0.36	0.72	0.73	0.68	0.39
time_class2018	0.85	0.9	0.39	1.33	1.00	1.00	0.98

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	estimate	lower.HPD	upper.HPD
BL - 2013	-0.20	-0.85	0.44
BL - 2018	-0.85	-1.43	-0.29
2013 - 2018	-0.65	-1.20	-0.13

Pairwise contrasts
Results are given on the response scale.
contrast	ratio	lower.HPD	upper.HPD
BL / 2013	0.82	0.39	1.45
BL / 2018	0.43	0.21	0.69
2013 / 2018	0.52	0.28	0.84

## emmip
## fails b/c no groups to contrast
pl.emmip.anc.sm <- plotF.emmip_TC_response_bw(stmod_stally2.anc.sm) +
  labs(y = "Stem count", title="Short saplings") +
  theme(legend.position = "none") +
  ylim(0,8)

pl.emmip.anc.lg <- plotF.emmip_TC_response_bw(stmod_stally2.anc.lg) +
  labs(y="", title="Tall saplings") +
  theme(legend.position = c(.15, .85)) +
  ylim(0,8)

pl.emmip.anc.sm + pl.emmip.anc.lg +
  plot_annotation(title = "Non-core winter range")

Replacement Figure 13 (#3)

Revised figure 13 - 3

Tall Saplings (1.5-2.5m Tall) - Kawuneeche Valley (AK)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.ak.lg' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.25,5.25])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	-5.523	0.021	2.050	-8.309	-5.258	-3.106	9246	1
time_class2013	-0.115	0.029	2.921	-3.853	-0.055	3.502	10382	1
time_class2018	-0.148	0.028	2.939	-3.930	-0.051	3.488	10640	1
mean_PPD	0.018	0.000	0.026	0.000	0.000	0.042	15306	1
log-posterior	-21.672	0.033	2.596	-25.135	-21.351	-18.631	6129	1

Panel plot of above MCMC area plots

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	estimate	lower.HPD	upper.HPD
BL - 2013	0.06	-5.73	6.01
BL - 2018	0.05	-5.67	5.97
2013 - 2018	0.02	-6.64	6.95

Pairwise contrasts
Results are given on the response scale.
contrast	ratio	lower.HPD	upper.HPD
BL / 2013	1.056940	1.096282e-05	158.2514
BL / 2018	1.052294	2.363998e-06	162.9399
2013 / 2018	1.015152	1.966738e-07	296.5660

Fig 21 revised 2022

Effects of Burning on Tall Saplings (1.5-2.5m Tall) - Combined Winter Range (AC+ANC)

Summary statistics
time_class	range_type	burned	fenced	variable	min	q1	med	q3	max
BL	core winter range	Unburned	Unfenced	percent_tot	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
BL	core winter range	Unburned	Unfenced	stem_den_ac	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
BL	core winter range	Unburned	Unfenced	stem_den_ha	1.000000e-07	0.0000001	0.0000001	0.0000001	1.000000e-07
BL	core winter range	Unburned	Unfenced	stem_tally	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
BL	core winter range	Unburned	Unfenced	sum_stem	0.000000e+00	6.0000000	14.0000000	26.0000000	5.900000e+01
BL	core winter range	Unburned	Unfenced	year	2.006000e+03	2006.0000000	2007.0000000	2009.0000000	2.009000e+03
BL	core winter range	Unburned	Fenced	percent_tot	0.000000e+00	0.0000000	0.0000000	0.0000000	5.000000e-02
BL	core winter range	Unburned	Fenced	stem_den_ac	0.000000e+00	0.0000000	0.0000000	0.0000000	1.618744e+02
BL	core winter range	Unburned	Fenced	stem_den_ha	1.000000e-07	0.0000001	0.0000001	0.0000001	4.000000e+02
BL	core winter range	Unburned	Fenced	stem_tally	0.000000e+00	0.0000000	0.0000000	0.0000000	1.000000e+00
BL	core winter range	Unburned	Fenced	sum_stem	0.000000e+00	13.0000000	24.5000000	46.0000000	1.390000e+02
BL	core winter range	Unburned	Fenced	year	2.007000e+03	2007.0000000	2008.0000000	2009.0000000	2.009000e+03
BL	core winter range	Burned	Unfenced	percent_tot	0.000000e+00	0.0000000	0.1666667	0.3333333	3.333333e-01
BL	core winter range	Burned	Unfenced	stem_den_ac	0.000000e+00	0.0000000	80.9372000	161.8744000	1.618744e+02
BL	core winter range	Burned	Unfenced	stem_den_ha	1.000000e-07	0.0000001	200.0000001	400.0000001	4.000000e+02
BL	core winter range	Burned	Unfenced	stem_tally	0.000000e+00	0.0000000	0.5000000	1.0000000	1.000000e+00
BL	core winter range	Burned	Unfenced	sum_stem	3.000000e+00	3.0000000	3.0000000	3.0000000	3.000000e+00
BL	core winter range	Burned	Unfenced	year	2.009000e+03	2009.0000000	2009.0000000	2009.0000000	2.009000e+03
BL	core winter range	Burned	Fenced	percent_tot	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
BL	core winter range	Burned	Fenced	stem_den_ac	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
BL	core winter range	Burned	Fenced	stem_den_ha	1.000000e-07	0.0000001	0.0000001	0.0000001	1.000000e-07
BL	core winter range	Burned	Fenced	stem_tally	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
BL	core winter range	Burned	Fenced	sum_stem	2.700000e+01	27.0000000	29.5000000	32.0000000	3.200000e+01
BL	core winter range	Burned	Fenced	year	2.006000e+03	2006.0000000	2006.0000000	2006.0000000	2.006000e+03
2013	core winter range	Unburned	Unfenced	percent_tot	0.000000e+00	0.0000000	0.0000000	0.0000000	5.000000e-01
2013	core winter range	Unburned	Unfenced	stem_den_ac	0.000000e+00	0.0000000	0.0000000	0.0000000	1.618744e+02
2013	core winter range	Unburned	Unfenced	stem_den_ha	1.000000e-07	0.0000001	0.0000001	0.0000001	4.000000e+02
2013	core winter range	Unburned	Unfenced	stem_tally	0.000000e+00	0.0000000	0.0000000	0.0000000	1.000000e+00
2013	core winter range	Unburned	Unfenced	sum_stem	0.000000e+00	7.0000000	22.0000000	34.0000000	1.440000e+02
2013	core winter range	Unburned	Unfenced	year	2.013000e+03	2013.0000000	2013.0000000	2013.0000000	2.013000e+03
2013	core winter range	Unburned	Fenced	percent_tot	0.000000e+00	0.0000000	0.1846591	0.3142857	8.000000e-01
2013	core winter range	Unburned	Fenced	stem_den_ac	0.000000e+00	0.0000000	728.4348000	1456.8696000	4.694358e+03
2013	core winter range	Unburned	Fenced	stem_den_ha	1.000000e-07	0.0000001	1800.0000001	3600.0000001	1.160000e+04
2013	core winter range	Unburned	Fenced	stem_tally	0.000000e+00	0.0000000	4.5000000	9.0000000	2.900000e+01
2013	core winter range	Unburned	Fenced	sum_stem	4.000000e+00	10.0000000	28.0000000	70.0000000	1.040000e+02
2013	core winter range	Unburned	Fenced	year	2.013000e+03	2013.0000000	2013.0000000	2013.0000000	2.013000e+03
2013	core winter range	Burned	Unfenced	percent_tot	0.000000e+00	0.0000000	0.0000000	0.0000000	1.908397e-03
2013	core winter range	Burned	Unfenced	stem_den_ac	0.000000e+00	0.0000000	0.0000000	0.0000000	1.618744e+02
2013	core winter range	Burned	Unfenced	stem_den_ha	1.000000e-07	0.0000001	0.0000001	0.0000001	4.000000e+02
2013	core winter range	Burned	Unfenced	stem_tally	0.000000e+00	0.0000000	0.0000000	0.0000000	1.000000e+00
2013	core winter range	Burned	Unfenced	sum_stem	1.700000e+01	92.5000000	168.5000000	346.5000000	5.240000e+02
2013	core winter range	Burned	Unfenced	year	2.013000e+03	2013.0000000	2013.0000000	2013.0000000	2.013000e+03
2013	core winter range	Burned	Fenced	percent_tot	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
2013	core winter range	Burned	Fenced	stem_den_ac	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
2013	core winter range	Burned	Fenced	stem_den_ha	1.000000e-07	0.0000001	0.0000001	0.0000001	1.000000e-07
2013	core winter range	Burned	Fenced	stem_tally	0.000000e+00	0.0000000	0.0000000	0.0000000	0.000000e+00
2013	core winter range	Burned	Fenced	sum_stem	4.400000e+01	44.0000000	50.0000000	52.0000000	5.200000e+01
2013	core winter range	Burned	Fenced	year	2.013000e+03	2013.0000000	2013.0000000	2013.0000000	2.013000e+03
2018	core winter range	Unburned	Unfenced	percent_tot	0.000000e+00	0.0000000	0.0000000	0.0000000	2.295082e-01
2018	core winter range	Unburned	Unfenced	stem_den_ac	0.000000e+00	0.0000000	0.0000000	0.0000000	2.266242e+03
2018	core winter range	Unburned	Unfenced	stem_den_ha	1.000000e-07	0.0000001	0.0000001	0.0000001	5.600000e+03
2018	core winter range	Unburned	Unfenced	stem_tally	0.000000e+00	0.0000000	0.0000000	0.0000000	1.400000e+01
2018	core winter range	Unburned	Unfenced	sum_stem	0.000000e+00	6.5000000	18.0000000	46.5000000	9.500000e+01
2018	core winter range	Unburned	Unfenced	year	2.018000e+03	2018.0000000	2018.0000000	2018.0000000	2.018000e+03
2018	core winter range	Unburned	Fenced	percent_tot	0.000000e+00	0.0000000	0.2222222	0.3333333	6.111111e-01
2018	core winter range	Unburned	Fenced	stem_den_ac	0.000000e+00	0.0000000	323.7488000	1133.1208000	3.723111e+03
2018	core winter range	Unburned	Fenced	stem_den_ha	1.000000e-07	0.0000001	800.0000001	2800.0000001	9.200000e+03
2018	core winter range	Unburned	Fenced	stem_tally	0.000000e+00	0.0000000	2.0000000	7.0000000	2.300000e+01
2018	core winter range	Unburned	Fenced	sum_stem	3.000000e+00	4.0000000	9.0000000	21.0000000	5.700000e+01
2018	core winter range	Unburned	Fenced	year	2.018000e+03	2018.0000000	2018.0000000	2018.0000000	2.018000e+03
2018	core winter range	Burned	Unfenced	percent_tot	0.000000e+00	0.0000000	0.0000000	0.1114394	2.303371e-01
2018	core winter range	Burned	Unfenced	stem_den_ac	0.000000e+00	0.0000000	0.0000000	2994.6764000	6.636850e+03
2018	core winter range	Burned	Unfenced	stem_den_ha	1.000000e-07	0.0000001	0.0000001	7400.0000001	1.640000e+04
2018	core winter range	Burned	Unfenced	stem_tally	0.000000e+00	0.0000000	0.0000000	18.5000000	4.100000e+01
2018	core winter range	Burned	Unfenced	sum_stem	2.000000e+00	39.5000000	96.5000000	147.0000000	1.780000e+02
2018	core winter range	Burned	Unfenced	year	2.018000e+03	2018.0000000	2018.0000000	2018.0000000	2.018000e+03
2018	core winter range	Burned	Fenced	percent_tot	7.142857e-02	0.1428571	0.2169231	0.3333333	4.400000e-01
2018	core winter range	Burned	Fenced	stem_den_ac	1.618744e+02	323.7488000	1052.1836000	1780.6184000	2.104367e+03
2018	core winter range	Burned	Fenced	stem_den_ha	4.000000e+02	800.0000001	2600.0000001	4400.0000001	5.200000e+03
2018	core winter range	Burned	Fenced	stem_tally	1.000000e+00	2.0000000	6.5000000	11.0000000	1.300000e+01
2018	core winter range	Burned	Fenced	sum_stem	1.400000e+01	14.0000000	25.0000000	39.0000000	3.900000e+01
2018	core winter range	Burned	Fenced	year	2.018000e+03	2018.0000000	2018.0000000	2018.0000000	2.018000e+03

Pairwise contrasts
Results are given on the response scale
contrast	burned	fenced	estimate	lower.HPD	upper.HPD
2013 - BL	Unburned	Unfenced	1.347698e+00	7.257894e-01	2.507768e+00
2018 - BL	Unburned	Unfenced	4.217323e+00	2.562921e+00	7.303658e+00
2018 - 2013	Unburned	Unfenced	3.123271e+00	1.991909e+00	5.006520e+00
2013 - BL	Burned	Unfenced	1.314597e-01	7.629043e-03	3.272656e+00
2018 - BL	Burned	Unfenced	6.497319e+00	6.208916e-01	1.064042e+02
2018 - 2013	Burned	Unfenced	4.828052e+01	1.317339e+01	2.702049e+02
2013 - BL	Unburned	Fenced	9.300851e+01	2.492602e+01	4.712025e+02
2018 - BL	Unburned	Fenced	5.217879e+01	1.440035e+01	2.750313e+02
2018 - 2013	Unburned	Fenced	5.622801e-01	4.324993e-01	7.278152e-01
2013 - BL	Burned	Fenced	1.137793e-03	1.850106e-16	2.614599e+06
2018 - BL	Burned	Fenced	1.599496e+05	2.306616e+01	1.511368e+11
2018 - 2013	Burned	Fenced	2.190260e+08	2.755728e+01	2.063195e+20

EMMIP

Replacement Figure 17

# small sap
plot.emmip.smsap.acanc.burned.bw <- emmip(stmod_smsap_TCxFxB_acanc, fenced ~ time_class | burned, type = "response",
  breaks = seq(0.50, 1, by = 0.10),
  CIarg = list(lwd = 1.85, alpha = 1), dotarg = list(size = 3.5, alpha = 1)) +
  theme_minimal() +
  scale_color_manual(values = c("black","grey 75")) +
  labs(color="", x = "Year", y = "Stem count", subtitle = "Short saplings")  +
  # scale_color_manual(values = c("black","grey 75")) + 
  labs(x = "Year") +
  theme(legend.position = "none")

plot.emmip.smsap.acanc.burned.bw

## tall sap
plot.emmip.lgsap.acanc.burned.bw <- emmip(stmod_lgsap_TCxFxB_acanc, fenced ~ time_class | burned, type = "response",
  breaks = seq(0.50, 1, by = 0.10),
  CIarg = list(lwd = 1.85, alpha = 1), dotarg = list(size = 3.5, alpha = 1)) +
  theme_minimal() +
  scale_color_manual(values = c("black","grey 75")) +
  labs(color="", x = "Year", y = "Stem count", subtitle = "Tall saplings")  +
  scale_color_manual(values = c("black","grey 75")) + labs(x = "Year") +
  theme(legend.position = "bottom")

plot.emmip.lgsap.acanc.burned.bw

All Stem Diameters - Combined Winter Range (AC+ANC)

## [[1]]

## 
## [[2]]

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2.acanc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.27,5.25,6.16,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	0.932	0.005	0.145	0.744	0.934	1.118	825	1.003
time_class2013	0.489	0.001	0.045	0.432	0.489	0.547	7615	1.000
time_class2018	0.294	0.001	0.046	0.235	0.294	0.353	7248	1.000
fencedFenced	0.918	0.009	0.313	0.515	0.919	1.317	1174	1.002
time_class2013:fencedFenced	-0.582	0.001	0.078	-0.682	-0.583	-0.482	8729	1.000
time_class2018:fencedFenced	-1.194	0.001	0.091	-1.312	-1.193	-1.076	9005	1.000
mean_PPD	5.868	0.001	0.110	5.728	5.866	6.010	14023	1.000
log-posterior	-6003.994	0.210	8.598	-6015.143	-6003.626	-5993.249	1683	1.001

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	fenced	estimate	lower.HPD	upper.HPD
2013 - BL	Unfenced	0.49	0.40	0.58
2018 - BL	Unfenced	0.29	0.20	0.38
2018 - 2013	Unfenced	-0.20	-0.27	-0.13
2013 - BL	Fenced	-0.09	-0.22	0.03
2018 - BL	Fenced	-0.90	-1.06	-0.75
2018 - 2013	Fenced	-0.81	-0.95	-0.66

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
2013 / BL	Unfenced	1.63	1.49	1.78
2018 / BL	Unfenced	1.34	1.22	1.46
2018 / 2013	Unfenced	0.82	0.77	0.88
2013 / BL	Fenced	0.91	0.81	1.03
2018 / BL	Fenced	0.41	0.35	0.47
2018 / 2013	Fenced	0.45	0.38	0.51

All Stem Diameters - Core Winter Range (AC)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_stally2' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.26,5.23,5.42,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	1.108	0.007	0.210	0.842	1.108	1.375	1010	1.001
time_class2013	0.550	0.001	0.053	0.483	0.551	0.618	5452	1.000
time_class2018	0.177	0.001	0.055	0.107	0.178	0.247	5554	1.001
fencedFenced	0.730	0.011	0.373	0.264	0.728	1.208	1186	1.001
time_class2013:fencedFenced	-0.643	0.001	0.081	-0.747	-0.643	-0.538	5551	1.001
time_class2018:fencedFenced	-1.077	0.001	0.096	-1.199	-1.076	-0.956	6104	1.000
mean_PPD	7.027	0.001	0.148	6.838	7.026	7.218	13908	1.000
log-posterior	-4871.545	0.161	7.052	-4880.799	-4871.206	-4862.852	1929	1.000

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	fenced	estimate	lower.HPD	upper.HPD
2013 - BL	Unfenced	0.55	0.44	0.65
2018 - BL	Unfenced	0.18	0.07	0.29
2018 - 2013	Unfenced	-0.37	-0.45	-0.30
2013 - BL	Fenced	-0.09	-0.22	0.03
2018 - BL	Fenced	-0.90	-1.05	-0.74
2018 - 2013	Fenced	-0.81	-0.96	-0.66

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
2013 / BL	Unfenced	1.73	1.56	1.92
2018 / BL	Unfenced	1.19	1.07	1.33
2018 / 2013	Unfenced	0.69	0.64	0.74
2013 / BL	Fenced	0.91	0.80	1.03
2018 / BL	Fenced	0.41	0.35	0.47
2018 / 2013	Fenced	0.45	0.38	0.52

Effects of Burning on All Stems - Core Winter Range (AC)

Summary statistics
time_class	range_type	burned	fenced	variable	min	q1	med	q3	max
BL	core winter range	Unburned	Unfenced	percent_tot	0	0.000000e+00	0.000000e+00	3.333333e-01	1.000000e+00
BL	core winter range	Unburned	Unfenced	site_number	2	9.000000e+00	2.000000e+01	3.000000e+01	4.000000e+01
BL	core winter range	Unburned	Unfenced	stem_den_ac	0	0.000000e+00	0.000000e+00	3.237488e+02	8.741218e+03
BL	core winter range	Unburned	Unfenced	stem_den_ha	0	0.000000e+00	0.000000e+00	8.000000e+02	2.160000e+04
BL	core winter range	Unburned	Unfenced	stem_tally	0	0.000000e+00	0.000000e+00	2.000000e+00	5.400000e+01
BL	core winter range	Unburned	Unfenced	sum_stem	0	6.000000e+00	1.400000e+01	2.600000e+01	5.900000e+01
BL	core winter range	Unburned	Unfenced	utm_e_nad83	445348	4.464030e+05	4.465360e+05	4.473970e+05	4.505480e+05
BL	core winter range	Unburned	Unfenced	utm_n_nad83	4467366	4.468377e+06	4.472409e+06	4.473033e+06	4.473650e+06
BL	core winter range	Unburned	Unfenced	year	2006	2.006000e+03	2.007000e+03	2.009000e+03	2.009000e+03
BL	core winter range	Unburned	Fenced	percent_tot	0	0.000000e+00	3.597122e-02	3.333333e-01	1.000000e+00
BL	core winter range	Unburned	Fenced	site_number	8	3.100000e+01	3.950000e+01	6.000000e+01	6.200000e+01
BL	core winter range	Unburned	Fenced	stem_den_ac	0	0.000000e+00	0.000000e+00	9.712464e+02	2.152930e+04
BL	core winter range	Unburned	Fenced	stem_den_ha	0	0.000000e+00	0.000000e+00	2.400000e+03	5.320000e+04
BL	core winter range	Unburned	Fenced	stem_tally	0	0.000000e+00	0.000000e+00	6.000000e+00	1.330000e+02
BL	core winter range	Unburned	Fenced	sum_stem	0	1.300000e+01	2.450000e+01	4.600000e+01	1.390000e+02
BL	core winter range	Unburned	Fenced	utm_e_nad83	446665	4.481500e+05	4.483695e+05	4.495780e+05	4.500990e+05
BL	core winter range	Unburned	Fenced	utm_n_nad83	4467680	4.468965e+06	4.471327e+06	4.472678e+06	4.472916e+06
BL	core winter range	Unburned	Fenced	year	2007	2.007000e+03	2.008000e+03	2.009000e+03	2.009000e+03
BL	core winter range	Burned	Unfenced	percent_tot	0	0.000000e+00	3.333333e-01	3.333333e-01	3.333333e-01
BL	core winter range	Burned	Unfenced	site_number	27	2.700000e+01	2.700000e+01	2.700000e+01	2.700000e+01
BL	core winter range	Burned	Unfenced	stem_den_ac	0	0.000000e+00	1.618744e+02	1.618744e+02	1.618744e+02
BL	core winter range	Burned	Unfenced	stem_den_ha	0	0.000000e+00	4.000000e+02	4.000000e+02	4.000000e+02
BL	core winter range	Burned	Unfenced	stem_tally	0	0.000000e+00	1.000000e+00	1.000000e+00	1.000000e+00
BL	core winter range	Burned	Unfenced	sum_stem	3	3.000000e+00	3.000000e+00	3.000000e+00	3.000000e+00
BL	core winter range	Burned	Unfenced	utm_e_nad83	446615	4.466150e+05	4.466150e+05	4.466150e+05	4.466150e+05
BL	core winter range	Burned	Unfenced	utm_n_nad83	4467326	4.467326e+06	4.467326e+06	4.467326e+06	4.467326e+06
BL	core winter range	Burned	Unfenced	year	2009	2.009000e+03	2.009000e+03	2.009000e+03	2.009000e+03
BL	core winter range	Burned	Fenced	percent_tot	0	0.000000e+00	0.000000e+00	2.187500e-01	1.000000e+00
BL	core winter range	Burned	Fenced	site_number	10	1.000000e+01	1.450000e+01	1.900000e+01	1.900000e+01
BL	core winter range	Burned	Fenced	stem_den_ac	0	0.000000e+00	0.000000e+00	1.133121e+03	4.370609e+03
BL	core winter range	Burned	Fenced	stem_den_ha	0	0.000000e+00	0.000000e+00	2.800000e+03	1.080000e+04
BL	core winter range	Burned	Fenced	stem_tally	0	0.000000e+00	0.000000e+00	7.000000e+00	2.700000e+01
BL	core winter range	Burned	Fenced	sum_stem	27	2.700000e+01	2.950000e+01	3.200000e+01	3.200000e+01
BL	core winter range	Burned	Fenced	utm_e_nad83	446929	4.469290e+05	4.470120e+05	4.470950e+05	4.470950e+05
BL	core winter range	Burned	Fenced	utm_n_nad83	4467463	4.467463e+06	4.467476e+06	4.467489e+06	4.467489e+06
BL	core winter range	Burned	Fenced	year	2006	2.006000e+03	2.006000e+03	2.006000e+03	2.006000e+03
2013	core winter range	Unburned	Unfenced	percent_tot	0	0.000000e+00	0.000000e+00	3.636364e-01	1.000000e+00
2013	core winter range	Unburned	Unfenced	site_number	2	9.000000e+00	2.000000e+01	3.000000e+01	4.000000e+01
2013	core winter range	Unburned	Unfenced	stem_den_ac	0	0.000000e+00	0.000000e+00	6.474976e+02	2.023430e+04
2013	core winter range	Unburned	Unfenced	stem_den_ha	0	0.000000e+00	0.000000e+00	1.600000e+03	5.000000e+04
2013	core winter range	Unburned	Unfenced	stem_tally	0	0.000000e+00	0.000000e+00	4.000000e+00	1.250000e+02
2013	core winter range	Unburned	Unfenced	sum_stem	0	7.000000e+00	2.200000e+01	3.400000e+01	1.440000e+02
2013	core winter range	Unburned	Unfenced	utm_e_nad83	445348	4.464030e+05	4.465360e+05	4.473970e+05	4.505480e+05
2013	core winter range	Unburned	Unfenced	utm_n_nad83	4467366	4.468377e+06	4.472409e+06	4.473033e+06	4.473650e+06
2013	core winter range	Unburned	Unfenced	year	2013	2.013000e+03	2.013000e+03	2.013000e+03	2.013000e+03
2013	core winter range	Unburned	Fenced	percent_tot	0	2.884615e-02	1.666667e-01	2.788462e-01	9.411765e-01
2013	core winter range	Unburned	Fenced	site_number	8	3.100000e+01	4.400000e+01	6.100000e+01	6.700000e+01
2013	core winter range	Unburned	Fenced	stem_den_ac	0	1.618744e+02	4.856232e+02	1.456870e+03	1.230245e+04
2013	core winter range	Unburned	Fenced	stem_den_ha	0	4.000000e+02	1.200000e+03	3.600000e+03	3.040000e+04
2013	core winter range	Unburned	Fenced	stem_tally	0	1.000000e+00	3.000000e+00	9.000000e+00	7.600000e+01
2013	core winter range	Unburned	Fenced	sum_stem	4	1.000000e+01	2.800000e+01	7.000000e+01	1.040000e+02
2013	core winter range	Unburned	Fenced	utm_e_nad83	446665	4.477950e+05	4.483260e+05	4.495780e+05	4.500990e+05
2013	core winter range	Unburned	Fenced	utm_n_nad83	4467489	4.468593e+06	4.471308e+06	4.472678e+06	4.472916e+06
2013	core winter range	Unburned	Fenced	year	2013	2.013000e+03	2.013000e+03	2.013000e+03	2.013000e+03
2013	core winter range	Burned	Unfenced	percent_tot	0	0.000000e+00	3.721374e-02	4.279044e-01	8.928571e-01
2013	core winter range	Burned	Unfenced	site_number	27	4.500000e+01	6.350000e+01	6.450000e+01	6.500000e+01
2013	core winter range	Burned	Unfenced	stem_den_ac	0	0.000000e+00	2.428116e+02	9.064966e+03	4.305859e+04
2013	core winter range	Burned	Unfenced	stem_den_ha	0	0.000000e+00	6.000000e+02	2.240000e+04	1.064000e+05
2013	core winter range	Burned	Unfenced	stem_tally	0	0.000000e+00	1.500000e+00	5.600000e+01	2.660000e+02
2013	core winter range	Burned	Unfenced	sum_stem	17	9.250000e+01	1.685000e+02	3.465000e+02	5.240000e+02
2013	core winter range	Burned	Unfenced	utm_e_nad83	446545	4.465800e+05	4.468505e+05	4.471050e+05	4.471240e+05
2013	core winter range	Burned	Unfenced	utm_n_nad83	4467326	4.467333e+06	4.467378e+06	4.467427e+06	4.467439e+06
2013	core winter range	Burned	Unfenced	year	2013	2.013000e+03	2.013000e+03	2.013000e+03	2.013000e+03
2013	core winter range	Burned	Fenced	percent_tot	0	0.000000e+00	0.000000e+00	4.800000e-01	6.346154e-01
2013	core winter range	Burned	Fenced	site_number	10	1.000000e+01	1.900000e+01	6.600000e+01	6.600000e+01
2013	core winter range	Burned	Fenced	stem_den_ac	0	0.000000e+00	0.000000e+00	3.723111e+03	5.341855e+03
2013	core winter range	Burned	Fenced	stem_den_ha	0	0.000000e+00	0.000000e+00	9.200000e+03	1.320000e+04
2013	core winter range	Burned	Fenced	stem_tally	0	0.000000e+00	0.000000e+00	2.300000e+01	3.300000e+01
2013	core winter range	Burned	Fenced	sum_stem	44	4.400000e+01	5.000000e+01	5.200000e+01	5.200000e+01
2013	core winter range	Burned	Fenced	utm_e_nad83	446929	4.469290e+05	4.470050e+05	4.470950e+05	4.470950e+05
2013	core winter range	Burned	Fenced	utm_n_nad83	4467463	4.467463e+06	4.467489e+06	4.467494e+06	4.467494e+06
2013	core winter range	Burned	Fenced	year	2013	2.013000e+03	2.013000e+03	2.013000e+03	2.013000e+03
2018	core winter range	Unburned	Unfenced	percent_tot	0	0.000000e+00	5.681818e-03	3.809524e-01	1.000000e+00
2018	core winter range	Unburned	Unfenced	site_number	2	1.000000e+01	2.100000e+01	3.100000e+01	6.800000e+01
2018	core winter range	Unburned	Unfenced	stem_den_ac	0	0.000000e+00	0.000000e+00	9.712464e+02	1.214058e+04
2018	core winter range	Unburned	Unfenced	stem_den_ha	0	0.000000e+00	0.000000e+00	2.400000e+03	3.000000e+04
2018	core winter range	Unburned	Unfenced	stem_tally	0	0.000000e+00	0.000000e+00	6.000000e+00	7.500000e+01
2018	core winter range	Unburned	Unfenced	sum_stem	0	6.500000e+00	1.800000e+01	4.650000e+01	9.500000e+01
2018	core winter range	Unburned	Unfenced	utm_e_nad83	445348	4.464100e+05	4.465395e+05	4.474950e+05	4.505480e+05
2018	core winter range	Unburned	Unfenced	utm_n_nad83	4467366	4.468958e+06	4.472412e+06	4.473022e+06	4.473650e+06
2018	core winter range	Unburned	Unfenced	year	2018	2.018000e+03	2.018000e+03	2.018000e+03	2.018000e+03
2018	core winter range	Unburned	Fenced	percent_tot	0	0.000000e+00	2.000000e-01	3.333333e-01	6.666667e-01
2018	core winter range	Unburned	Fenced	site_number	8	3.100000e+01	4.400000e+01	6.100000e+01	6.700000e+01
2018	core winter range	Unburned	Fenced	stem_den_ac	0	0.000000e+00	3.237488e+02	6.474976e+02	3.723111e+03
2018	core winter range	Unburned	Fenced	stem_den_ha	0	0.000000e+00	8.000000e+02	1.600000e+03	9.200000e+03
2018	core winter range	Unburned	Fenced	stem_tally	0	0.000000e+00	2.000000e+00	4.000000e+00	2.300000e+01
2018	core winter range	Unburned	Fenced	sum_stem	3	4.000000e+00	9.000000e+00	2.100000e+01	5.700000e+01
2018	core winter range	Unburned	Fenced	utm_e_nad83	446665	4.477950e+05	4.483260e+05	4.495780e+05	4.500990e+05
2018	core winter range	Unburned	Fenced	utm_n_nad83	4467489	4.468593e+06	4.471308e+06	4.472678e+06	4.472916e+06
2018	core winter range	Unburned	Fenced	year	2018	2.018000e+03	2.018000e+03	2.018000e+03	2.018000e+03
2018	core winter range	Burned	Unfenced	percent_tot	0	0.000000e+00	1.329911e-01	4.083072e-01	5.086207e-01
2018	core winter range	Burned	Unfenced	site_number	27	4.500000e+01	6.350000e+01	6.450000e+01	6.500000e+01
2018	core winter range	Burned	Unfenced	stem_den_ac	0	0.000000e+00	7.284348e+02	6.394039e+03	9.712464e+03
2018	core winter range	Burned	Unfenced	stem_den_ha	0	0.000000e+00	1.800000e+03	1.580000e+04	2.400000e+04
2018	core winter range	Burned	Unfenced	stem_tally	0	0.000000e+00	4.500000e+00	3.950000e+01	6.000000e+01
2018	core winter range	Burned	Unfenced	sum_stem	2	3.950000e+01	9.650000e+01	1.470000e+02	1.780000e+02
2018	core winter range	Burned	Unfenced	utm_e_nad83	446545	4.465800e+05	4.468505e+05	4.471050e+05	4.471240e+05
2018	core winter range	Burned	Unfenced	utm_n_nad83	4467326	4.467333e+06	4.467378e+06	4.467427e+06	4.467439e+06
2018	core winter range	Burned	Unfenced	year	2018	2.018000e+03	2.018000e+03	2.018000e+03	2.018000e+03
2018	core winter range	Burned	Fenced	percent_tot	0	7.692308e-02	2.000000e-01	2.857143e-01	5.000000e-01
2018	core winter range	Burned	Fenced	site_number	10	1.000000e+01	1.900000e+01	6.600000e+01	6.600000e+01
2018	core winter range	Burned	Fenced	stem_den_ac	0	3.237488e+02	8.093720e+02	1.294995e+03	2.104367e+03
2018	core winter range	Burned	Fenced	stem_den_ha	0	8.000000e+02	2.000000e+03	3.200000e+03	5.200000e+03
2018	core winter range	Burned	Fenced	stem_tally	0	2.000000e+00	5.000000e+00	8.000000e+00	1.300000e+01
2018	core winter range	Burned	Fenced	sum_stem	14	1.400000e+01	2.500000e+01	3.900000e+01	3.900000e+01
2018	core winter range	Burned	Fenced	utm_e_nad83	446929	4.469290e+05	4.470050e+05	4.470950e+05	4.470950e+05
2018	core winter range	Burned	Fenced	utm_n_nad83	4467463	4.467463e+06	4.467489e+06	4.467494e+06	4.467494e+06
2018	core winter range	Burned	Fenced	year	2018	2.018000e+03	2.018000e+03	2.018000e+03	2.018000e+03

Pairwise contrasts
Results are given on the response scale
contrast	burned	fenced	estimate	lower.HPD	upper.HPD
2013 - BL	Unburned	Unfenced	1.4477160	1.2918801	1.6172958
2018 - BL	Unburned	Unfenced	1.4501550	1.2964982	1.6266971
2018 - 2013	Unburned	Unfenced	1.0015036	0.9050915	1.1108198
2013 - BL	Burned	Unfenced	5.3623899	1.6517603	20.5791732
2018 - BL	Burned	Unfenced	2.2802934	0.6827064	8.6486092
2018 - 2013	Burned	Unfenced	0.4247683	0.3764684	0.4805162
2013 - BL	Unburned	Fenced	0.7912445	0.6862436	0.9046157
2018 - BL	Unburned	Fenced	0.3294917	0.2735920	0.3908882
2018 - 2013	Unburned	Fenced	0.4162960	0.3487131	0.4960705
2013 - BL	Burned	Fenced	1.8359081	1.3606438	2.5357842
2018 - BL	Burned	Fenced	0.9804953	0.6969817	1.4013976
2018 - 2013	Burned	Fenced	0.5334445	0.4072155	0.7056563

All Stem Diameters - Noncore Winter Range (ANC)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class + (1 | site_id)

## Priors for model 'stmod_anc_stally1' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.27,5.27])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	0.668	0.004	0.175	0.444	0.669	0.888	1700	1.005
time_class2013	0.204	0.001	0.084	0.095	0.204	0.311	10588	1.000
time_class2018	0.605	0.001	0.077	0.505	0.604	0.703	10520	1.000
mean_PPD	3.421	0.001	0.148	3.232	3.419	3.610	13629	1.000
log-posterior	-1093.668	0.097	4.799	-1099.894	-1093.304	-1087.825	2455	1.001

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	estimate	lower.HPD	upper.HPD
2013 - BL	0.2	0.04	0.37
2018 - BL	0.6	0.46	0.76
2018 - 2013	0.4	0.27	0.55

Pairwise contrasts
Results are given on the response scale.
contrast	ratio	lower.HPD	upper.HPD
2013 / BL	1.23	1.03	1.44
2018 / BL	1.83	1.56	2.12
2018 / 2013	1.49	1.29	1.71

# comment 2.2 "provide analyses on single graphic...
## Plot combinations
contr.plot.ac <- emmeans(stmod_stally2, ~ time_class | fenced) %>% 
  pairs(reverse = TRUE) %>% 
  # plot(type = "response") +
  plot() +
  theme_minimal() +
  geom_vline(aes(xintercept = 0), color= "black", size=1, lty="dashed") +
  labs(title= "Core winter range", caption = "", x = "Estimate", y = "Contrast")

contr.plot.anc <- emmeans(stmod_anc_stally1, ~ time_class) %>% 
  pairs(reverse = TRUE) %>% 
  plot(type = "response") +
  theme_minimal() +
  geom_vline(aes(xintercept = 0), color= "black", size=1, lty="dashed") +
  labs(title= "Noncore winter range", caption = "Point estimate displayed: median 
HPD interval probability: 0.95, Results on log scale", x = "Estimate", y = "")

# patchwork
contr.plot.ac + contr.plot.anc 

ggsave("./output/figures_exported/revised_figs_aspen/emm_stemcnt_AC_ANC.png", width = 8, height = 3.75, dpi=300)
ggsave("./output/figures_exported/revised_figs_aspen/emm_stemcnt_AC_ANC.pdf", width = 8, height = 3.75)

All Stem Diameters - Kawuneeche Valley (AK)

Modeling and Posterior Description

Prior summary: stem_tally ~ time_class + (1 | site_id)

## Priors for model 'stmod_ak_stally1' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.28,5.28])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Parameter	mean	mcse	sd	10%	50%	90%	n_eff	Rhat
(Intercept)	0.567	0.004	0.246	0.257	0.571	0.866	3931	1
time_class2013	0.191	0.002	0.152	-0.004	0.192	0.384	10000	1
time_class2018	0.377	0.001	0.144	0.193	0.376	0.561	9863	1
mean_PPD	2.475	0.002	0.202	2.217	2.467	2.733	14246	1
log-posterior	-504.307	0.049	3.034	-508.326	-503.985	-500.694	3768	1

## Priors for model 'stmod_ak_stally1' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.28,5.28])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	estimate	lower.HPD	upper.HPD
2013 - BL	0.19	-0.11	0.49
2018 - BL	0.38	0.09	0.66
2018 - 2013	0.18	-0.07	0.47

Pairwise contrasts
Results are given on the response scale.
contrast	ratio	lower.HPD	upper.HPD
2013 / BL	1.21	0.88	1.60
2018 / BL	1.46	1.06	1.89
2018 / 2013	1.20	0.90	1.55

Willow Macroplot Data

Willow Height - Combined Winter Range Macroplot Data

Histograms of Willow Height

Modeling and Posterior Description

Prior summary: plant_ht_cm ~ time_class*fenced + (1 | site_id)

# Get the prior summary
prior_summary(wc.wnc.stmod_TCxF)

## Priors for model 'wc.wnc.stmod_TCxF' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.70,5.00,5.75,...])
## 
## Auxiliary (shape)
##  ~ exponential(rate = 1)
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.04	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.00	fixed	conditional
fencedFenced	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013:fencedFenced	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2018:fencedFenced	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

## fenced = Unfenced:
##  contrast    estimate lower.HPD upper.HPD
##  2013 - BL     0.1364    0.0825     0.186
##  2018 - BL     0.2166    0.1715     0.263
##  2018 - 2013   0.0801    0.0367     0.127
## 
## fenced = Fenced:
##  contrast    estimate lower.HPD upper.HPD
##  2013 - BL     0.4719    0.3745     0.566
##  2018 - BL     0.9345    0.8500     1.023
##  2018 - 2013   0.4627    0.3857     0.540
## 
## Point estimate displayed: median 
## Results are given on the log (not the response) scale. 
## HPD interval probability: 0.95

Pairwise contrasts
Results are given on the log scale
contrast	fenced	estimate	lower.HPD	upper.HPD
2013 - BL	Unfenced	0.14	0.08	0.19
2018 - BL	Unfenced	0.22	0.17	0.26
2018 - 2013	Unfenced	0.08	0.04	0.13
2013 - BL	Fenced	0.47	0.37	0.57
2018 - BL	Fenced	0.93	0.85	1.02
2018 - 2013	Fenced	0.46	0.39	0.54

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
2013 / BL	Unfenced	1.15	1.09	1.20
2018 / BL	Unfenced	1.24	1.19	1.30
2018 / 2013	Unfenced	1.08	1.03	1.13
2013 / BL	Fenced	1.60	1.45	1.76
2018 / BL	Fenced	2.55	2.33	2.77
2018 / 2013	Fenced	1.59	1.47	1.71

EMMIP

#### Figure 29 (#1) revised 2022
plotF.emmip_FxTC_response_bw(wc.wnc.stmod_TCxF) +
  labs(y="Willow height (cm)") +
  theme(legend.position = "right")

# ggsave("./output/figures_202202/emmip/Fig29_1_emmip_TCxF_wcwnc_allwillowheight.png", width = 7.25, height = 3.75) 
# 
# ggsave("./output/figures_202202/emmip/Fig29_1_emmip_TCxF_wcwnc_allwillowheight.pdf", width = 7.25, height = 3.75) 

# bw
plotF.emmip_FxTC_response_bw(wc.wnc.stmod_TCxF) +
  labs(title="Combined winter range") +
  theme(legend.position = "right") +
  labs(x="Year")

ggsave("./output/figures_exported/emmip_TCxF_wcwnc_allwillowheight_bw.png", width = 5, height = 4.25) 

Effects of Burning

Prior summary height ~ time_class * burned * fenced + (1 | site_id)

## Priors for model 'stmod_ht.wc.wnc.TCxFxB' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.70,5.00,5.75,...])
## 
## Auxiliary (shape)
##  ~ exponential(rate = 1)
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

## fenced = Unfenced, burned = Burned:
##  time_class emmean lower.HPD upper.HPD
##  BL           4.26      3.88      4.64
##  2013         3.85      3.47      4.20
##  2018         4.07      3.72      4.45
## 
## fenced = Fenced, burned = Burned:
##  time_class emmean lower.HPD upper.HPD
##  BL           3.79      3.48      4.13
##  2013         4.08      3.79      4.40
##  2018         4.81      4.51      5.10
## 
## fenced = Unfenced, burned = Unburned:
##  time_class emmean lower.HPD upper.HPD
##  BL           4.70      4.55      4.84
##  2013         4.85      4.71      5.00
##  2018         4.93      4.77      5.06
## 
## fenced = Fenced, burned = Unburned:
##  time_class emmean lower.HPD upper.HPD
##  BL           4.25      3.96      4.54
##  2013         4.85      4.55      5.13
##  2018         5.11      4.83      5.39
## 
## Point estimate displayed: median 
## Results are given on the log (not the response) scale. 
## HPD interval probability: 0.95

Pairwise contrasts
Results are given on the log scale
contrast	fenced	estimate	lower.HPD	upper.HPD
2013 - BL	Unfenced	-0.13	-0.27	0.01
2018 - BL	Unfenced	0.02	-0.12	0.16
2018 - 2013	Unfenced	0.15	0.04	0.27
2013 - BL	Fenced	0.45	0.34	0.55
2018 - BL	Fenced	0.94	0.84	1.04
2018 - 2013	Fenced	0.50	0.42	0.57

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
2013 / BL	Unfenced	0.88	0.76	1.01
2018 / BL	Unfenced	1.02	0.88	1.17
2018 / 2013	Unfenced	1.16	1.03	1.30
2013 / BL	Fenced	1.56	1.40	1.73
2018 / BL	Fenced	2.56	2.31	2.82
2018 / 2013	Fenced	1.64	1.52	1.77

EMMIP

#### Figure 40 (#1)

Willow Height - Core Winter Range Macroplot Data

Histograms of Willow Height

Modeling and Posterior Description

Prior summary: plant_ht_cm ~ time_class*fenced + (1 | site_id)

## Priors for model 'stmod_ht2' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.65,5.01,5.06,...])
## 
## Auxiliary (shape)
##  ~ exponential(rate = 1)
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

## model comparison.
## Run originally and code retained. stmod_ht2 selected
loo1 <- loo(stmod_ht1,
            k_threshold = 0.7) 
loo2 <- loo(stmod_ht2,
            k_threshold = 0.7) 
comp <- loo_compare(loo2, loo1)
print(comp, simplify = FALSE, digits = 2)

Core Winter Range
Willow Height Posterior Description
Parameter	Component	Median	MAD	CI	CI_low	CI_high	pd	ps	ROPE_CI	ROPE_low	ROPE_high	ROPE_Percentage	Rhat	ESS
(Intercept)	conditional	4.4	0.1	0.9	4.3	4.6	1.0	1.0	0.9	-0.1	0.1	0.0	1	997.1
time_class2013	conditional	0.0	0.0	0.9	-0.1	0.1	0.6	0.0	0.9	-0.1	0.1	1.0	1	5191.9
time_class2018	conditional	0.1	0.0	0.9	0.0	0.2	1.0	0.6	0.9	-0.1	0.1	0.4	1	5152.4
fencedFenced	conditional	-0.4	0.1	0.9	-0.6	-0.2	1.0	1.0	0.9	-0.1	0.1	0.0	1	1015.9
time_class2013:fencedFenced	conditional	0.5	0.1	0.9	0.4	0.6	1.0	1.0	0.9	-0.1	0.1	0.0	1	5281.0
time_class2018:fencedFenced	conditional	0.8	0.1	0.9	0.7	0.9	1.0	1.0	0.9	-0.1	0.1	0.0	1	4898.9
shape	distributional	4.7	0.1	0.9	4.4	4.9	1.0	1.0	0.9	-0.1	0.1	0.0	NA	NA

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2013	0.9	-0.1	0.1	1.0	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.4	fixed	conditional
fencedFenced	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2013:fencedFenced	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2018:fencedFenced	0.9	-0.1	0.1	0.0	fixed	conditional

Contrasts

Results are given on the response scale.

Pairwise contrasts
Results are given on the log scale
contrast	fenced	estimate	lower.HPD	upper.HPD
2013 - BL	Unfenced	0.01	-0.07	0.08
2018 - BL	Unfenced	0.11	0.03	0.18
2018 - 2013	Unfenced	0.10	0.04	0.16
2013 - BL	Fenced	0.47	0.38	0.56
2018 - BL	Fenced	0.93	0.85	1.01
2018 - 2013	Fenced	0.46	0.39	0.53

Pairwise contrasts
Results are given on the response scale.
contrast	fenced	ratio	lower.HPD	upper.HPD
2013 / BL	Unfenced	1.01	0.93	1.08
2018 / BL	Unfenced	1.11	1.03	1.19
2018 / 2013	Unfenced	1.10	1.04	1.17
2013 / BL	Fenced	1.60	1.46	1.75
2018 / BL	Fenced	2.55	2.33	2.75
2018 / 2013	Fenced	1.59	1.48	1.70

plotF.emmip_FxTC_response(stmod_ht2) +
  labs(title="Core winter range") +
  theme(legend.position = "right")

ggsave("./output/figures_exported/emmip_TCxF_wc_willowht.png", width = 5, height = 4.25) # save plot

Figure 29 (#2) revised 2022

plotF.emmip_FxTC_response_bw(stmod_ht2) +
  labs(y="Willow height (cm)") +
  theme(legend.position = "right")

ggsave("./output/figures_202202/emmip/Fig29_2_emmip_TCxF_wc_allwillowheight.png", width = 7.25, height = 3.75) 

ggsave("./output/figures_202202/emmip/Fig29_2_emmip_TCxF_wc_allwillowheight.pdf", width = 7.25, height = 3.75)

Willow Height - Noncore Winter Range Macroplot Data

Histograms of Willow Height

Modeling and Posterior Description

Prior summary: plant_ht_cm ~ time_class + (1 | site_id)

## Priors for model 'nc.stmod_ht1' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.77,5.02])
## 
## Auxiliary (shape)
##  ~ exponential(rate = 1)
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.0	fixed	conditional

Contrasts

Results are given on the the response scale.

Pairwise contrasts
Results are given on the log scale
contrast	estimate	lower.HPD	upper.HPD
2013 - BL	0.22	0.15	0.29
2018 - BL	0.28	0.22	0.34
2018 - 2013	0.06	-0.01	0.13

Pairwise contrasts
Results are given on the response scale.
contrast	ratio	lower.HPD	upper.HPD
2013 / BL	1.25	1.16	1.34
2018 / BL	1.32	1.24	1.41
2018 / 2013	1.06	0.99	1.13

Figure 30 revised 2022

plotF.emmip_TC_response_bw(nc.stmod_ht1) +
  labs(y="Willow height (cm)") +
  theme(legend.position = "right")

ggsave("./output/figures_202202/emmip/Fig30_emmip_TC_wnc_allwillowheight.png", width = 7.25, height = 3.75) 

ggsave("./output/figures_202202/emmip/Fig30_emmip_TC_wnc_allwillowheight.pdf", width = 7.25, height = 3.75) 

Willow Height - Kawuneeche Valley Macroplot Data

Histograms of Willow Height

time_class	variable	mean	sd	min	q1	med	q3	max	mad	iqr	cv	skewness	se.skewness	kurtosis	n.valid	pct.valid	sem
2011	plant_ht_cm	47.65323	23.42608	14	30	43	59.5	168	20.7564	29.25	0.4915948	1.6734970	0.2173713	5.3253116	124	100	2.103724
2016	plant_ht_cm	56.40741	25.88473	10	40	50	75.0	155	22.2390	32.50	0.4588889	0.8201672	0.2085259	0.8220138	135	100	2.227803
2017	plant_ht_cm	55.39568	29.36567	10	35	50	70.0	145	22.2390	35.00	0.5301075	0.9429452	0.2055673	0.5114885	139	100	2.490763
2018	plant_ht_cm	60.51383	34.24835	5	35	55	85.0	150	37.0650	50.00	0.5659590	0.4952152	0.1530955	-0.6332637	253	100	2.153175

Modeling and Posterior Description

Prior summary: plant_ht_cm ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_ht_wk_TCxF' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [6.16,6.10,5.12,...])
## 
## Auxiliary (shape)
##  ~ exponential(rate = 1)
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

#### PP check plot
ppc.plot(stmod_ht_wk_TCxF) +
  labs(caption = "WK: willow_ht~TCxF")

ggsave("./output/figures_202108/willow_height_stats/ppc_willowht_wk.png", dpi = 300, width = 4.75, height = 3.75)

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2016	0.9	-0.1	0.1	0.3	fixed	conditional
time_class2017	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.0	fixed	conditional
fencedFenced	0.9	-0.1	0.1	0.3	fixed	conditional
time_class2016:fencedFenced	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2017:fencedFenced	0.9	-0.1	0.1	0.0	fixed	conditional
time_class2018:fencedFenced	0.9	-0.1	0.1	0.0	fixed	conditional

Contrasts

## time_class = 2011:
##  fenced   emmean lower.HPD upper.HPD
##  Unfenced   4.01      3.68      4.37
##  Fenced     3.87      3.43      4.29
## 
## time_class = 2016:
##  fenced   emmean lower.HPD upper.HPD
##  Unfenced   3.87      3.53      4.22
##  Fenced     4.31      3.89      4.73
## 
## time_class = 2017:
##  fenced   emmean lower.HPD upper.HPD
##  Unfenced   3.80      3.46      4.15
##  Fenced     4.42      4.00      4.85
## 
## time_class = 2018:
##  fenced   emmean lower.HPD upper.HPD
##  Unfenced   3.56      3.24      3.92
##  Fenced     4.55      4.14      4.97
## 
## Point estimate displayed: median 
## Results are given on the log (not the response) scale. 
## HPD interval probability: 0.95

Pairwise contrasts
Results are given on the log (not the response) scale.
contrast	fenced	estimate	lower.HPD	upper.HPD
2016 - 2011	Unfenced	-0.14	-0.27	0.00
2017 - 2011	Unfenced	-0.21	-0.34	-0.08
2017 - 2016	Unfenced	-0.07	-0.20	0.06
2018 - 2011	Unfenced	-0.45	-0.58	-0.32
2018 - 2016	Unfenced	-0.31	-0.43	-0.19
2018 - 2017	Unfenced	-0.24	-0.36	-0.13
2016 - 2011	Fenced	0.44	0.27	0.61
2017 - 2011	Fenced	0.55	0.38	0.73
2017 - 2016	Fenced	0.11	-0.05	0.26
2018 - 2011	Fenced	0.69	0.52	0.85
2018 - 2016	Fenced	0.24	0.10	0.39
2018 - 2017	Fenced	0.13	-0.02	0.28

Pairwise contrasts
Results are given on the response scale.
contrast	ratio	lower.HPD	upper.HPD
2016 / 2011	1.16	1.05	1.30
2017 / 2011	1.19	1.06	1.32
2017 / 2016	1.02	0.92	1.13
2018 / 2011	1.12	1.01	1.24
2018 / 2016	0.97	0.87	1.05
2018 / 2017	0.95	0.86	1.04

#### Figure 41

COVER

Willow Cover - Combined Winter Range Macroplot Data

# munge
canopy.all.willow.sum <- canopy.all.willow.sum %>% 
  mutate(burned = case_when(burned == "Y" ~ "Burned",
                            burned == "N" ~ "Unburned",
                            is.na(burned) ~ "Unburned",
                            TRUE ~ burned)) %>% 
  mutate(fenced = case_when(site_id == "WK01" ~ "Fenced",
                            site_id == "WK02" ~ "Fenced",
                            site_id == "WK03" ~ "Unfenced",
                            site_id == "WK04" ~ "Unfenced",
                            site_id == "WK05" ~ "Fenced",
                            site_id == "WK06" ~ "Unfenced",
                            site_id == "WK07" ~ "Unfenced",
                            site_id == "WK08" ~ "Unfenced",
                            site_id == "WK09" ~ "Unfenced",
                            site_id == "WK10" ~ "Fenced",
                            TRUE ~ fenced)) %>% 
  mutate(time_class = as_factor(time_class)) %>%
  mutate(fenced = as_factor(fenced)) %>%
  mutate(burned = as_factor(burned)) %>% 
  mutate(time_class = fct_rev(time_class))

Combined WC and WNC macroplot data: willow cover - fenced + time_class
time_class	fenced	variable	mean	sd	min	med	max	n.valid	pct.valid
2013	Unfenced	cover_allwillow	0.31	0.30	0.00	0.25	1	60	100
2013	Fenced	cover_allwillow	0.29	0.30	0.00	0.18	1	28	100
2018	Unfenced	cover_allwillow	0.45	0.35	0.00	0.41	1	66	100
2018	Fenced	cover_allwillow	0.55	0.35	0.05	0.56	1	28	100
BL	Unfenced	cover_allwillow	0.33	0.27	0.00	0.27	1	51	100
BL	Fenced	cover_allwillow	0.17	0.24	0.00	0.11	1	17	100

Combined WC and WNC macroplot data: willow cover ~ burned + fenced + time_class
time_class	fenced	burned	variable	mean	sd	min	med	max	n.valid	pct.valid
2013	Unfenced	Unburned	cover_allwillow	0.34	0.29	0.00	0.27	1.00	51	100
2013	Unfenced	Burned	cover_allwillow	0.11	0.23	0.01	0.02	0.72	9	100
2013	Fenced	Unburned	cover_allwillow	0.43	0.30	0.02	0.48	1.00	15	100
2013	Fenced	Burned	cover_allwillow	0.12	0.21	0.00	0.04	0.79	13	100
2018	Unfenced	Unburned	cover_allwillow	0.49	0.33	0.00	0.47	1.00	58	100
2018	Unfenced	Burned	cover_allwillow	0.15	0.34	0.00	0.03	1.00	8	100
2018	Fenced	Unburned	cover_allwillow	0.69	0.33	0.14	0.76	1.00	15	100
2018	Fenced	Burned	cover_allwillow	0.39	0.29	0.05	0.26	1.00	13	100
BL	Unfenced	Unburned	cover_allwillow	0.34	0.27	0.00	0.27	1.00	50	100
BL	Unfenced	Burned	cover_allwillow	0.03	NA	0.03	0.03	0.03	1	100
BL	Fenced	Unburned	cover_allwillow	0.27	0.27	0.01	0.20	1.00	10	100
BL	Fenced	Burned	cover_allwillow	0.03	0.03	0.00	0.02	0.09	7	100

WK macroplot data: willow cover ~ fenced * time_class
time_class	fenced	variable	mean	sd	min	med	max	n.valid	pct.valid
2016	Unfenced	cover_allwillow	0.23	0.23	0.04	0.18	0.54	4	100
2016	Fenced	cover_allwillow	0.31	0.20	0.19	0.20	0.54	3	100
2018	Unfenced	cover_allwillow	0.26	0.25	0.02	0.17	0.53	5	100
2018	Fenced	cover_allwillow	0.46	0.27	0.14	0.48	0.75	4	100
2017	Unfenced	cover_allwillow	0.23	0.24	0.02	0.16	0.58	4	100
2017	Fenced	cover_allwillow	0.34	0.25	0.17	0.23	0.62	3	100
BL	Unfenced	cover_allwillow	0.20	0.17	0.04	0.17	0.42	5	100
BL	Fenced	cover_allwillow	0.13	0.13	0.04	0.06	0.28	3	100

Modeling and Posterior Description

Prior summary: cover_allwillow ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_cov.wc.wnc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.15,5.61,5.49,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log odds ratio, not the response scale
contrast	fenced	estimate	lower.HPD	upper.HPD
2018 - 2013	Unfenced	0.47	0.13	0.81
BL - 2013	Unfenced	-0.34	-0.74	0.06
BL - 2018	Unfenced	-0.81	-1.20	-0.40
2018 - 2013	Fenced	1.97	1.44	2.48
BL - 2013	Fenced	-0.66	-1.32	0.01
BL - 2018	Fenced	-2.63	-3.30	-1.95

Pairwise contrasts
Results are given on the response scale
contrast	fenced	odds.ratio	lower.HPD	upper.HPD
2018 / 2013	Unfenced	1.59	1.11	2.22
BL / 2013	Unfenced	0.71	0.46	1.03
BL / 2018	Unfenced	0.45	0.29	0.65
2018 / 2013	Fenced	7.20	3.93	11.56
BL / 2013	Fenced	0.52	0.23	0.93
BL / 2018	Fenced	0.07	0.03	0.13

EMMIP

#### Fig 39 - 1

Effects of Burning - Combined Range

Prior summary cover_allwillow ~ time_class * burned * fenced + (1 | site_id)

## Priors for model 'stmod_cov.wc.wnc.TCxFxB' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.15,5.61,6.19,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

EMMIP

Figure 40 (#2)

Willow Cover - Core Winter Range Macroplot Data

WC macroplot data: willow cover ~ WC: all willow cover - time_class
time_class	variable	mean	sd	min	med	max	n.valid	pct.valid
2013	cover_allwillow	0.24	0.27	0	0.14	1	61	100
2018	cover_allwillow	0.44	0.35	0	0.33	1	62	100
BL	cover_allwillow	0.22	0.24	0	0.18	1	36	100

WC macroplot data: All willow cover - burned + time_class
time_class	burned	variable	mean	sd	min	med	max	n.valid	pct.valid
2013	Unburned	cover_allwillow	0.31	0.27	0	0.27	1.00	39	100
2013	Burned	cover_allwillow	0.12	0.22	0	0.03	0.79	22	100
2018	Unburned	cover_allwillow	0.51	0.35	0	0.46	1.00	41	100
2018	Burned	cover_allwillow	0.30	0.33	0	0.20	1.00	21	100
BL	Unburned	cover_allwillow	0.27	0.25	0	0.22	1.00	28	100
BL	Burned	cover_allwillow	0.03	0.03	0	0.02	0.09	8	100

Modeling and Posterior Description

Prior summary: cover_allwillow ~ time_class * fenced + (1 | site_id)

## Priors for model 'stmod_cov.wc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.11,5.95,5.00,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

#### Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Results are given on the the response scale.

Pairwise contrasts
Results are given on the log odds ratio, not the response scale
contrast	fenced	estimate	lower.HPD	upper.HPD
2018 - 2013	Unfenced	0.39	-0.08	0.84
BL - 2013	Unfenced	-0.10	-0.68	0.45
BL - 2018	Unfenced	-0.49	-1.05	0.04
2018 - 2013	Fenced	2.04	1.51	2.56
BL - 2013	Fenced	-0.71	-1.39	-0.08
BL - 2018	Fenced	-2.75	-3.43	-2.07

Pairwise contrasts
Results are given on the response scale
contrast	fenced	odds.ratio	lower.HPD	upper.HPD
2018 / 2013	Unfenced	1.48	0.89	2.27
BL / 2013	Unfenced	0.91	0.48	1.51
BL / 2018	Unfenced	0.61	0.31	0.99
2018 / 2013	Fenced	7.69	4.09	12.12
BL / 2013	Fenced	0.49	0.21	0.86
BL / 2018	Fenced	0.06	0.03	0.11

EMMIP

#### Figure 39 (#2)

Willow Cover - Noncore Winter Range Macroplot Data

WNC macroplot data: willow cover
time_class	fenced	variable	mean	sd	min	med	max	n.valid	pct.valid
2013	Unfenced	cover_allwillow	0.43	0.32	0.02	0.37	1	27	100
2018	Unfenced	cover_allwillow	0.56	0.33	0.01	0.53	1	32	100
BL	Unfenced	cover_allwillow	0.37	0.28	0.03	0.30	1	32	100

Modeling and Posterior Description

Prior summary: cover_allwillow ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.wnc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.21,5.21])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

Pairwise contrasts
Results are given on the log odds ratio, not the response scale
contrast	estimate	lower.HPD	upper.HPD
2018 - 2013	0.49	-0.04	1.05
BL - 2013	-0.54	-1.11	0.05
BL - 2018	-1.03	-1.61	-0.48

Pairwise contrasts
Results are given on the response scale
contrast	odds.ratio	lower.HPD	upper.HPD
2018 / 2013	1.63	0.84	2.63
BL / 2013	0.58	0.29	0.98
BL / 2018	0.36	0.18	0.59

Figure 39 (#3). Modeled mean percent willow canopy cover (points) and credible intervals (bars) for all willow species on the noncore elk winter range of Rocky Mountain National Park. All noncore plots were unfenced. BL = baseline data (collected 2006-2008).

Willow Cover - Kawuneeche Valley Macroplot Data

WNC macroplot data: willow cover
time_class	fenced	variable	mean	sd	min	med	max	n.valid	pct.valid
2016	Unfenced	cover_allwillow	0.23	0.23	0.04	0.18	0.54	4	100
2016	Fenced	cover_allwillow	0.31	0.20	0.19	0.20	0.54	3	100
2018	Unfenced	cover_allwillow	0.26	0.25	0.02	0.17	0.53	5	100
2018	Fenced	cover_allwillow	0.46	0.27	0.14	0.48	0.75	4	100
2017	Unfenced	cover_allwillow	0.23	0.24	0.02	0.16	0.58	4	100
2017	Fenced	cover_allwillow	0.34	0.25	0.17	0.23	0.62	3	100
BL	Unfenced	cover_allwillow	0.20	0.17	0.04	0.17	0.42	5	100
BL	Fenced	cover_allwillow	0.13	0.13	0.04	0.06	0.28	3	100

WNC macroplot data: willow cover
time_class	variable	mean	sd	min	med	max	n.valid	pct.valid
2016	cover_allwillow	0.27	0.20	0.04	0.20	0.54	7	100
2018	cover_allwillow	0.35	0.26	0.02	0.34	0.75	9	100
2017	cover_allwillow	0.28	0.23	0.02	0.21	0.62	7	100
BL	cover_allwillow	0.17	0.15	0.04	0.12	0.42	8	100

Modeling and Posterior Description

## # A tibble: 4 x 1
##   time_class
##   <fct>     
## 1 BL        
## 2 2018      
## 3 2016      
## 4 2017

Prior summary

## Priors for model 'stmod_cov.wk' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0,0,...], scale = [2.5,2.5,2.5,...])
##   Adjusted prior:
##     ~ normal(location = [0,0,0,...], scale = [5.42,5.88,5.62,...])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Contrasts

emmip

Figure 43

Wild Basin Comparisons

location	variable	mean	sd	min	med	max	n.valid
Fishermans Nook	mean	290.0	7.1	285.0	290.0	295.0	2
Hollowell Park	mean	223.6	124.0	83.5	194.9	508.8	12
Wild Basin	mean	216.6	112.5	10.0	210.0	413.3	27
Lower Hidden Valley	mean	216.4	NA	216.4	216.4	216.4	1
McGraw Ranch	mean	183.1	76.6	94.3	170.8	309.2	8
Hidden Valley	mean	160.5	104.8	67.2	116.7	377.5	9
Upper Beaver Meadows	mean	154.7	90.2	55.0	146.6	280.0	8
Horseshoe Park	mean	143.2	63.2	49.0	168.6	226.8	12
Moraine Park	mean	120.1	62.3	30.0	122.0	270.0	33
Endovalley	mean	112.0	69.8	57.0	90.7	287.2	9
Kawuneeche	mean	67.0	32.0	25.0	73.5	111.7	9

site_type	mean	sd	min	med	max	n.valid
WB	216.6	112.5	10.0	210.0	413.3	27
WNC	199.7	104.4	67.2	179.8	508.8	32
WC	127.9	67.4	30.0	118.7	287.2	62
WK	67.0	32.0	25.0	73.5	111.7	9

Modeling and Posterior Description

## stan_aov
##  family:       gaussian [identity]
##  formula:      mean ~ site_type
##  observations: 130
##  predictors:   4
## ------
##              Median MAD_SD
## (Intercept)   214.2   16.7
## site_typeWC   -85.0   19.3
## site_typeWK  -143.3   32.1
## site_typeWNC  -16.3   22.6
## 
## Auxiliary parameter(s):
##               Median MAD_SD
## R2             0.2    0.1  
## log-fit_ratio  0.0    0.1  
## sigma         86.7    5.4  
## 
## ANOVA-like table:
##                   Median  MAD_SD 
## Mean Sq site_type 89489.5 28094.6
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg

## 
## Model Info:
##  function:     stan_lmer
##  family:       gaussian [identity]
##  formula:      mean ~ 1 + (1 | site_type)
##  algorithm:    sampling
##  sample:       14000 (posterior sample size)
##  priors:       see help('prior_summary')
##  observations: 130
##  groups:       site_type (4)
## 
## Estimates:
##                                            mean     sd       10%      50%   
## (Intercept)                                  12.1     39.2     -5.0      0.8
## b[(Intercept) site_type:WB]                 202.4     42.8    168.8    211.2
## b[(Intercept) site_type:WC]                 115.3     40.5     86.0    125.0
## b[(Intercept) site_type:WK]                  53.8     46.4      4.5     60.6
## b[(Intercept) site_type:WNC]                186.3     42.1    153.5    195.2
## sigma                                        87.3      5.6     80.4     87.1
## Sigma[site_type:(Intercept),(Intercept)]  54338.8  59684.2  13351.1  36556.1
##                                            90%   
## (Intercept)                                  37.6
## b[(Intercept) site_type:WB]                 236.3
## b[(Intercept) site_type:WC]                 142.4
## b[(Intercept) site_type:WK]                 100.8
## b[(Intercept) site_type:WNC]                218.3
## sigma                                        94.6
## Sigma[site_type:(Intercept),(Intercept)] 112230.4
## 
## Fit Diagnostics:
##            mean   sd    10%   50%   90%
## mean_PPD 158.8   10.9 144.9 158.8 172.9
## 
## The mean_ppd is the sample average posterior predictive distribution of the outcome variable (for details see help('summary.stanreg')).
## 
## MCMC diagnostics
##                                          mcse   Rhat   n_eff
## (Intercept)                                 1.6    1.0   597
## b[(Intercept) site_type:WB]                 1.6    1.0   709
## b[(Intercept) site_type:WC]                 1.6    1.0   617
## b[(Intercept) site_type:WK]                 1.5    1.0   899
## b[(Intercept) site_type:WNC]                1.6    1.0   669
## sigma                                       0.1    1.0  6680
## Sigma[site_type:(Intercept),(Intercept)] 1172.7    1.0  2590
## mean_PPD                                    0.1    1.0 13129
## log-posterior                               0.1    1.0   887
## 
## For each parameter, mcse is Monte Carlo standard error, n_eff is a crude measure of effective sample size, and Rhat is the potential scale reduction factor on split chains (at convergence Rhat=1).

Parameter	Effects	Component	Median	Mean	MAP	CI	CI_low	CI_high	pd	ps	Rhat	ESS
(Intercept)	fixed	conditional	7.663099e-01	12.12726	3.875211e-01	0.95	-30.094037	122.02677	0.5865714	0.1807857	1.003529	596.6051
b[(Intercept) site_type:WB]	random	conditional	2.111739e+02	202.37520	2.124639e+02	0.95	92.706889	263.08201	0.9959286	0.9944286	1.002704	709.4914
b[(Intercept) site_type:WC]	random	conditional	1.249533e+02	115.33694	1.279367e+02	0.95	4.978789	162.88107	0.9618571	0.9547143	1.003050	616.9764
b[(Intercept) site_type:WK]	random	conditional	6.060418e+01	53.82350	6.819934e+01	0.95	-51.921102	132.75646	0.9077143	0.8863571	1.001695	899.2485
b[(Intercept) site_type:WNC]	random	conditional	1.952327e+02	186.32204	1.977933e+02	0.95	75.950557	243.44464	0.9937857	0.9923571	1.002482	669.2157
Sigma[site_type:(Intercept),(Intercept)]	random	conditional	3.655611e+04	54338.82343	1.931836e+04	0.95	585.611229	159846.31186	1.0000000	1.0000000	1.001544	2590.1329
sigma	fixed	sigma	8.708692e+01	87.34466	8.634923e+01	0.95	76.780007	98.54709	1.0000000	1.0000000	NA	NA

Upland Data

## [[1]]

## 
## [[2]]

Upland Shrub Cover, Pooled Shrub Species - Combined Winter Range

Upland cover - site type + time_class
time_class	site_type	variable	mean	sd	min	med	max	n.valid	pct.valid
BL	core range	perc_cover	13.32	8.76	2.00	11.00	31.00	21	100
BL	non-core range	perc_cover	21.21	15.03	1.33	19.50	47.33	20	100
2013	core range	perc_cover	10.75	8.11	0.33	8.83	24.33	20	100
2013	non-core range	perc_cover	24.32	19.75	1.67	19.00	66.33	22	100
2018	core range	perc_cover	13.61	8.04	1.67	12.67	34.33	18	100
2018	non-core range	perc_cover	20.06	15.11	0.33	18.67	56.67	22	100

Summary by Range Type

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.allshr' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.25,5.32])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.50	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.46	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.00	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Upland Shrub Cover, Pooled Shrub Species - Core Winter Range

Upland cover - site type + time_class
time_class	site_type	variable	mean	sd	min	med	max	n.valid	pct.valid
BL	core range	perc_cover	13.32	8.76	2.00	11.00	31.00	21	100
BL	non-core range	perc_cover	21.21	15.03	1.33	19.50	47.33	20	100
2013	core range	perc_cover	10.75	8.11	0.33	8.83	24.33	20	100
2013	non-core range	perc_cover	24.32	19.75	1.67	19.00	66.33	22	100
2018	core range	perc_cover	13.61	8.04	1.67	12.67	34.33	18	100
2018	non-core range	perc_cover	20.06	15.11	0.33	18.67	56.67	22	100

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.allshr.uc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.25,5.32])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.50	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.46	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.00	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Upland Shrub Cover, Pooled Shrub Species - Noncore Winter Range

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.allshr.unc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.25,5.32])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.50	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.46	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.00	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Individual Species

shrub_species	n	percent
PUTR	79	0.26072607
RICE	52	0.17161716
CHNA	39	0.12871287
ARTR	35	0.11551155
ROWO	24	0.07920792
JUCO	19	0.06270627
PRVI	13	0.04290429
CEFE	8	0.02640264
SYRO	5	0.01650165
ACGL	3	0.00990099
RIXX	3	0.00990099
RUDE	3	0.00990099
DAFR	2	0.00660066
JUSC	2	0.00660066
MARE	2	0.00660066
PHMO	2	0.00660066
SYOR	2	0.00660066
XXXX	2	0.00660066
ALIN	1	0.00330033
AMAL	1	0.00330033
CEXX	1	0.00330033
JAAM	1	0.00330033
PRXX	1	0.00330033
PSME	1	0.00330033
SHCA	1	0.00330033
UNKN	1	0.00330033

Upland Shrub Cover, Purshia tridentata - Combined Winter Range

PUTR cover ~ site type + time_class
time_class	site_type	variable	mean	sd	min	med	max	n.valid	pct.valid
BL	core range	perc_cover	7.23	4.40	1.00	6.67	14.00	13	100
BL	non-core range	perc_cover	12.41	5.09	5.33	12.67	20.67	13	100
2013	core range	perc_cover	9.12	4.37	3.67	7.67	18.33	11	100
2013	non-core range	perc_cover	15.31	10.12	2.00	14.00	30.67	15	100
2018	core range	perc_cover	8.95	6.30	1.67	6.00	20.33	13	100
2018	non-core range	perc_cover	13.02	8.26	2.00	13.83	25.00	14	100

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.putr' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.29,5.24])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.29	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.37	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.00	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Upland Shrub Cover, Purshia tridentata - Core Winter Range

Upland cover - site type + time_class
time_class	site_type	variable	mean	sd	min	med	max	n.valid	pct.valid
BL	core range	perc_cover	7.23	4.40	1.00	6.67	14.00	13	100
BL	non-core range	perc_cover	12.41	5.09	5.33	12.67	20.67	13	100
2013	core range	perc_cover	9.12	4.37	3.67	7.67	18.33	11	100
2013	non-core range	perc_cover	15.31	10.12	2.00	14.00	30.67	15	100
2018	core range	perc_cover	8.95	6.30	1.67	6.00	20.33	13	100
2018	non-core range	perc_cover	13.02	8.26	2.00	13.83	25.00	14	100

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.putr.uc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.29,5.24])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.29	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.37	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.00	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Upland Shrub Cover, Purshia tridentata - Noncore Winter Range

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.putr.unc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.29,5.24])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.29	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.37	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.00	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Upland Shrub Cover, Ribes cereum - Combined Winter Range

Upland cover ~ site type + time_class
time_class	site_type	variable	mean	sd	min	med	max	n.valid	pct.valid
BL	core range	perc_cover	7.19	6.28	0.33	10.00	15.67	9	100
BL	non-core range	perc_cover	3.96	2.97	1.00	3.50	8.67	8	100
2013	core range	perc_cover	3.67	3.94	0.33	2.00	13.33	9	100
2013	non-core range	perc_cover	6.03	6.85	0.33	3.00	20.00	10	100
2018	core range	perc_cover	2.38	3.54	0.33	1.33	11.00	8	100
2018	non-core range	perc_cover	3.67	4.52	0.33	1.33	12.33	8	100

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.rice' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.14,5.36])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.27	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.17	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.91	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Upland Shrub Cover, Ribes cereum - Core Winter Range

Upland cover - site type + time_class
time_class	site_type	variable	mean	sd	min	med	max	n.valid	pct.valid
BL	core range	perc_cover	7.19	6.28	0.33	10.00	15.67	9	100
BL	non-core range	perc_cover	3.96	2.97	1.00	3.50	8.67	8	100
2013	core range	perc_cover	3.67	3.94	0.33	2.00	13.33	9	100
2013	non-core range	perc_cover	6.03	6.85	0.33	3.00	20.00	10	100
2018	core range	perc_cover	2.38	3.54	0.33	1.33	11.00	8	100
2018	non-core range	perc_cover	3.67	4.52	0.33	1.33	12.33	8	100

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.rice.uc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.14,5.36])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.27	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.17	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.91	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Upland Shrub Cover, Ribes cereum - Noncore Winter Range

Modeling and Posterior Description

Prior summary: cover ~ time_class + (1 | site_id)

## Priors for model 'stmod_cov.rice.unc' 
## ------
## Intercept (after predictors centered)
##  ~ normal(location = 0, scale = 2.5)
## 
## Coefficients
##   Specified prior:
##     ~ normal(location = [0,0], scale = [2.5,2.5])
##   Adjusted prior:
##     ~ normal(location = [0,0], scale = [5.04,5.31])
## 
## Covariance
##  ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
## ------
## See help('prior_summary.stanreg') for more details

Effect Description

Probability of Direction

Region of Practical Equivalence (ROPE)

Percent in ROPE
Parameter	CI	ROPE_low	ROPE_high	ROPE_Percentage	Effects	Component
(Intercept)	0.9	-0.1	0.1	0.00	fixed	conditional
time_class2013	0.9	-0.1	0.1	0.18	fixed	conditional
time_class2018	0.9	-0.1	0.1	0.16	fixed	conditional
(phi)	0.9	-0.1	0.1	0.00	fixed	conditional
mean_PPD	0.9	-0.1	0.1	0.44	fixed	conditional
log-posterior	0.9	-0.1	0.1	0.00	fixed	conditional

Contrasts

EMM summaries with type = “response,” the tests and confidence intervals are done before back-transforming. The ratios estimated are ratios of geometric means. A model with a log response is in fact a model for relative effects of any of its linear predictors, and this back-transformation to ratios goes hand-in-hand with that.

Purshia tridenta

Ribes cereum

Combined range

Core range

Noncore range

References

Gabry J., Mahr T. 2021. “bayesplot: Plotting for Bayesian Models.” R package version 1.8.1, <URL: https://mc-stan.org/bayesplot/>.

Gabray, J., Goodrich, B. 2020a. Estimating Generalized Linear Models for Continuous Data with rstanarm. ver:2020-07-20 <URL:https://mc-stan.org/rstanarm/articles/continuous.html/>.

Gabray, J., Goodrich, B. 2020b. Estimating Generalized Linear Models for Count Data with rstanarm. ver:2020-07-20 <URL:https://mc-stan.org/rstanarm/articles/count.html/>.

Gabray, J., Goodrich, B. 2020c. Modeling Rates/Proportions using Beta Regression with rstanarm. version 2020-07-20 <URL:https://mc-stan.org/rstanarm/articles/betareg.html/>.

Gabry J., Simpson D., Vehtari A., Betancourt M., Gelman A. 2019. “Visualization in Bayesian workflow.” J.R. Stat. Soc. A, 182, 389-402. doi: 10.1111/rssa.12378 (URL: https://doi.org/10.1111/rssa.12378).

Gelman, A., Hill, J., Vehtari, A. 2020. Regression and other stories. Cambridge University Press.

Goodrich B, Gabry J, Ali I & Brilleman S. 2020. rstanarm: Bayesian applied regression modeling via Stan. R package version 2.21.1 https://mc-stan.org/rstanarm

Lenth, R.V. 2021. emmeans: Estimated Marginal Means, aka Least-Squares Means. R package version 1.6.3. https://CRAN.R-project.org/package=emmeans

Brilleman, SL, MJ Crowther, M Moreno-Betancur, J Buros Novik, and R Wolfe. 2018. “StatCon 2018. 10-12 Jan 2018.” In. Pacific Grove, CA, USA. https://github.com/stan-dev/stancon_talks/.

Gabry, Jonah, Daniel Simpson, Aki Vehtari, Michael Betancourt, and Andrew Gelman. 2019. “Visualization in Bayesian Workflow.” Journal of the Royal Statistical Society: Series A (Statistics in Society) 182 (2): 389–402. https://doi.org/10.1111/rssa.12378.

Goodrich, B, J Gabry, I Ali, and S Brilleman. 2020. “Rstanarm: Bayesian Applied Regression Modeling via Stan. R Package Version 2.19.3 Https://Mc-Stan.org/Rstanarm.”

Kruschke, John K. 2018. “Rejecting or Accepting Parameter Values in Bayesian Estimation.” Advances in Methods and Practices in Psychological Science 1 (2): 270–80. https://doi.org/10.1177/2515245918771304.

Makowski, Dominique, Mattan Ben-Shachar, and Daniel Lüdecke. 2019. “bayestestR: Describing Effects and Their Uncertainty, Existence and Significance Within the Bayesian Framework.” Journal of Open Source Software 4 (40): 1541.

Wickham, Hadley, Mara Averick, Jennifer Bryan, Winston Chang, Lucy D’Agostino McGowan, Romain François, Garrett Grolemund, Alex Hayes, Lionel Henry, and Jim Hester. 2019. “Welcome to the Tidyverse.” Journal of Open Source Software 4 (43): 1686.

Zeigenfuss, Linda C., and Therese L. Johnson. 2015. “Monitoring of Vegetation Response to Elk Population and Habitat Management in Rocky Mountain National Park, 200814.” Reston, VA. https://doi.org/10.3133/ofr20151216.

Zeigenfuss, Linda C, Therese L Johnson, and Zachary Wiebe. 2011. “Monitoring Plan for Vegetation Responses to Elk Management in Rocky Mountain National Park.”