Help with interpreting xtmixed output

[Jem Lane]

Hi,
I'm using Stata 12IC for Windows. I'm trying to perfom multilevel mixed effects modelling on a dataset where there are 2 groups of samples (id) defined by GRP (1/0). I have made repeated measurements on these samples to obtain the dependent variable OUT, at increasing concentrations of a drug (CONC 0-5) (CONC was increased sequentially at each timepoint, therefore this variable serves as both timepoint and concentration). I am treating CONC as a continuous variable here. I have chosen to have GRP as an explanatory variable/fixed effect as there are only 2 possible groups.
I have used the following command for a random slopes model:

xtmixed OUT i.GRP##c.CONC || id: CONC, mle variance

Please see attached pictures of the graph for the mean responses of the two groups, and the output from this command. I would be grateful for help in checking I have understood the output correctly. Am I correct in thinking the following (referring to coloured circles)?

1 (red): the intercept of GRP0 at CONC 0
2 (green): the p-value for null hypothesis that there is no difference between GRP0 (the reference) and GRP1. This indicates the null hypothesis should be rejected: there is a difference. But is this just at baseline, or across all concentrations?
3 (blue): I am not sure here. Is it the p-value for the null hypothesis that there is no difference in OUT across all concentrations, but for both groups combined? And as it is <0.05, this suggests OUT does vary with CONC in some way?
4 (purple): p-value for the interaction - highly significant, suggesting there is an interaction and that OUT varies with CONC differently for each group

Thanks
Jem

[Rich Goldstein]

the attached pictures do not show up (at least on my system); please read the FAQ which discussed how to show output, etc.

[Carlo Lazzaro]

jeml (as per FAQ, please note the preference for real full names and surnames in this forum. Just click on the Contact us button at the bottom of the screen and re-register, Thank you):
unfortuntely, your images are not readable (as it is often the case), for me at least.
Your chance to get helpful replies would increase if you post the results that Stata gave you, too.
You can copy all your Stata session on a .do file, copy and past it between the hash-shaped code delimiters (as you did for the Stata code in your message).
Usually this approach works out fine.

[Jem Lane]

Thanks Rich and Carlo - that's a shame the images don't show for you. I annotated the output to make my query clearer, but as you're not able to view it I will look into other ways to display and repost.
Will also amend my name on my profile.
regards

Jem Lane

[Jem Lane]

Hi Rich / Carlo,
I have contacted Statalist admin to request a change of username - awaiting response. To recap (and avoid need to scroll), I'm trying to perfom multilevel mixed effects modelling on a dataset where there are 2 groups of samples (id) defined by GRP (TOP/BOTTOM). I have made repeated measurements on these samples to obtain the dependent variable OUT, at increasing concentrations of a drug (CONC 0-5) (CONC was increased sequentially at each timepoint, therefore this variable serves as both timepoint and concentration). I am treating CONC as a continuous variable here. I have chosen to have GRP as an explanatory variable/fixed effect as there are only 2 possible groups. Below is the command I used for a random slopes model with the output:

. xtmixed OUT i.GRP##c.CONC || id: CONC, mle variance

Performing EM optimization: 

Performing gradient-based optimization: 

Iteration 0:   log likelihood = -389.52551  
Iteration 1:   log likelihood = -389.52551  

Computing standard errors:

Mixed-effects ML regression                     Number of obs      =       114
Group variable: id                              Number of groups   =        19

                                                Obs per group: min =         6
                                                               avg =       6.0
                                                               max =         6


                                                Wald chi2(3)       =     20.46
Log likelihood = -389.52551                     Prob > chi2        =    0.0001

------------------------------------------------------------------------------
         OUT |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
       1.GRP |   18.73228   8.295914     2.26   0.024     2.472583    34.99197
        CONC |  -1.495238   .5285768    -2.83   0.005     -2.53123   -.4592467
             |
  GRP#c.CONC |
          1  |   2.635238   .7285928     3.62   0.000     1.207222    4.063254
             |
       _cons |   38.20106   6.018489     6.35   0.000     26.40504    49.99708
------------------------------------------------------------------------------

