https://github.com/cran/spatstat
Tip revision: 6cceeeb7d80470e77e1363c4bcaa85ae6ccaa844 authored by Adrian Baddeley on 19 May 2011, 10:05:10 UTC
version 1.22-1
version 1.22-1
Tip revision: 6cceeeb
residppm.R
#
# residppm.R
#
# computes residuals for fitted point process model
#
#
# $Revision: 1.14 $ $Date: 2011/02/15 07:24:08 $
#
residuals.ppm <- function(object, type="raw", ..., check=TRUE, drop=FALSE,
fittedvalues = fitted.ppm(object, check=check, drop=drop)) {
verifyclass(object, "ppm")
if(check && missing(fittedvalues) && damaged.ppm(object))
stop("object format corrupted; try update(object, use.internal=TRUE)")
type <- pickoption("type", type,
c(inverse="inverse",
raw="raw",
pearson="pearson",
Pearson="pearson",
score="score"))
typenames <- c(inverse="inverse-lambda residuals",
raw="raw residuals",
pearson="Pearson residuals",
score="score residuals")
typename <- typenames[[type]]
# Extract quadrature points and weights
Q <- quad.ppm(object, drop=drop)
U <- union.quad(Q) # quadrature points
Z <- is.data(Q) # indicator data/dummy
# W <- w.quad(Q) # quadrature weights
# Compute fitted conditional intensity at quadrature points
lambda <- fittedvalues
# indicator is 1 if lambda > 0
# (adjusted for numerical behaviour of predict.glm)
indicator <- (lambda > .Machine$double.eps)
if(type == "score") {
# need the covariates
X <- model.matrix(object)
if(drop) {
gs <- getglmsubset(object)
ok <- !is.na(gs) && gs
X <- X[ok,]
}
}
# Evaluate residual measure components
discrete <- switch(type,
raw = rep(1, sum(Z)),
inverse = 1/lambda[Z],
pearson = 1/sqrt(lambda[Z]),
score = X[Z, ]
)
density <- switch(type,
raw = -lambda,
inverse = -indicator,
pearson = -indicator * sqrt(lambda),
score = -lambda * X)
# Residual measure (return value)
res <- msr(Q, discrete, density)
# name the residuals
attr(res, "type") <- type
attr(res, "typename") <- typename
return(res)
}