Linear regression of mixed effects and population attributable fraction (PAF) in randomized clinical trials

[Andressa Freire]

Hello Stata community

I have a database in long format (as follows in the dataex example) from a randomized clinical study carried out with children (parti) where the outcome is a continuous variable of nutritional status (cc_BMI_zscore_). The intervention group was exposed to lifestyle changes such as consuming a Mediterranean diet and promoting physical activity. In both groups, measures to calculate BMI and questionnaires to assess adherence to the Mediterranean diet (Dmed_score - continuous variable) and degree of physical activity (AF_score - continuous variable) at time 0 (baseline) and time 1 (12 months later).

I would like to use mixed effects linear regression to evaluate changes in BMI over times 0 and 1 in both groups and see if the Mediterranean diet and physical activity influenced this outcome. Furthermore, I would like to know if it is possible and how to calculate the population attributable fraction (PAF) in this type of study to see how much of this change in BMI can be attributable to the Mediterranean diet and/or physical activity. Would it also be interesting to create categories for the AF_score Dmed_score variables?

I would really like your help
Stata version 17

----------------------- copy starting from the next line -----------------------

* Example generated by -dataex-. To install: ssc install dataex
clear
input double parti byte tempo double(Grupo AF_score cc_BMI_ cc_BMI_zscore_ cc_BMI_cat_ Mean_WC_ Dmed_score) float edad_calc_
1003 0 2  5 15.86   .21  0               49.6 12 4.4462695
1003 1 2 11 15.45     0  0  52.56666666666666 11  5.483915
1005 0 1  9  14.8  -.46  0               49.7 10  5.667351
1005 1 1 12 16.41     0  0  56.56666666666666 14  6.814511
1006 0 1  6 16.34   .48  0                 56 11  4.073922
1006 1 1  7 17.19     1  1  57.53333333333334 11  5.204654
1010 0 2 11 14.28  -1.1  0               45.1 11  5.456537
1010 1 2 11 14.91     0  0 53.866666666666674 13  6.494182
1021 0 1  8 17.21  1.07  1              55.45  7  4.895277
1021 1 1 12 16.91     0  0  58.46666666666667 12   6.02601
1023 0 2 10 18.04  1.18  1              58.25  9  6.995209
1023 1 2  8 18.43     1  1               63.5  8   8.03833
1024 0 2  9 15.08  -.29  0              47.45 14  6.294319
1024 1 2 12 16.15     0  0               54.6 17  7.315537
1025 0 1 10 14.54 -1.03  0               45.1 18  4.273785
1025 1 1  9 13.67    -1 -1 48.666666666666664 16  5.412731
1028 0 1 10 15.57  -.09  0               51.6 11  3.852156
1028 1 1 11 15.22     0  0  52.23333333333333 12  5.037645
1037 0 1 11 17.22  1.04  0               53.3 12 4.0985627
1037 1 1 10 17.16     1  1               54.9 14  5.190965
1038 0 2  7  21.7  2.78  2               63.4  7  7.033539
1038 1 2  7 24.35     2  2                 71 13  8.123203
1039 0 1  7 16.81   .81  0              57.75 13  5.642711
1039 1 1  7 17.54     0  0               59.9 14  6.691308
1040 0 1  8 16.27   .39  0              51.25  8 3.6769335
1040 1 1  7 16.42     0  0               55.4 12 4.6844625
1041 0 2 11  14.8  -.58  0               56.7 13  5.845312
1041 1 2 11    19     1  1               62.5 11  6.913073
1042 0 2  9 14.88  -.62  0               54.8 12 4.6735115
1042 1 2 11 28.05     2  2 54.333333333333336 10  5.670089
1046 0 1 13 14.85  -.43  0              50.25 15  5.262149
1046 1 1 10 15.25     0  0               53.5 13  6.351814
1047 0 2  9 25.45  4.61  2               58.2  7 4.3340178
1047 1 2  8 19.85     2  2              66.65  6  5.418207
1048 0 2  9 15.13  -.66  0               49.5 12   3.59206
1048 1 2  9 16.37     0  0  52.03333333333333 12  4.643395
1049 0 2 11 16.43   .59  0               54.8 11  5.771389
1049 1 2  9 17.39     0  0  56.76666666666667 12  6.822724
1053 0 2  8 14.03 -1.29 -1               52.4 11  4.219028
1053 1 2 11  16.1     0  0               52.6  7  5.289528
1054 0 1 11 13.34 -1.96 -1               55.5  5  7.301848
1054 1 1 10 13.67    -1 -1               54.6 10  8.314853
1055 0 1 10 14.47  -.91  0               53.4  5  5.494866
1055 1 1  9 14.17     0  0                 58 11  6.507871
1056 0 1  6 15.35  -.29  0                  .  8  3.649555
1056 1 1  7 17.21     1  1               54.2 10  4.605065
1058 0 2  9 14.17  -.96  0              51.05 11  5.316906
1058 1 2  7 15.04     0  0               50.5 10  6.502396
1059 0 1  9 14.44  -.88  0               51.2 10  6.650239
1059 1 1 10 14.28     0  0              53.35 12   7.66872
end
label values Grupo Grupo
label def Grupo 1 "Intervención", modify
label def Grupo 2 "Control", modify