------------------------------------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
id: Independent              |
                   var(CONC) |   1.172252   .8348628      .2902727    4.734082
                  var(_cons) |   313.6955   104.3665      163.4232    602.1478
-----------------------------+------------------------------------------------
               var(Residual) |   23.49005   3.831319      17.06277    32.33838
------------------------------------------------------------------------------
LR test vs. linear regression:       chi2(2) =   208.38   Prob > chi2 = 0.0000

Note: LR test is conservative and provided only for reference.

Could I check my understanding of the following 4 things in the output?

1 (red): the intercept of GRP0 at CONC 0
2 (green): the p-value for null hypothesis that there is no difference between GRP0 (the reference) and GRP1. This indicates the null hypothesis should be rejected: there is a difference. But is this just at baseline, or across all concentrations?
3 (blue): I am not sure here. Is it the p-value for the null hypothesis that there is no difference in OUT across all concentrations, for both groups combined? And as it is <0.05, this suggests OUT does vary with CONC in some way?
4 (purple): p-value for the interaction - highly significant, suggesting there is an interaction and that OUT varies with CONC differently for each group

thanks

Jem Lane


P.S. To give you a visual impression of the results, here (hopefully) is a graph of mean values for each group, measured repeatedly at several drug concentrations/timepoints (mCONC). OUT is the ordinate (in png format this time):

[Carlo Lazzaro]

Jem:
thanks for re-registering.
About your questions concernong the fixed effect part of your multilevel mixed-effects linear regression
1 (red): the intercept of GRP0 at CONC 0. Your interpretation is right; however, it is probably less meaningful than expected, if you do not have CONC=0 in your dataset. A suitable option is to center that variable around a meningful value of CONC.
2 (green): there is a difference between the two groups across all concentrations.
3 (blue): OUT does vary with CONC, no matter which group you consider.
4 (purple): OUT varies with CONC differently for each group.

As an aside, you can add -baselevels- option at the end of your code to make things clearer:

xtmixed OUT i.GRP##c.CONC || id: CONC, mle variance baselevels

[Clyde Schechter]

I disagree with two of Carlo's interpretations. The use of an equation with a dichotomous grouping variable (GRP, coded 0 vs1) and a continuous predictor (CONC) and their interaction is similar to having two equations:

Group 0: OUT = b0 + b1*CONC + error
Group 1: OUT = c0 + c1*CONC + error

