Question A

Plot 1: Visit Number and Percent Onchylosis

This plot shows the cross-sectional effect of time on percent of observed moderate-to-severe toe infection. It shows that at Visit 1, nearly 40% were observed to have toe infection, and at visit 7, there was about 10%. In other words, toe infections appear to reduce over time.

Plot 2: Visit Number and Percent Onchylosis By Treatment

Similar to the first plot, this plot also shows the cross-sectional effect treatment over time on percent of observed moderate-to-severe toe infection. The trends by treatment over time are somewhat similar. There was a slighly lower percentage of toe infections among those treated with Terbinafine compared to those treated Itraconazole. In addition, the diminishing effect of time on infections was slightly reduced in those treated itraconazole compared to terbinafine.

Question B:

Consider a marginal model for the log odds of moderate or severe onycholysis.

Question B.1:

Assume “exchangeable” correlations (standard Pearson correlation based specification or logor based specification) for the association among the repeated binary response.

## 
## Call:
## geeglm(formula = Y ~ Month + trt.month, family = "binomial", 
##     data = data, id = ID, corstr = "exchange")
## 
##  Coefficients:
##             Estimate  Std.err   Wald Pr(>|W|)    
## (Intercept) -0.57822  0.13041 19.661 9.25e-06 ***
## Month       -0.17132  0.02957 33.574 6.86e-09 ***
## trt.month   -0.07770  0.05379  2.086    0.149    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Estimated Scale Parameters:
##             Estimate Std.err
## (Intercept)    1.088  0.5265
## 
## Correlation: Structure = exchangeable  Link = identity 
## 
## Estimated Correlation Parameters:
##       Estimate Std.err
## alpha   0.4217  0.2203
## Number of clusters:   294   Maximum cluster size: 7

Question B.2:

What is the interpretation of \(\beta_1\) in this model?

The estimate of \(\beta_1\) = -0.17 (which is the log odds ratio). The OR is then $e^{-0.17} or 0.84. This can be interpreted as the odds ratio of experiencing moderate-to-severe onycholysis for each one month increase for the target population in average. The estimated effect of month is significant at the 0.05 level, and these results provide evidence that the odds of disease severity decreases with each increasing month (at least of the months that were measured in the study) for the population on average.

Question B.3:

What is the interpretation of \(\beta_2\) in this model?

The estimate of \(\beta_2\) = -0.08 (which is the log odds ratio). The OR is then $e^{-0.08} or 0.93. This can be interpreted as the odds ratio of the treatment by month interaction, or the treatment effect on onycholysis with time, for the target population on average. The estimated difference in the slopes is not statistically significant at the 0.05 level, indicating that the outcomes are not changing over time (in other words, the within-subjects effect - the development over time - is not significant) for each of the treatment groups, for the population on average.

Question B.4:

What conclusions can you draw about the effect of treatment on changes in the log odds of moderate or severe onycholysis over time? Provide results that support your conclusions.

I would conclude that the effect of treatment on changes in the log odds of moderate or severe onycholysis over time are not statistically signiciant at the 0.05 level for the target population on average. The \(\beta_2\) for the treatment by month interaction is small, with a p-value > 0.05. Thus, it appears that a marginal logistic regression model without this term is defensible. My interpretation for these results can be reiterated from my answers to (B.3) above: the effect of treatment on onycholysis severity does not significantly vary across the time of the study period (for the target population on average).

Question C:

Consider a generalized linear mixed model, with randomly varying intercepts, for the patient- specific log odds of moderate or severe onycholysis.

Question C.1:

Using Laplace approximation, fit a model with

\[ \begin{aligned} logit \{E(Y_{ij})\} &= \beta_0 + \beta_1 Month_{ij} + \beta_2 Treatment_i \times Month_{ij} \end{aligned} \]

where, given \(b_i\), \(Y_{ij}\) is assumed to have a Bernoulli distribution. Assume that \(b_1 \sim N(0,\sigma^2_b)\).

## Generalized linear mixed model fit by maximum likelihood (Laplace
##   Approximation) [glmerMod]
##  Family: binomial  ( logit )
## Formula: Y ~ Month + trt.month + (1 | ID)
##    Data: data
## 
##      AIC      BIC   logLik deviance df.resid 
##   1263.8   1286.0   -627.9   1255.8     1904 
## 
## Scaled residuals: 
##    Min     1Q Median     3Q    Max 
##  -3.25  -0.15  -0.07  -0.01  48.20 
## 
## Random effects:
##  Groups Name        Variance Std.Dev.
##  ID     (Intercept) 20.6     4.54    
## Number of obs: 1908, groups:  ID, 294
## 
## Fixed effects:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.6491     0.7087   -3.74  0.00019 ***
## Month        -0.3956     0.0457   -8.66  < 2e-16 ***
## trt.month    -0.1455     0.0669   -2.18  0.02952 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##           (Intr) Month 
## Month      0.248       
## trt.month  0.073 -0.513

Question C.2:

What is the estimate of \(\sigma^2_{b}\)? Give an interpretation to the magnitude of the estimated variance?

The estimate of \(\sigma^2_{b}\) is 20.6.Therefore, \(\hat{\sigma_b}\) is \(\sqrt{20.6}\) = 4.54.

The \(\sigma^2_{b}\) is large, indicating that there is substantial variability to experience toe infection.

The esimated propensity of toe infection at baseline is wide: the 95% range of propensity of experiencing moderate-to-severe disease at baseline is:

\[\left( \frac{exp(\hat{\beta}_{0}-1.96*\hat{\sigma}_b)}{1+exp(\hat{\beta}_{0}-1.96*\hat{\sigma}_b)}, \frac{exp(\hat{\beta}_{0}+1.96*\hat{\sigma}_b)}{1+exp(\hat{\beta}_{0}+1.96*\hat{\sigma}_b)} \right)\]

