modelling total diversity with random effects

Source: R/model_betta_random.R

This function extends betta() to permit random effects modelling.

betta_random(
  chats = NULL,
  ses,
  X = NULL,
  groups = NULL,
  formula = NULL,
  data = NULL
)

Arguments

chats

A vector of estimates of total diversity at different sampling locations. Required with the groups argument, and optionally with the X argument.
ses

The standard errors in chats, the diversity estimates. This can either be a vector of standard errors (with the arguments chats and X), or the name of the variable in the dataframe data that contains the standard errors (with the arguments formula and data).
X

A numeric matrix of covariates corresponding to fixed effects. If not supplied, an intercept-only model will be fit. Optional with the chats and groups arguments.
groups

A categorical variable representing the random-effects groups that each of the estimates belong to. Required with the chats argument and optionally with the X argument.
formula

A formula object of the form \(y ~ x | group\). Required with the data argument.
data

A dataframe containing the response, response standard errors, covariates, and grouping variable. Required with the formula argument.

Value

table

A coefficient table for the model parameters. The columns give the parameter estimates, standard errors, and p-values, respectively. This model is only as effective as your diversity estimation procedure; for this reason please confirm that your estimates are appropriate and that your model is not misspecified. betta_pic may be useful for this purpose.

cov

Estimated covariance matrix of the parameter estimates.

ssq_u

The estimate of the heterogeneity variance.

ssq_g

Estimates of within-group variance. The estimate will be zero for groups with only one observation.

homogeneity

The test statistic and p-value for the test of homogeneity.

global

The test statistic and p-value for the test of model explanatory power.

blups

The conditional expected values of the diversity estimates (conditional on the random effects). Estimates of variability for the random effects case are unavailable at this time; please contact the maintainer if needed.

function.args

A list containing values initially passed to betta_random.

References

Willis, A., Bunge, J., and Whitman, T. (2015). Inference for changes in biodiversity. arXiv preprint.

See also

betta;

Author

Amy Willis

Examples


df <- data.frame(chats = c(2000, 3000, 4000, 3000),
ses = c(100, 200, 150, 180),
Cont_var = c(100, 150, 100, 50),
groups = c("a", "a", "a", "b"))

# formula notation
betta_random(formula = chats ~ Cont_var| groups,
ses = ses,
 data = df)
#> $table
#>                 Estimates Standard Errors p-values
#> (Intercept) 3000.32999173       384.97695    0.000
#> Cont_var      -0.02907528         4.35074    0.995
#> 
#> $cov
#>             (Intercept)   Cont_var
#> (Intercept)  860819.418 -8851.3807
#> Cont_var      -8851.381   109.9433
#> 
#> $ssq_u
#>        
#> 241802 
#> 
#> $ssq_group
#>        a        b 
#> 701816.2      0.0 
#> 
#> $homogeneity
#> [1] 144.1595   0.0000
#> 
#> $global
#> [1] 6.064440e+01 6.883383e-15
#> 
#> $blups
#> [1] 2744.513 2996.960 3248.350 2999.867
#> 
#> $function.args
#> $function.args$chats
#>    1    2    3    4 
#> 2000 3000 4000 3000 
#> 
#> $function.args$ses
#> [1] 100 200 150 180
#> 
#> $function.args$X
#>   (Intercept) Cont_var
#> 1           1      100
#> 2           1      150
#> 3           1      100
#> 4           1       50
#> attr(,"assign")
#> [1] 0 1
#> 
#> $function.args$groups
#> [1] "a" "a" "a" "b"
#> 
#> $function.args$formula
#> chats ~ Cont_var | groups
#> <environment: 0x7f87ad347350>
#> 
#> $function.args$data
#>   chats ses Cont_var groups
#> 1  2000 100      100      a
#> 2  3000 200      150      a
#> 3  4000 150      100      a
#> 4  3000 180       50      b
#> 
#> $function.args$model_type
#> [1] "mixed"
#> 
#> 

# direct input
betta_random(c(2000, 3000, 4000, 3000), c(100, 200, 150, 180),
X = cbind(Int = 1, Cont_var = c(100, 150, 100, 50)),
groups = c("a", "a", "a", "b"))
#> $table
#>              Estimates Standard Errors p-values
#> Int      3000.32999173       384.97695    0.000
#> Cont_var   -0.02907528         4.35074    0.995
#> 
#> $cov
#>                 Int   Cont_var
#> Int      860819.418 -8851.3807
#> Cont_var  -8851.381   109.9433
#> 
#> $ssq_u
#>        
#> 241802 
#> 
#> $ssq_group
#>        a        b 
#> 701816.2      0.0 
#> 
#> $homogeneity
#> [1] 144.1595   0.0000
#> 
#> $global
#> [1] 6.064440e+01 6.883383e-15
#> 
#> $blups
#> [1] 2744.513 2996.960 3248.350 2999.867
#> 
#> $function.args
#> $function.args$chats
#> [1] 2000 3000 4000 3000
#> 
#> $function.args$ses
#> [1] 100 200 150 180
#> 
#> $function.args$X
#>      Int Cont_var
#> [1,]   1      100
#> [2,]   1      150
#> [3,]   1      100
#> [4,]   1       50
#> 
#> $function.args$groups
#> [1] "a" "a" "a" "b"
#> 
#> $function.args$formula
#> NULL
#> 
#> $function.args$data
#> NULL
#> 
#> $function.args$model_type
#> [1] "mixed"
#> 
#> 

## handles missing data
betta_random(c(2000, 3000, 4000, 3000), c(100, 200, 150, NA),
groups = c("a", NA,
    "b", "b"))
#> $table
#>      Estimates Standard Errors p-values
#> [1,]      3000        1004.049    0.003
#> 
#> $cov
#>         [,1]
#> [1,] 1008115
#> 
#> $ssq_u
#>       
#> 2e+06 
#> 
#> $ssq_group
#> a b 
#> 0 0 
#> 
#> $homogeneity
#> [1] 144.4444   0.0000
#> 
#> $global
#> [1] 8.92755 0.00000
#> 
#> $blups
#> [1] 2004.975       NA 3988.875       NA
#> 
#> $function.args
#> $function.args$chats
#> [1] 2000 3000 4000 3000
#> 
#> $function.args$ses
#> [1] 100 200 150  NA
#> 
#> $function.args$X
#>      [,1]
#> [1,]    1
#> [2,]    1
#> [3,]    1
#> [4,]    1
#> 
#> $function.args$groups
#> [1] "a" NA  "b" "b"
#> 
#> $function.args$formula
#> NULL
#> 
#> $function.args$data
#> NULL
#> 
#> $function.args$model_type
#> [1] "mixed"
#> 
#>