Retrieve summary statistics for seasonal metrics and trend results

Usage

anlz_sumstats(
  mod,
  metfun = mean,
  doystr = 1,
  doyend = 364,
  yrstr = 2000,
  yrend = 2019,
  yromit = NULL,
  nsim = 10000,
  confint = 0.95,
  useave = FALSE,
  ...
)

Arguments

mod: input model object as returned by anlz_gam
metfun: function input for metric to calculate, e.g., mean, var, max, etc
doystr: numeric indicating start Julian day for extracting averages
doyend: numeric indicating ending Julian day for extracting averages
yrstr: numeric for starting year for trend model, see details
yrend: numeric for ending year for trend model, see details
yromit: optional numeric vector for years to omit from the plot, see details
nsim: numeric indicating number of random draws for simulating uncertainty
confint: numeric from zero to one indicating confidence interval level for summarizing the mixed-effects meta-analysis model, see details
useave: logical indicating if anlz_avgseason is used for the seasonal metric calculation, see details
...: additional arguments passed to metfun, e.g., na.rm = TRUE

Value

A list object with named elements:

mixmet: mixmeta object of the fitted mixed-effects meta-analysis trend model
metseason: tibble object of the fitted seasonal metrics as returned by anlz_metseason or anlz_avgseason
summary: summary of the mixmet object
coeffs: tibble object of the slope estimate coefficients from the mixmet model. An approximately linear slope estimate will be included as slope.approx if trans = 'log10' for the GAM used in mod.

Details

This function is primarily for convenience to return summary statistics of a fitted GAM from anlz_gam.

Note that confint only applies to the summary and coeffs list outputs. It does not apply to the metseason list element output that is default set to 95

Set useave = T to speed up calculations if metfun = mean. This will use anlz_avgseason to estimate the seasonal summary metrics using a non-stochastic equation.

Examples

library(dplyr)

# data to model
tomod <- rawdat %>%
  filter(station %in% 34) %>%
  filter(param %in% 'chl') %>% 
  filter(yr > 2015)

mod <- anlz_gam(tomod, trans = 'log10')

anlz_sumstats(mod, metfun = mean, doystr = 90, doyend = 180, yrstr = 2016, 
  yrend = 2019, nsim = 100)
#> $mixmet
#> Call:  mixmeta::mixmeta(formula = met ~ yr, S = S, data = totrnd, random = ~1 | 
#>     yr, method = "reml")
#> 
#> Fixed-effects coefficients:
#> (Intercept)           yr  
#>    -74.8387       0.0375  
#> 
#> 4 units, 1 outcome, 4 observations, 2 fixed and 1 random-effects parameters
#> logLik     AIC     BIC  
#> 1.0379  3.9243  0.0037  
#> 
#> 
#> $metseason
#> # A tibble: 4 × 7
#>      yr   met     se bt_lwr bt_upr bt_met dispersion
#>   <dbl> <dbl>  <dbl>  <dbl>  <dbl>  <dbl>      <dbl>
#> 1  2016 0.766 0.0856   4.53   9.81   6.67     0.0504
#> 2  2017 0.772 0.0846   4.62   9.92   6.77     0.0504
#> 3  2018 1.03  0.0624   9.19  16.1   12.2      0.0504
#> 4  2019 0.808 0.0694   5.37  10.0    7.34     0.0504
#> 
#> $summary
#> Call:  mixmeta::mixmeta(formula = met ~ yr, S = S, data = totrnd, random = ~1 | 
#>     yr, method = "reml")
#> 
#> Univariate extended random-effects meta-regression
#> Dimension: 1
#> Estimation method: REML
#> 
#> Fixed-effects coefficients
#>              Estimate  Std. Error        z  Pr(>|z|)   95%ci.lb  95%ci.ub   
#> (Intercept)  -74.8387    134.7875  -0.5552    0.5787  -339.0172  189.3399   
#> yr             0.0375      0.0668   0.5615    0.5744    -0.0934    0.1685   
#> ---
#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
#> 
#> Random-effects (co)variance components
#>  Formula: ~1 | yr
#>  Structure: General positive-definite
#>   Std. Dev
#>     0.1277
#> 
#> Univariate Cochran Q-test for residual heterogeneity:
#> Q = 8.9722 (df = 2), p-value = 0.0113
#> I-square statistic = 77.7%
#> 
#> 4 units, 1 outcome, 4 observations, 2 fixed and 1 random-effects parameters
#> logLik     AIC     BIC  
#> 1.0379  3.9243  0.0037  
#> 
#> 
#> $coeffs
#> # A tibble: 1 × 8
#>   slope.approx  slope slope.se     z     p likelihood   ci.lb ci.ub
#>          <dbl>  <dbl>    <dbl> <dbl> <dbl>      <dbl>   <dbl> <dbl>
#> 1        0.696 0.0375   0.0668 0.562 0.574      0.713 -0.0934 0.168
#>