\[\left( \frac{exp(-2.64911-1.96*4.5387223)}{1+exp(-2.64911-1.96*4.5387223)}, \frac{exp(-2.64911+1.96*4.5387223)}{1+exp(-2.64911+1.96*4.5387223)} \right)\]

=9.684195210^{-6},0.9980671

Question C.3:

What is the interpretation of \(\beta_1\) in this model?

The estimate of \(\beta_1\) = -0.4 (which is the subject-specific log odds ratio). The OR is then $e^{-0.4} or 0.67. This can be interpreted as the subject-specific odds ratio of experiencing moderate-to-severe onycholysis for each one month increase. The estimated within-subject effect of month is significant at the 0.05 level, and these results provide evidence that the subject-specific odds of disease severity decreases with each increasing month (at least of the months that were measured in the study).

Question C.4:

What is the interpretation of \(\beta_2\) in this model?

The estimate of \(\beta_2\) = -0.15 (which is the log odds ratio). The OR is then $e^{-0.15} or 0.86. This can be interpreted as odds ratio of experiencing the disease between two individuals who happen to have the same value for the unobserved random effects (i.e., \(b_i=b_{i'}\)), but who differ by one unit in the covariate (e.g., one using the Terbinafine treatment, and one using the Itraconazole treatment). This estimated effect is significant at the 0.05 level.

Question C.5:

Suppose for id=18, the predicted random effect \(\hat{b}_{id=18}\) = 1.141, calculate the predicted subject-specific (conditional) probability of having toe infection at week 10 using his/her assigned group.

NOTE: I had a different predicted random effect estimate for subject 18, 2.061876! Therefore, my other model estimates might be wrong.

\(\hat{\beta}_{0}\) = -2.6491067

\(\hat{b}_{id=18}\) = 1.14

\(\hat{\beta}_{1}\) = -0.395591

\(\hat{\beta}_{10 weeks}\) = -0.395591 * (10 weeks/4 weeks per month)

\(\hat{\beta}_{2}\) = 0 (since the treatment for this group has been assigned to 0 - Itraconazole)

\[ \begin{aligned} \hat{Y}_{18} = (\hat{\beta}_{0}+\hat{b}_{id=18}) + \hat{\beta}_{1}*(10/4) \end{aligned} \]

\[ \begin{aligned} \hat{Y}_{18} = -2.6491067 + 1.14 + -0.9889776 = -2.4970843 \end{aligned} \]

The odds ratio = exp(-2.4970843) = 0.0823247

For subject 18, the predicted subject-specific (conditional) probability of having a toe infection at week 10 is

\[ \frac{0.0823247}{1+0.0823247} = 0.0760588,~i.e~7.61\% \]

Question D:

Compare and contrast the estimates of \(\beta_2\) from the marginal and mixed models. Why might they differ?

The \(\beta_2\) for the marginal model was -0.08, which was not statistically significant at the 0.05 level. The \(\beta_2\) for the mixed model was -0.15, which was statistically significant at the 0.05 level. These estimates could differ, because the mixed model controls for the random effects when estimating the fixed effects, and in contrast in the marginal model, the coefficients are directly estimated among each participant’s residuals. As such, the interpretations for \(\beta_2\) are different for the marginal and mixed models. In the marginal model the beta estimates the log odds ratio of the two drugs over time for the target population on average. In the mixed model, the beta estimates the log odds change by one unit difference in the covariate for a given subject over time (in other words, for subjects who have the same value for unobserved random effects but take different treatments).

Question E:

Repeat the analysis from problem (c) sequentially increasing the number of quadrature points used. Compare the estimates and standard errors of the model parameters when the number of quadrature points is 3, 9, 15, 21, 49. Do the results depend on the number of quadrature points?

GLMM with increasing number of quadrature points
Dependent variable:
Y
nAGQ=1 nAGQ=3 nAGQ=9 nAGQ=15 nAGQ=21 nAGQ=49
(1) (2) (3) (4) (5) (6)
Month -0.396*** -0.402*** -0.385*** -0.390*** -0.388*** -0.388***
(0.046) (0.046) (0.043) (0.044) (0.043) (0.043)
trt.month -0.146** -0.145** -0.141** -0.143** -0.142** -0.142**
(0.067) (0.068) (0.064) (0.065) (0.065) (0.065)
Constant -2.649*** -2.127*** -1.651*** -1.729*** -1.695*** -1.697***
(0.709) (0.503) (0.312) (0.341) (0.330) (0.330)
Observations 1,908 1,908 1,908 1,908 1,908 1,908
Log Likelihood -627.915 -631.119 -625.601 -625.250 -625.465 -625.436
Akaike Inf. Crit. 1,263.829 1,270.238 1,259.202 1,258.500 1,258.929 1,258.871
Bayesian Inf. Crit. 1,286.044 1,292.453 1,281.417 1,280.715 1,281.144 1,281.087
Note: p<0.1; p<0.05; p<0.01

The beta estimates and standard errors for the model parameters vary very slightly with increasing number of quadrature points (see above table). The intercept/constat values, however, decreased from -2.650 (SE=0.709) to -1.700 (SE=0.330). The estimates and standard errors seemed to settle (and decrease somewhat) from nAGQ=1 to nAGQ=21,49. The random effect standard deviation decreases from 20.6 to 16 from nAGQ’s 1 to 21. However, the fixed effect and random effect standard deviation is smallest when there are 9 quadrature points. If there is high random effect variance, then it would appear that a higher number of quadrature points would better approximate the estimates. But I am not confident to say whether or not nAGQ=9 should be used, although the AIC is very slightly smaller in the models with a greater number of quadrature points.