#' Asymptotic variance of population average potential outcome vector for
#' different alphas.
#'
#' @param ygroup An matrix including the group average potential outcome
#' estimates where rows correspond to group, and columns to values of alpha.
#' @param ps String. Can take values 'true', or 'estimated' for known or
#' estimated propensity score. Defaults to 'true'.
#' @param scores A matrix with rows corresponding to the parameters of the
#' propensity score model and columns for groups. Includes the score of the
#' propensity score evaluated for the variables of each group. Can be left
#' NULL for ps = 'true'.
#'
#' @export
IEvar <- function(ygroup, ps = c('true', 'estimated'), scores = NULL) {
ps <- match.arg(ps)
n_neigh <- dim(ygroup)[1]
alpha <- as.numeric(dimnames(ygroup)[[2]])
ie_var <- cov(ygroup)
ie_var <- ie_var * (n_neigh - 1) / (n_neigh ^ 2)
if (ps == 'true') {
return(ie_var)
}
# Else, we are based on the estimated propensity score.
if (is.null(scores)) {
stop('Provide score matrix.')
}
var_est_ps <- array(0, dim = dim(ie_var))
num_gamma <- dim(scores)[1]
# --- Calculating B11, the information matrix of the cluster ps.
B11 <- matrix(0, nrow = num_gamma, ncol = num_gamma)
for (nn in 1 : n_neigh) {
scores_nn <- scores[, nn, drop = FALSE]
B11 <- B11 + scores_nn %*% t(scores_nn)
}
B11 <- B11 / n_neigh
B11_inv <- chol2inv(chol(B11))
# Calculating C21, and D12.
ypop <- apply(ygroup, 2, mean)
C21 <- array(0, dim = c(length(alpha), num_gamma))
D12 <- array(0, dim = c(num_gamma, length(alpha)))
for (nn in 1 : n_neigh) {
C21 <- C21 - t(ygroup[nn, , drop = FALSE]) %*% t(scores[, nn, drop = FALSE])
D12 <- D12 + scores[, nn, drop = FALSE] %*% (ygroup[nn, , drop = FALSE] - ypop)
}
C21 <- C21 / n_neigh
D12 <- D12 / n_neigh
chol_B11_inv <- chol(B11_inv)
mat1 <- C21 %*% t(chol_B11_inv) %*% t(C21 %*% t(chol_B11_inv))
mat <- C21 %*% B11_inv %*% D12
var_est_ps <- mat1 + mat + t(mat)
var_est_ps <- ie_var + var_est_ps / n_neigh
return(var_est_ps)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.