5.8 Comparing IPW with Aggregated Estimate

We know that the probablity of the treatment increases with \(X1\), as shown below.

plot(X[, 1], prob)
par(new = T)
abline(v = 0, col = "red", lty = "dashed")

Next, we would want to divide \(X1\) into segments such that the probability of treatment remains more or less similar in each segment. Since we generated the data ourselves, we know that the probability of treatment increases for observations with \(X1 > 0\). However, lets take 10 segments, which cuts the distribution of \(X1\) in decile.

quant  <- quantile(X[, 1], p = seq(0, 1, 0.1))
group  <- cut(X[, 1], quant, include.lowest = TRUE)
dat$group  <- group

Next, we want to estimate the treatment effects for each segment and then aggregate it.

dat_sum   <- data.frame(dat  %>% 
                mutate(status = ifelse(W == 1, "treatment", "control"))  %>% 
                group_by(group)  %>% 
                summarize(num_treat = sum(W), n_group = n()))  %>% 
                mutate(prop_treat = num_treat / n_group) %>%  
                dplyr::select(c(group, prop_treat))


dat   <- dat  %>%  
            merge(dat_sum, by = "group", all.x = T)  %>%  
            mutate(score = Y * ((W / prop_treat) - ((1 - W) / (1 - prop_treat))) 
            )


paste("aggregated treatment effect is: ", round(mean(dat$score), 3))

## [1] "aggregated treatment effect is:  15.019"

paste("se of aggregated treatment effect is: ", round(sd(dat$score) / sqrt(length(W)), 3))

## [1] "se of aggregated treatment effect is:  3.122"

Now, let’s estimate the oracle IPW with known propensity score.

score_oracle  <- ((W * Y) / (prob)) - 
                ((1 - W) * Y / (1 - prob))
ipw.oracle  <- mean(score_oracle)
se.oracle  <- sd(score_oracle) / sqrt(length(score_oracle))

paste("oracle IPW estimate: ", round(ipw.oracle, 3))

## [1] "oracle IPW estimate:  16.496"

paste("standard error: ", round(se.oracle, 3))

## [1] "standard error:  3.165"