-
Notifications
You must be signed in to change notification settings - Fork 5
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Dev #26
Dev #26
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Small comments on residuals.nonprobsvy
@@ -346,19 +357,32 @@ cooks.distance.nonprobsvy <- function(model, | |||
hatvalues.nonprobsvy <- function(model, | |||
...) { # TODO reduce execution time and glm.fit object and customise to variable selection | |||
if (any(c("nonprobsvy_dr", "nonprobsvy_ipw") %in% class(model))) { | |||
X <- model$X | |||
propensity_scores <- model$prob | |||
W <- Matrix::Diagonal(x = propensity_scores * (1 - propensity_scores)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Preety sure you can just leave W
as a vector and this will work faster
} | ||
#hats <- Matrix::Diagonal(x = W %*% object$X %*% XWX_inv %*% t(object$X)) | ||
#names(hat_values) <- row.names(model$parameters) | ||
H <- as.matrix(sqrt(W) %*% X %*% solve(t(X) %*% W %*% X) %*% t(X) %*% sqrt(W)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Also can you not just call hatvalues
method on original engine like hatvalues.glm
on glm
's?
"pearsonSTD" = r/sqrt( (1 - hatvalues(object)$selection) * object$selection$variance) | ||
) # TODO studentized_pearson, studentized_deviance | ||
"pearsonSTD" = r/sqrt( (1 - hatvalues(object)$selection) * object$selection$variance), | ||
stop("Invalid type of residual. Choose from 'pearson', 'working', 'deviance', 'response', 'pearsonSTD'.") |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
this switch will result in base R error if type
is NULL
and will return NULL
if type == NA
probably doesn't matter here but its always good to remember that stop
in switch
can break if end user tries to break it :)
class = "family") | ||
} | ||
|
||
gaussian_nonprobsvy <- function(link = "identity") { | ||
mu <- function(eta) eta | ||
variance <- function(mu, y = NULL) rep.int(1, length(mu)) #mean((y - mu)^2) rep.int(1, length(mu)) | ||
variance <- function(mu, y = NULL) mean((y - mu)^2) #rep.int(1, length(mu)) rep.int(1, length(mu)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Just as a reminder you don't have to write all those for every exponential family you want to use you can also just do something simpler like:
gaussian_nonprobsvy <- function(link = "identity") {
x <- stats::gaussian(link = link)
x$mu_der <- function(mu) 1
# Sooner or later you will also probably need a class for those as well
class(x) <- c(class(x), "nonprobsvy_family")
}
No description provided.