residuals.vlmc: Compute Residuals of a Fitted VLMC Object

Description

Compute residuals of a fitted vlmc object.

This is yet a matter of research and may change in the future.

Usage

# S3 method for vlmc
residuals(object,
        type = c("classwise",
                 "deviance", "pearson", "working", "response", "partial"),
        y = object$y, ...)

Value

If type = "classwise" (the default), a numeric matrix of dimension

\(n \times m\) of values \(I_{i,j} - p_{i,j}\) where the indicator \(I_{i,j}\) is 1 iff

y[i] == a[j] and a is the alphabet (or levels) of

y, and \(p_{i,j}\) are the elements of the estimated (1-step ahead) predicted probabilities, p <- predict(object). Hence, for each \(i\), the only positive residual stands for the observed class.

For all other types, the result is a numeric vector of the length of the original time-series (with first element NA).

For type = "deviance",

\(r_i = \pm\sqrt{-2\log(P_i)}\)

where \(P_i\) is the predicted probability for the i-th observation which is the same as \(p_{i,y_i}\) above (now assuming \(y_i \in \{1,2,\dots,m\)). The sum of the squared deviance residuals is the deviance of the fitted model.

Arguments

object: typically the result of vlmc(..).
type: The type of residuals to compute, defaults to "classwise" which returns an \(n \times m\) matrix, see below. The other types only make sense when the discrete values of y are ordered which always includes the binary case (\(m=2\)).
The "deviance" residuals \(r\) are defined similarly as for logistic regression, see below.
"pearson", "working" and "response" are currently identical and give the difference of the underlying integer code (of the discrete data). Note that "partial" residuals are not yet defined!
y: discrete time series with respect to which the residuals are to be computed.
...: possibly further arguments (none at the moment).

Author

Martin Maechler

Examples

Run this code

example(vlmc)
rp <- residuals(vlmc.pres)
stopifnot(all(abs(apply(rp[-1,],1,sum)) < 1e-15))
matplot(seq(presidents), rp, ylab = "residuals", type="l")
## ``Tukey-Anscombe'' (the following is first stab at plot method):
matplot(fitted(vlmc.pres), rp, ylab = "residuals", xaxt = "n",
        type="b", pch=vlmc.pres$alpha)
axis(1, at = 0:(vlmc.pres$alpha.len-1),
     labels = strsplit(vlmc.pres$alpha,"")[[1]])


summary(rd <- residuals(vlmc.pres, type = "dev"))
rd[1:7]
## sum of squared dev.residuals === deviance :
all.equal(sum(rd[-1] ^ 2),
          deviance(vlmc.pres))
stopifnot(all.equal(sum(rd[-1] ^ 2),
                    deviance(vlmc.pres),tol=1e-12))