-
Notifications
You must be signed in to change notification settings - Fork 0
/
05nPP.Rmd
55 lines (42 loc) · 2.64 KB
/
05nPP.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
# Naive Per Protocol Estimation {-}
The next effect estimate we calculate is the naive Per Protocol (naive PP) estimate. This corresponds to comparing all those who were randomized to receive the treatment, and truly did receive it against those who were randomized to receive the control and truly did receive the control.
The estimation of this treatment effect corresponds to:
$$logit(Y) = \gamma_{0} + \gamma_{A} A + \gamma_{t} t$$
For this estimate, we start by artificially censoring the dataset so that we only retain observations from individuals at time points where they were adherent to their randomized treatment.
```{r, cache=TRUE}
#### First artificially censor the original data within each arm
simulated.data.treated <- simulated.data[simulated.data$Z==1 &
simulated.data$Alag1==1,]
simulated.data.control <- simulated.data[simulated.data$Z==0 &
simulated.data$Alag1==0,]
#### Join these to form the complete artificially censored dataset
simulated.data.combined <- rbind(simulated.data.treated,
simulated.data.control)
head(simulated.data.combined)
```
Next, we fit the pooled logistic regression model:
```{r, cache=TRUE}
naivePPFit <- glm(Y ~ t0 + A,
data = simulated.data.combined,
family = binomial(link="logit"))
naivePPFit
```
Lastly, we calculate the appropriate standard error, p-values, and confidence intervals using the clustered sandwich standard errors:
```{r, cache=TRUE}
cov.m1 <- vcovCL(naivePPFit, type = "HC0",
cluster = simulated.data.combined$id)
co.test <- coeftest(naivePPFit, vcov = cov.m1)
co.CI <- coefci(naivePPFit, vcov = cov.m1, level = 0.95)
cat("Naive PP effect estimate:", "\n",
"log(OR)", round(coefficients(naivePPFit)[["A"]],3), "\n",
"95% CI:", round(co.CI["A","2.5 %"],3), ",",round(co.CI["A","97.5 %"],3), "\n",
"p-value", round(co.test["A","Pr(>|z|)"],3), "\n",
"robust standard error", round(sqrt(diag(cov.m1))[["A"]],3))
```
The naive PP estimate concludes that the effect of treatment is significant, with a log odds ratio of `r coefficients(naivePPFit)[["A"]]`. This result would indicate that the receipt of treatment is associated with an increase in the likelihood of the outcome occurring. However, naive PP estimates may be subject to bias if non-adherence is associated with the outcome of interest.
Alternatively, we can again use the `summ()` function to obtain the same estimates:
```{r, cache=TRUE}
summ(naivePPFit, robust = "HC0",
confint = TRUE, digits = 3, cluster="id",
model.info = FALSE, model.fit = FALSE)
```