anovir (version 0.1.0)

nll_frailty_correlated: Negative log-likelihood function: correlated frailty model

Description

Function calculating the negative log-likelihood (nll) for patterns of mortality in infected and uninfected treatments when there is unobserved variation in background mortality and mortality due to infection, where these two sources of variation are positively correlated.

Usage

nll_frailty_correlated(
  a1 = a1,
  b1 = b1,
  a2 = a2,
  b2 = b2,
  theta01 = theta01,
  theta02 = theta02,
  rho = rho,
  data = data,
  time = time,
  censor = censor,
  infected_treatment = infected_treatment,
  d1 = "Weibull",
  d2 = "Weibull"
)

Arguments

a1, b1

location and scale parameters describing background mortality

a2, b2

location and scale parameters describing mortality due to infection

theta01

parameter describing variance of unobserved variation in background rate of mortality

theta02

parameter describing variance of unobserved variation in rate of mortality due to infection

rho

parameter for strength of positive correlation between theta01 and theta02; rho > 0

data

name of data frame containing survival data

time

name of data frame column identifying time of event; time > 0

censor

name of data frame column idenifying if event was death (0) or right-censoring (1)

infected_treatment

name of data frame column identifying if data are from an infected (1) or uninfected (0) treatment

d1, d2

names of probability distributions chosen to describe background mortality and mortality due to infection, respectively; both default to the Weibull distribution

Details

This function assumes both sources of unobserved variation follow the Gamma distribution where both have a mean = 1.0, and variances 'theta01' and 'theta02', respectively. The nll function is calculated based on seven parameters; location and scale parameters for background mortality and mortality due to infection, respectively, the unobserved variance of in each source of mortality, and the strength of the positive correlation between them.

See Also

nll_frailty

Examples

# NOT RUN {
# step #1: parameterise nll function to be passed to 'mle2'
    m01_prep_function <- function(
      a1 = a1, b1 = b1, a2 = a2, b2 = b2,
      theta01 = theta01, theta02 = theta02, rho = rho){
        nll_frailty_correlated(
          a1 = a1, b1 = b1, a2 = a2, b2 = b2,
          theta01 = theta01, theta02 = theta02, rho = rho,
          data = data_lorenz,
          time = t,
          censor = censored,
          infected_treatment = g,
          d1 = "Gumbel",
          d2 = "Gumbel"
          )}

# step #2: send 'prep_function' to 'mle2' for maximum likelihood estimation
#          lower bounds of estimates set to 1e-6
  m01 <- mle2(m01_prep_function,
    start = list(a1 = 20, b1 = 5, a2 = 20, b2 = 4,
                 theta01 = 1, theta02 = 1, rho = 1),
    method = "L-BFGS-B",
    lower = list(a1 = 1e-6, b1 = 1e-6, a2 = 1e-6, b2 = 1e-6,
                 theta01 = 1e-6, theta02 = 1e-6, rho = 1e-6)
                 )

  summary(m01)

# NB no std. errors estimated and estimate of 'theta01' at lower boundary

# rerun model with theta01 set at lower boundary
    m02 <- mle2(m01_prep_function,
      start = list(a1 = 20, b1 = 5, a2 = 20, b2 = 4,
                   theta01 = 1, theta02 = 1, rho = 1),
      fixed = list(theta01 = 1e-6),
      method = "L-BFGS-B",
      lower = list(a1 = 1e-6, b1 = 1e-6, a2 = 1e-6,
                   b2 = 1e-6, theta02 = 1e-6, rho = 1e-6)
      )

summary(m02)

# NB std. error of 'rho' crosses lower boundary

# rerun model with theta01 and rho set at lower limits
    m03 <- mle2(m01_prep_function,
      start = list(a1 = 20, b1 = 5, a2 = 20, b2 = 4,
                   theta01 = 1, theta02 = 1, rho = 1),
      fixed = list(theta01 = 1e-6, rho = 1e-6),
      method = "L-BFGS-B",
      lower = list(a1 = 1e-6, b1 = 1e-6, a2 = 1e-6,
                   b2 = 1e-6, theta02 = 1e-6)
      )

    summary(m03)

# result of m03 corresponds with estimates from 'nll_frailty' model,
  # i.e., where it is assumed there is no unobserved variation in the
  # rate of background mortality and where the gamma distribution describes
  # the unobserved variation in virulence

# step #1: parameterise nll function to be passed to 'mle2'
    m04_prep_function <- function(a1, b1, a2, b2, theta){
      nll_frailty(a1, b1, a2, b2, theta,
        data = data_lorenz,
        time = t,
        censor = censored,
        infected_treatment = g,
        d1 = "Gumbel",
        d2 = "Gumbel",
        d3 = "Gamma"
        )}

# step #2: send 'prep_function' to 'mle2' for maximum likelihood estimation
    m04 <- mle2(m04_prep_function,
      start = list(a1 = 20, b1 = 5, a2 = 20, b2 = 5, theta = 2)
      )

    summary(m04)

# compare estimated values m03 vs. m04
    coef(m03) ; coef(m04)


# }