------------------ copy up to and including the previous line ------------------

[Carlo Lazzaro]

Andressa:
do you mean something along the following lines?
[CODE][. mixed cc_BMI_zscore_ i.tempo c.Dmed_score##c.AF_score || Grupo:

Performing EM optimization ...

Performing gradient-based optimization:
Iteration 0: log likelihood = -70.643888
Iteration 1: log likelihood = -70.643835
Iteration 2: log likelihood = -70.643835

Computing standard errors ...

Mixed-effects ML regression Number of obs = 50
Group variable: Grupo Number of groups = 2
Obs per group:
min = 24
avg = 25.0
max = 26
Wald chi2(4) = 13.36
Log likelihood = -70.643835 Prob > chi2 = 0.0096

-----------------------------------------------------------------------------------------
cc_BMI_zscore_ | Coefficient Std. err. z P>|z| [95% conf. interval]
------------------------+----------------------------------------------------------------
1.tempo | .5010304 .285549 1.75 0.079 -.0586355 1.060696
Dmed_score | -.724814 .297876 -2.43 0.015 -1.30864 -.1409877
AF_score | -.8618048 .3574404 -2.41 0.016 -1.562375 -.1612345
|
c.Dmed_score#c.AF_score | .0646116 .0304226 2.12 0.034 .0049845 .1242387
|
_cons | 9.273637 3.416893 2.71 0.007 2.57665 15.97062
-----------------------------------------------------------------------------------------

------------------------------------------------------------------------------
Random-effects parameters | Estimate Std. err. [95% conf. interval]
-----------------------------+------------------------------------------------
Grupo: Identity |
var(_cons) | .0249118 .064622 .0001543 4.022078
-----------------------------+------------------------------------------------
var(Residual) | .9685254 .197756 .6490983 1.445145
------------------------------------------------------------------------------
LR test vs. linear model: chibar2(01) = 0.27 Prob >= chibar2 = 0.3013

.
/CODE]

[Andressa Freire]

Hi Carlo, thank you very much for your help here. In fact, I feel a little lost where to allocate each variable in the mixed code. But I'll try to explain here what I thought... and if you can tell me something, I'd be very grateful.

First I thought I could analyze physical activity and the Mediterranean diet separately in the code. Taking the Mediterranean diet as an example, I created a dichotomous variable (dmed_cat) from the numeric variable and tested the following code:

gen Dmed_cat=.
replace Dmed_cat=1 if Dmed_score<9
replace Dmed_cat=2 if Dmed_score>=9

mixed cc_BMI_zscore_ i.Grupo##i.tempo##i.Dmed_cat || parti:


my question, for example, is whether I should put II Group: or || parti: in the code...
if it were ||Grupo:, then it would be:

mixed cc_BMI_zscore_ i.tempo##i.Dmed_cat || Grupo:

[Jeff Wooldridge]

Andressa: Can you explain why you want to use a mixed model? If you have randomized treatments, you can simply regress cc_BMI_zscore on your continuous treatment measures, probably interact them, maybe put them in squared, and also include the baseline BMI measure. I don't see what using a mixed model buys you. You'll estimate the marginal effects of each treatment on BMI, and I would think that's what you want to think about effectiveness of the program.

[Erik Ruzek]

If you want to use a mixed model (see post #4), then Grupo should be in the fixed part of the model (everything before the || in the code) and the childid (parti) should go in the random part of the model (after the || in the code).

** Data management
recode Grupo (2 = 0), gen (grupo)     // control == 0, treatment == 1
gen Dmed_cat=.
replace Dmed_cat=1 if Dmed_score<9
replace Dmed_cat=2 if Dmed_score>=9

** Mixed model
mixed cc_BMI_zscore_ i.Dmed_cat#i.grupo AF_score || parti: , reml 

** Check residuals 
predict sresidz, rstandard 
predict xbz, xb
scatter sresidz xbz           // looks a bit off, but maybe more data will resolve this
qnorm sresidz                  // looks ok

** Graph the interaction
margins Dmed_cat#grupo
marginsplot

I am not familiar with the population attributable fraction (PAF), but noticed that Roger Newson has written a program for this punaf. I'm not sure if that will work with mixed effects models.

[Andressa Freire]

Thank you Carlo, Jeff e Erik for the great insights that helped me develop the analyzes I needed