theta_CI: Compute Confidence Interval for Theta

In DImodels: Diversity-Interactions (DI) Models

Compute Confidence Interval for Theta

Description

This function allows the computation of a confidence interval for theta from a model object created from DI that includes the argument estimate_theta = TRUE, or from certain model objects created from autoDI.

A description of the non-linear parameter theta is available in Connolly et al 2013.

Usage

theta_CI(obj, conf = .95, n =100)

Arguments

obj	DI model object.
conf	Confidence level of the interval. The default is 0.95.
n	Number of subintervals in which the root is sought. The default is 100.

Details

The confidence interval calculated here is based on the values obtained when profiling the log-likelihood function for different values of theta. It is obtained in four steps:

1. define a grid of values for theta ranging from 0.01 to 2.5 of length 101 including the profile likelihood estimate of theta

2. fit the DI model setting theta equal to each value in the grid and obtain the log-likelihood value corresponding to each value of theta

3. obtain linear interpolations between the log-likelihood (l) values (here we use approxfun)

4. calculate the lower and upper values of the CI by obtaining the values of theta corresponding to a log-likelihood value of max_\theta(l) - 0.5*\chi^2_(1-\alpha;1), where max_\theta(l) is the maximum value of the profile log-likelihood obtained in the grid and \chi^2_(conf;1) is the conf*100% percentile of the chi-squared distribution with 1 d.f.

When fitting any DI model setting estimate_theta = TRUE, steps 1 and 2 are automatically done within the DI function call. The theta_CI function performs steps 3 and 4 above to return the CI.

Note that when maximising the log-likelihood to find the estimate for theta, the parametric space is limited between 0.01 and 1.5. The larger grid (up to 2.5) is constructed to allow for obtaining the upper bound of the confidence interval in case the estimate of theta is close to 1.5.

Value

The function returns a named numeric vector with two values: the lower and upper limits of the conf*100% CI for theta.

Author(s)

Rafael A. Moral, John Connolly and Caroline Brophy

References

Connolly J, T Bell, T Bolger, C Brophy, T Carnus, JA Finn, L Kirwan, F Isbell, J Levine, A Lüscher, V Picasso, C Roscher, MT Sebastia, M Suter and A Weigelt (2013) An improved model to predict the effects of changing biodiversity levels on ecosystem function. Journal of Ecology, 101, 344-355.

See Also

DI autoDI

Other examples using the theta_CI function:
The Bell dataset examples.
The sim2 dataset examples.
The sim5 dataset examples.

Examples

## Load the sim5 data
  data(sim5)
## View the first five entries
  head(sim5)
## Explore the variables in sim5
  str(sim5)
  
## Fit the functional group model, with theta, using DI and the FG tag
  m1 <- DI(y = "response", prop = 3:11, 
           FG = c("FG1","FG1","FG1","FG1","FG1","FG2","FG2","FG3","FG3"), DImodel = "FG", 
           estimate_theta = TRUE, data = sim5)
  summary(m1)
  CI_95 <- theta_CI(m1, conf = .95)
  CI_95
## Graph the profile likelihood
  library(ggplot2)
  ggplot(m1$profile_loglik, aes(x = grid, y = prof)) +
    theme_bw() +
    geom_line() +
    xlim(0,1.5) +
    xlab(expression(theta)) +
    ylab("Log-likelihood") + 
    geom_vline(xintercept = CI_95, lty = 3) + 
    labs(title = "   Log-likelihood versus theta", 
         caption = "dotted vertical lines are upper and lower bounds of 95% CI for theta")

DImodels documentation built on May 29, 2024, 7:05 a.m.