where
b0 = _cons, b1 = _b[CONS]
c0 = _cons + _b[1.GRP], c1 = _b[CONS] +_b[1.GRP#CONS]

(Note: I'm abstracting away from the random slopes and intercepts by just writing an undifferentiated "error" term. Nothing I say here would be different were I to explicitly represent the error structure.)

Keeping these two equations in mind, the green p-value refers to the coefficient of 1.GRP, which, in the two equation framework is the expected difference between the group outcomes when CONC = 0; it does not follow, in general, that the group differences persist at other levels of CONC. (As can be seen in your graphs, the differences do persist and even increase, but the green p-value does not tell you that. If your interaction term coefficient had turned out to be a large negative, then the OUT vs CONC slope in group 1 would be negative, and the two group trajectories could have crossed each other at some point, meaning that at that point there would be no group difference in expected outcome OUT.)

The blue p-value refers to the coefficient of CONC, which is b1 in the two equation framework. That means that that the (expected) slope of the OUT to CONC relationship is both large enough and precisely measured enough that we can be confident it is not zero in GROUP 0. In the general situation, it does not follow that the slope of the OUT to CONC relationship is different from zero in group 1. (Although, again, your data do show a non-zero relationship in group 1 as well. But my point is that the blue p-value doesn't tell you that.) If your data were different, it could have been the case that the interaction coefficient _b[1.GRP#CONC] could have been negative and roughly equal in magnitude to the positive _b[CONC]. Had that happened, you would have a roughly 0 relationship between CONC and OUT in group 1 in spite of what happened in group 0.

You can probably get a better understanding of what your model is telling you (and how well it fits your data) by running the -margins- command and -marginsplot- thereafter. And if you would like direct output of the expected OUT vs CONC slopes in each group, -margins- with the -dydx- option will give it to you.

[Roman Mostazir]

For interpretation convenience, I take "OUT" to be the "Blood Pressure". The red coefficient tells you the mean level of blood pressure for the Top group (Group==0). The green tells you your Bottom group (Group==1) has a higher blood pressure by a magnitude of 18.73 compare to the Top (your graph says that).

Now, I think these two margins command will help you to understand your model better:

Margins command: 1

margins,dydx(CONC)at(GRP=(0 1))

You should be presented with two coefficients with associated other estimates:

1. GRP_0= -1.495 (The blue coefficient)

2. GRP_1= -4.13

The interaction effect in the model gives you the difference between this two (-1.495)-(-4.13) = 2.64. Thus, the blue coefficient -1.495 is actually giving you the effect of Concentration on Group=0 (Top) when Concentration is 'zero', meaning you expect the Concentration to decreases the blood pressure on Top group. HANG ON FOR A MOMENT !!! DON'T INTERPRET YET what will happen if the Concentration increases, as the increase of Concentration has a different effect at two levels of the groups (significant interaction). Therefore we interpret the effect of increasing Concentration with the maroon coefficient.

The maroon (interaction) coefficient tells you, as you increase the concentration, the blood pressure on Bottom group (GROUP-1) decreases by a greater magnitude of 2.64 compare to GROUP-0 and that happens linearly for each level of Concentration. Therefore, the 2nd margins command you need to understand this is below:

Margins command: 2

margins i.GRP, at(CONC=(0(1)6))

marginsplot, xdim(CONC)

This plot will tell you a better picture, however, for group comparison tests i.e. whether the groups at each concentration are different or not, you need to use 'contrast'. That can be tested later, first try those margins command and understand your model.

[Carlo Lazzaro]

Jem:
Clyde is right.
I left my desk and forgot to correct this part of my previous reply.
Sorry for any inconvenience my delay could have caused to Jem.

[Jem Lane]

Dear Carlo, Clyde and Roman,

thank you all for your input - I think I'm almost there with your help.

Clyde - your detailed explanation was really useful, particularly in terms of explaining the 'what if' scenarios regarding what the output would mean should the data have been different. For completeness, could I ask what is meant by "_b" and "_b1" in your two equations:

b0 = _cons, b1 = _b[CONS]
c0 = _cons + _b1[1.GRP], c1 = _b[CONS] + _b[1.GRP.CONS]

Do they have a value, and if so, where are they in the output? Is it analogous to u0_ij and u1_ij in alternative statistical language?


Carlo - thanks for the margins code. I ran it, but got different values to those you quoted. I got

GRP 0 = -1.495
GRP 1 = 1.14 (from the original output, CONC + GRP#c.CONC = -1.495 + 2.635 = 1.14. Perhaps this is what you meant?)

So these values are the average slope for each group across all values of CONC, right (after taking into account interaction effect)?

I ran

margins i.GRP, at (CONC=(0(1)6))

as you suggested. This gives what I presume are the values predicted by the model for each group at each concentration, right? But for some reason, although I have 6 concentrations (0-5), the output gave 7 'margins' for each group (output omitted here for brevity). Any ideas?

The marginsplot command plots these predicted values. So, I have a couple of questions re: margins:
i. the confidence intervals on marginsplot overlap for CONC 0 and 1 (see graph below). This seems strange given that in my original output, the green p-value (for 1.GRP) was 0.024. Is it because the overlap is small, and therefore the difference at CONC = 0 is still significant?
ii. the values predicted do not exactly match the mean values obtained from the original data at CONC =0 and 5 for example - is a small amount of error to be expected?

Finally, yes, I do want to compare groups at each concentration to look for significant differences esp. at the final CONC (similar to a post-hoc test with ANOVA?). So is the contrast command the best for this? The handbook explanation if pretty brief... I tried running

contrast GRP##CONC

but I got
"factor CONC not found in e(b)
r(111);"

If 'contrast' is the correct command for what I want, what code should I use?

Apologies for the long post. Many thanks again

Jem

marginsplot showing overlap of CIs at CONC 0 & 1

[Roman Mostazir]

Carlo - thanks for the margins code. I ran it, but got different values to those you quoted. I got

GRP 0 = -1.495
GRP 1 = 1.14 (from the original output, CONC + GRP#c.CONC = -1.495 + 2.635 = 1.14. Perhaps this is what you meant?)

So these values are the average slope for each group across all values of CONC, right (after taking into account interaction effect)?

I think you meant my code rather Carlo. So I am taking the responsibility of the reply. Yes, that's right. I wasn't sure, as I don't have the data, that which way Concentration is affecting the Goup-1. I assumed the same direction as Group-0, but clearly they are being affected in opposite direction. Your recent margins plot says that. Thus, yes, 1.14 is the right coefficient for the Group-1.

I ran

margins i.GRP, at (CONC=(0(1)6))

as you suggested. This gives what I presume are the values predicted by the model for each group at each concentration, right? But for some reason, although I have 6 concentrations (0-5), the output gave 7 'margins' for each group (output omitted here for brevity). Any ideas?

Since you don't have 6 in your Concentration observation your code should be:

margins i.GRP, at (CONC=(0(1)5))

This will give you 6 levels of estimations for both groups.

The marginsplot command plots these predicted values. So, I have a couple of questions re: margins:
i. the confidence intervals on marginsplot overlap for CONC 0 and 1 (see graph below). This seems strange given that in my original output, the green p-value (for 1.GRP) was 0.024. Is it because the overlap is small, and therefore the difference at CONC = 0 is still significant?

ii. the values predicted do not exactly match the mean values obtained from the original data at CONC =0 and 5 for example - is a small amount of error to be expected?

Finally, yes, I do want to compare groups at each concentration to look for significant differences esp. at the final CONC (similar to a post-hoc test with ANOVA?). So is the contrast command the best for this? The handbook explanation if pretty brief... I tried running

contrast GRP##CONC

but I got
"factor CONC not found in e(b)
r(111);"


If 'contrast' is the correct command for what I want, what code should I use?

For confidence interval overlaps, you are certain that groups are different if the bars don't overlap, but the opposite is not always true. Therefore, you need to test the groups at each level of concentration. Though it is certain at other levels (from you graph) that they are different, yet, to get the p' values for each level of Concentration for two groups:

margins r.GRP, at (CONC=(0(1)5)) contrast

Note the r. notation instead i. before GRP. This will produce the test group difference at each level of Concentration and their S,E and CI around the differences.

For understanding the slight variation in the margins prediction, you need to understand the several competing theories behind margins and other prediction types. Richard Williams has a splendid documentation on this, and I highly recommend this https://www3.nd.edu/~rwilliam/stats/Margins01.pdf . Also Maarten Buis has another superb one http://www.stata-journal.com./article.html?article=st0127. Finally, to understand the Contrast tests, please do look at the manual. Help 'contrast' should take you there. In brief, yes contrast produces test results equivalent to ANOVA.

[Clyde Schechter]

For completeness, could I ask what is meant by "_b" and "_b1" in your two equations:

b0 = _cons, b1 = _b[CONS]
c0 = _cons + _b1[1.GRP], c1 = _b[CONS] + _b[1.GRP#CONS]

OK, after any Stata estimation command, there is a vector called _b which contains all the coefficients estimated in the estimation command. It remains available until another estimation command is run (after which it is overwritten by the coefficients of the new estimation command). So any time you want to write a command that makes use of the value of a particular coefficient from the most recent estimation, _b is your friend. (There are other ways to get the coefficients, but referring to _b is by far the easiest. By the way, there is also a vector _se available that contains the standard errors of the coefficients should you need to refer to them.)

As for _b1, that was a typographical error: _b was what I intended. Sorry for the confusion there. Also, although it doesn't seem to have tripped you up, I also made the repeated error of typing CONS when I meant CONC (your variable name). In the places where I wrote, in lower case, _cons, I was referring to the constant term of the regression output. (I would not use CONC as a variable name for exactly this reason: it invites confusion. I would probably use the full word concentration instead. Matter of programming style.) So the two equations should read:

b0 = _cons, b1 = _b[CONC]
c0 = _cons + _b[1.GRP], c1 = _b[CONC] + _b[1.GRP#CONC]

Concerning you question about the confidence intervals overlapping at CONC = 0, Roman is correct. When they don't overlap, the difference is statistically significant, but the inverse is not, in general, true. If you would like a non-technical reference that explains this and also shows a "workaround," have a look at Julious SA, Using confidence intervals around individual means to assess statistical significance between two means. Pharmaceutical Statistics 2004; 3:217-222.

[Carlo Lazzaro]

Dear All:
Roman, not me, deserves all the credits for the margins code.

[Jem Lane]

Dear Roman and Carlo,

Yes, I meant Roman - sorry for the mixup with names, and thanks Carlo for acknowledging.
Thanks a lot Roman for the 'contrast' command, and for explaining why I got too many margins - I should have spotted that! Also for the references re: margins.

Clyde - thank you for explaining the notation and how Stata works, plus the reference. And once again, for all your help on this and my previous post.

I think I have what I need now - most grateful for all your time and explanations.

Best wishes

Jem (still waiting for Statalist to update my profile name...)

[Jem Lane]

Sorry, a couple of other questions have come up since I have been running the commands discussed above, and I wonder if I could trouble you again...?

My model is:
xtmixed OUT i.GRP##c.CONCEN || id: CONCEN, mle variance /// I have renamed CONC as CONCEN

Firstly, if I run

margins r.GRP, at (CONCEN=(0(1)5)) contrast

I get the following output:

Contrasts of adjusted predictions
Expression   : Linear prediction, fixed portion, predict()

1._at        : CONCEN          =           0
2._at        : CONCEN          =           1
3._at        : CONCEN          =           2
4._at        : CONCEN          =           3
5._at        : CONCEN          =           4
6._at        : CONCEN          =           5

------------------------------------------------  
             |         df        chi2     P>chi2
-------------+----------------------------------
     GRP@_at |
 (1 vs 0) 1  |          1        5.10     0.0239
 (1 vs 0) 2  |          1        6.72     0.0095
 (1 vs 0) 3  |          1        8.46     0.0036
 (1 vs 0) 4  |          1       10.23     0.0014
 (1 vs 0) 5  |          1       11.96     0.0005
 (1 vs 0) 6  |          1       13.57     0.0002
      Joint  |          2       20.37     0.0000
------------------------------------------------
--------------------------------------------------------------
             |            Delta-method
             |   Contrast   Std. Err.     [95% Conf. Interval]
-------------+------------------------------------------------
     GRP@_at |
 (1 vs 0) 1  |   18.73228   8.295914      2.472583    34.99197
 (1 vs 0) 2  |   21.36751    8.24234      5.212823     37.5222
 (1 vs 0) 3  |   24.00275   8.252991      7.827186    40.17832
 (1 vs 0) 4  |   26.63799   8.327619      10.31616    42.95982
 (1 vs 0) 5  |   29.27323   8.464532      12.68305    45.86341
 (1 vs 0) 6  |   31.90847   8.660778      14.93365    48.88328
--------------------------------------------------------------

Could I check I have understood correctly:
i. the line highlighted in red is comparing groups(GRP) 1 & 0 at CONCEN=0, and there is a significant difference?
ii. But what is the line in blue? Is this some sort of overall comparison of groups statistic?

Secondly, although I am interested in comparing groups at each concentration (CONCEN), I am also interested to know if in GRP 0 for example, there is a significant change between the value of OUT at CONCEN = 0 (i.e. baseline), and at CONCEN=3 for example (where OUT deviates the most from baseline). I have tried altering the margins command above, and have read the contrast and margins entries, but they seem to require factor variables, whereas I am treating CONCEN as continuous. I get an error message

'CONCEN' not found in list of covariates

Related to this, whilst reading about the margins and contrast commands I noticed

mcomp - the multiple comparisons adjustment, which as default is set to noadjust. Should I incorporate this somehow if I am comparing OUT at CONCEN = 0 and 1,2,3 etc?

thanks

Jem