I have a dataset which I think requires a multivariate multilevel analysis. I am unsure both of the appropriate model and of how to fit it with
R. I have come up with a tentative model, but my understanding of the math is so superficial that I cannot tell whether my analysis is "right" or whether it includes blatant errors. I would appreciate any insight on the model design or the model specification in
The study design
The question is whether the architectural design of a clinic will influence the outcome of a pathology for permanent residents in this clinic. We have collected data on 13 symptoms for 8 patients per clinic in 21 clinics.
There is a patient-level IV
medication and a clinic-level IV
architecture. All variables are continuous-ish.
The 13 symptoms are correlated +.20 on average, which I think indicates a multivariate multilevel analysis is appropriate.
To run the multivariate analysis with
nlme I have standardized my DVs, stacked these 13 DVs in a single column, and added a categorical dummy variable to flag which row corresponds to which symptom.
It looks like this:
Clinic Patient Symptom Score Medication Architecture 1 1 EP1 0.12 1 3.2 1 1 EP2 0.11 1 3.2 1 1 EP3 0.13 1 3.2 1 2 EP1 0.56 4 3.2 1 2 EP2 0.67 4 3.2 1 2 EP3 0.23 4 3.2 2 3 EP1 0.22 3 5.1 2 3 EP2 0.25 3 5.1 2 3 EP3 0.14 3 5.1 2 4 EP1 0.78 6 5.1 2 4 EP2 0.89 6 5.1 2 4 EP3 0.11 6 5.1
The model design
- To run the analysis as multivariate, I use both
symptom:architectureas IVs and I remove the intercept in both the fixed and random parts of the model. I do not include the main effect of
architectureas an IV.
- The effect of
medicationshould be the same within all clinics, so there is no random effect for this variable.
- I do not want to constrain equality between the effect of
architectureon the different
- Due to the multivariate nature of the analysis, I expect the residuals to be correlated, with different correlations between the 13 different symptoms; therefore I specify the covariance structure of residuals as
corSymm(non-zero but unstructured, if I get this correctly).
- I also expect heteroscedasticity between the different symptoms (there should be more variance on certain symptoms), so I add the option
The end result
This is the model I come up with:
model1 = lme(fixed = Score ~ symptom + medication:symptom + architecture:symptom + medication:architecture:symptom - 1, + random = ~ symptom - 1 | patient/clinic, + correlation = corSymm, + weights=varIdent(form= ~ 1|symptoms) + method = "ML")
In order to test the effect of the architectural variables, I would then compare this model to the following constrained model, dropping all the terms related to architecture:
model2 = lme(fixed = Score ~ symptom + medication:symptom - 1, + random = ~ symptom - 1 | patient/clinic, + correlation = corSymm, + weights=varIdent(form= ~ 1|symptoms) + method = "ML")
I would then run this comparison with the command
anova(model1, model2) and compare the log-likelihood of the two models.
Overall, do these model design and r specification look correct to you? Thank you so much for your help!