R: Calculates GoF statistics based on cumulative residual...
cumres.glm
R Documentation
Calculates GoF statistics based on cumulative residual processes
Description
Given the generalized linear models model
g(E(Y_i|X_{i1},...,X_{ik})) = ∑_{i=1}^k
β_jX_{ij}
the cumres-function calculates the
the observed cumulative sum of residual process,
cumulating the residuals, e_i, by the jth
covariate:
W_j(t) = n^{-1/2}∑_{i=1}^n
1_{{X_{ij}<t}}e_i
and Kolmogorov-Smirnov and
Cramer-von-Mises test statistics are calculated via
simulation from the asymptotic distribution of the
cumulative residual process under the null (Lin et al.,
2002).
Usage
## S3 method for class 'glm'
cumres(model,
variable = c("predicted", colnames(model.matrix(model))),
data = data.frame(model.matrix(model)), R = 1000,
b = 0, plots = min(R, 50), breakties = 1e-12,
seed = round(runif(1, 1, 1e+09)), ...)
Arguments
model
Model object (lm or glm)
variable
List of variable to order the residuals
after
data
data.frame used to fit model (complete
cases)
R
Number of samples used in simulation
b
Moving average bandwidth (0 corresponds to
infinity = standard cumulated residuals)
plots
Number of realizations to save for use in
the plot-routine
breakties
Add unif[0,breakties] to observations
seed
Random seed
...
additional arguments
Value
Returns an object of class 'cumres'.
Note
Currently linear (normal), logistic and poisson
regression models with canonical links are supported.
Author(s)
Klaus K. Holst
References
D.Y. Lin and L.J. Wei and Z. Ying (2002)
Model-Checking Techniques Based on Cumulative
Residuals. Biometrics, Volume 58, pp 1-12.
John Q. Su and L.J. Wei (1991) A lack-of-fit test
for the mean function in a generalized linear model.
Journal. Amer. Statist. Assoc., Volume 86, Number 414, pp
420-426.
See Also
cox.aalen in the
timereg-package for similar GoF-methods for
survival-data.
Examples
sim1 <- function(n=100, f=function(x1,x2) {10+x1+x2^2}, sd=1, seed=1) {
if (!is.null(seed))
set.seed(seed)
x1 <- rnorm(n);
x2 <- rnorm(n)
X <- cbind(1,x1,x2)
y <- f(x1,x2) + rnorm(n,sd=sd)
d <- data.frame(y,x1,x2)
return(d)
}
d <- sim1(100); l <- lm(y ~ x1 + x2,d)
system.time(g <- cumres(l, R=100, plots=50))
g
plot(g)
g1 <- cumres(l, c("y"), R=100, plots=50)
g1
g2 <- cumres(l, c("y"), R=100, plots=50, b=0.5)
g2