residuals.maxlogL | R Documentation |
maxlogL
model.residuals.maxlogL
is the maxlogLreg
specific method for the
generic function residuals which extracts the residuals from a fitted model.
## S3 method for class 'maxlogL'
residuals(object, parameter = NULL, type = "rqres", routine, ...)
object |
an object of |
parameter |
a character which specifies residuals for a specific parameter. |
type |
a character with the type of residuals to be computed.
The default value is |
routine |
a character specifying the integration routine.
|
... |
further arguments for the integration routine. |
For type = "deviance"
, the residuals are computed using the following
expression
r^D_i = \mbox{sign}(y_i - \hat{\mu}_i) d_i^{1/2},
where d_i
is the residual deviance of each data point. In this context,
\hat{\mu}
is the estimated mean, computed as the expected value using
the estimated parameters of a fitted maxlogLreg
model.
On the other hand, for type = "response"
the computation is simpler
r^R_i = (y_i - \hat{\mu}_i).
a vector with the specified residuals of a maxlogLreg
model.
Jaime Mosquera GutiƩrrez, jmosquerag@unal.edu.co
library(EstimationTools)
#----------------------------------------------------------------------------
# Example 1: Test deviance residuals
set.seed(123)
n <- 500
x <- runif(n = n, min = 0, max = 1)
y <- rweibull(n = n, shape = 1, scale = exp(4.5 + 0.5*x))
status <- rep(1, n) # sample(0:1, size = n, replace = TRUE)
distribution <- Surv(y, status) ~ dweibull
formulas <- list(
scale.fo = ~ x
)
fixed <- list(shape = 1)
links <- list(
over = "scale",
fun = "log_link"
)
model <- maxlogLreg(
formulas = formulas,
y_dist = distribution,
fixed = fixed,
link = links
)
# Using `residuals` method
cox_snell_residuals_test <- residuals(model, type = "cox-snell")
martingale_residuals_test <- residuals(model, type = "martingale")
deviance_residuals_test <- residuals(model, type = "right-censored-deviance")
# From scratch
cox_snell_residuals_ref <- -pweibull(
q = y,
shape = 1,
scale = exp(cbind(rep(1, n), x) %*% cbind(coef(model))),
lower.tail = FALSE,
log.p = TRUE
)
martingale_residuals_ref <- status - cox_snell_residuals_ref
deviance_residuals_ref <- sign(martingale_residuals_ref) * (
-2 * (martingale_residuals_ref + status*log(status - martingale_residuals_ref))
)^ 0.5
plot(cox_snell_residuals_test, cox_snell_residuals_ref)
plot(martingale_residuals_test, martingale_residuals_ref)
plot(deviance_residuals_test, deviance_residuals_ref)
#----------------------------------------------------------------------------
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.