VGAM (version 0.9-4)

posbernoulli.tb: Positive Bernoulli Family Function with Time and Behavioural Effects

Description

Fits a GLM/GAM-like model to multiple Bernoulli responses where each row in the capture history matrix response has at least one success (capture). Sampling occasion effects and behavioural effects are accommodated.

Usage

posbernoulli.tb(link = "logit",
                parallel.t = FALSE ~ 1,
                parallel.b = FALSE ~ 0,
                drop.b     = FALSE ~ 1,
                type.fitted = c("likelihood.cond", "mean.uncond"),
                imethod = 1, iprob = NULL,
                p.small = 1e-4, no.warning = FALSE,
                ridge.constant = 0.01, ridge.power = -4)

Arguments

link, imethod, iprob

See CommonVGAMffArguments for information.

parallel.t, parallel.b, drop.b

A logical, or formula with a logical as the response. See CommonVGAMffArguments for information. The parallel.-type arguments specify whether the constraint matrices have a

type.fitted

Character, one of the choices for the type of fitted value returned. The default is the first one. Partial matching is okay. For "likelihood.cond": the probability defined by the conditional likelihood. For "mean.uncond"

ridge.constant, ridge.power

Determines the ridge parameters at each IRLS iteration. They are the constant and power (exponent) for the ridge adjustment for the working weight matrices (the capture probability block matrix, hence the first $\tau$ diagonal values). At iteratio

p.small, no.warning

See posbernoulli.t.

Value

• An object of class "vglmff" (see vglmff-class). The object is used by modelling functions such as vglm, and vgam.

Details

This model (commonly known as $M_{tb}$/$M_{tbh}$ in the capture--recapture literature) operates on a response matrix of 0s and 1s ($n \times \tau$). See posbernoulli.t for information that is in common. It allows time and behavioural effects to be modelled.

Evidently, the expected information matrix (EIM) seems not of full rank (especially in early iterations), so ridge.constant and ridge.power are used to try fix up the problem. The default link functions are $(logit \,p_{c1},\ldots,logit \,p_{c\tau},logit \,p_{r2},\ldots,logit \,p_{r\tau})^T$ where the subscript $c$ denotes capture, the subscript $r$ denotes recapture, and it is not possible to recapture the animal at sampling occasion 1. Thus $M = 2\tau - 1$. The parameters are currently prefixed by pcapture and precapture for the capture and recapture probabilities. This VGAM family function may be further modified in the future.

References

See posbernoulli.t.

See Also

posbernoulli.b (including N.hat), posbernoulli.t, posbinomial, Select, Huggins89table1, Huggins89.t1, deermice, prinia.

Examples

# Example 1: simulated data
nTimePts <- 5  # (aka tau == # of sampling occasions)
nnn <- 1000   # Number of animals
pdata <- rposbern(n = nnn, nTimePts = nTimePts, pvars = 2)
dim(pdata); head(pdata)

M_tbh.1 <- vglm(cbind(y1, y2, y3, y4, y5) ~ x2,
                posbernoulli.tb, data = pdata, trace = TRUE)
coef(M_tbh.1)  # First element is the behavioural effect
coef(M_tbh.1, matrix = TRUE)
constraints(M_tbh.1, matrix = TRUE)
summary(M_tbh.1, presid = FALSE)  # Standard errors are approximate
head(fitted(M_tbh.1))
head(model.matrix(M_tbh.1, type = "vlm"), 21)
dim(depvar(M_tbh.1))

M_tbh.2 <- vglm(cbind(y1, y2, y3, y4, y5) ~ x2,
                posbernoulli.tb(parallel.t = FALSE ~  0),
                data = pdata, trace = TRUE)
coef(M_tbh.2)  # First element is the behavioural effect
coef(M_tbh.2, matrix = TRUE)
constraints(M_tbh.2, matrix = TRUE)
summary(M_tbh.2, presid = FALSE)  # Standard errors are approximate
head(fitted(M_tbh.2))
head(model.matrix(M_tbh.2, type = "vlm"), 21)
dim(depvar(M_tbh.2))

# Example 2: deermice subset data
fit1 <- vglm(cbind(y1, y2, y3, y4, y5, y6) ~ sex + weight,
             posbernoulli.t, data = deermice, trace = TRUE)
coef(fit1)
coef(fit1, matrix = TRUE)
constraints(fit1, matrix = TRUE)
summary(fit1, presid = FALSE)  # Standard errors are approximate

# fit1 is the same as Fit1 (a M_{th} model):
Fit1 <- vglm(cbind(y1, y2, y3, y4, y5, y6) ~ sex + weight,
             posbernoulli.tb(drop.b = TRUE ~ sex + weight,
                             parallel.t = TRUE),  # No parallelism for the intercept
             data = deermice, trace = TRUE)
constraints(Fit1)