VGAM (version 0.7-8)

fill: Creates a Matrix of Appropriate Dimension

Description

A support function for the argument xij, it generates a matrix of an appropriate dimension.

Usage

fill(x, values = 0, ncolx = ncol(x))

Arguments

x

A vector or matrix which is used to determine the dimension of the answer, in particular, the number of rows. After converting x to a matrix if necessary, the answer is a matrix of values values, of dimension n

values

Numeric. The answer contains these values which are recycled if necessary.

ncolx

The number of columns of the returned matrix. The default is the number of columns of x.

Value

• matrix(values, nrow=nrow(x), ncol=ncolx), i.e., a matrix consisting of values values, with the number of rows matching x, and the default number of columns is the number of columns of x.

Warning

The use of the xij argument overrides other arguments such as exchangeable and zero. Care is needed in such cases. See the examples below.

Details

The xij argument for vglm allows the user to input variables specific to each linear predictor. For example, consider the bivariate logit model where the first/second linear/additive predictor is the logistic regression of the first/second binary response respectively. The third linear/additive predictor is log(OR) = eta3, where OR is the odds ratio. If one has ocular pressure as a covariate in this model then xij is required to handle the ocular pressure for each eye, since these will be different in general. [This contrasts with a variable such as age, the age of the person, which has a common value for both eyes.] In order to input these data into vglm one often finds that functions fill, fill1, etc. are useful.

All terms in the xij argument must appear in the main formula argument in vglm.

References

More information can be found at http://www.stat.auckland.ac.nz/~yee.

See Also

vglm, vglm.control.

Examples

fill(runif(5))
fill(runif(5), ncol=3)
fill(runif(5), val=1, ncol=3)

# Generate eyes data for the examples below. Eyes are independent (OR=1).
set.seed(123)
n = 2000  # Number of people
eyes = data.frame(lop = round(runif(n), 2),
                  rop = round(runif(n), 2),
                  age = round(rnorm(n, 40, 10)))
eyes = transform(eyes,
    mop = (lop + rop) / 2, # mean ocular pressure
    eta1 = 0 - 2*lop + 0.04*age, # Linear predictor for left eye
    eta2 = 0 - 2*rop + 0.04*age) # Linear predictor for right eye
eyes = transform(eyes,
    leye = rbinom(n, size=1, prob=exp(eta1)/(1+exp(eta1))),
    reye = rbinom(n, size=1, prob=exp(eta2)/(1+exp(eta2))))

# Example 1
# Non-exchangeable errors (misspecified model)
fit1 = vglm(cbind(leye,reye) ~ lop + rop + fill(lop) + age,
            family = binom2.or(exchangeable=FALSE, zero=NULL),
            xij = op ~ lop + rop + fill(lop), data=eyes)
model.matrix(fit1, type="lm")[1:7,]   # LM model matrix
model.matrix(fit1, type="vlm")[1:7,]  # Big VLM model matrix
coef(fit1)
coef(fit1, matrix=TRUE)  # Looks wrong but is correct
coef(fit1, matrix=TRUE, compress=FALSE)  # Looks wrong but is correct
constraints(fit1)
max(abs(predict(fit1)-predict(fit1, new=eyes))) # Predicts correctly
summary(fit1)


# Example 2
# Nonexchangeable errors (misspecified model), OR is a function of mop
fit2 = vglm(cbind(leye,reye) ~ lop + rop + mop + age,
            family = binom2.or(exchangeable=FALSE, zero=NULL),
            xij = op ~ lop + rop + mop, data=eyes)
model.matrix(fit2, type="lm")[1:7,]   # LM model matrix
model.matrix(fit2, type="vlm")[1:7,]  # Big VLM model matrix
coef(fit2)
coef(fit2, matrix=TRUE)  # correct
coef(fit2, matrix=TRUE, compress=FALSE)  # correct
max(abs(predict(fit2)-predict(fit2, new=eyes))) # Predicts correctly
summary(fit2)


# Example 3. This model is correctly specified.
# Exchangeable errors
fit3 = vglm(cbind(leye,reye) ~ lop + rop + fill(lop) + age,
            family = binom2.or(exchangeable=TRUE, zero=3),
            xij = op ~ lop + rop + fill(lop), data=eyes)
model.matrix(fit3, type="lm")[1:7,]   # LM model matrix
model.matrix(fit3, type="vlm")[1:7,]  # Big VLM model matrix
coef(fit3)
coef(fit3, matrix=TRUE) # Looks wrong but is correct
coef(fit3, matrix=TRUE, compress=FALSE) # Looks wrong but is correct
predict(fit3, new=eyes[1:4,])  # Note the 'scalar' OR, i.e., zero=3
max(abs(predict(fit3)-predict(fit3, new=eyes))) # Predicts correctly
summary(fit3)


# Example 4. This model uses regression splines on ocular pressure.
# It assumes exchangeable errors.
BS = function(x, ...) bs(c(x,...), df=3)[1:length(x),]
fit4 = vglm(cbind(leye,reye) ~ BS(lop,rop,mop) + BS(rop,lop,mop) +
            fill(BS(lop,rop,mop)) + age,
            family = binom2.or(exchangeable=TRUE, zero=3),
            xij = BS(op) ~ BS(lop,rop,mop) + BS(rop,lop,mop) +
                  fill(BS(lop,rop,mop)), data=eyes)
model.matrix(fit4, type="lm")[1:7,]   # LM model matrix
model.matrix(fit4, type="vlm")[1:7,]  # Big VLM model matrix
coef(fit4)
coef(fit4, matrix=TRUE) # Looks wrong but is correct
coef(fit4, matrix=TRUE, compress=FALSE) # Looks wrong but is correct
predict(fit4, new=eyes[1:4,])  # Note the 'scalar' OR, i.e., zero=3
max(abs(predict(fit4)-predict(fit4, new=eyes))) # Predicts correctly
summary(fit4)