Announcement

Collapse
No announcement yet.
X
  • Filter
  • Time
  • Show
Clear All
new posts

  • Checking assumtions for repeated measures analysis

    HI,

    I am trying to analyse repeated measures data I have for two groups of samples using Stata 12IC. The data are defined as follows: each sample (id), group (GRP = 0/1), repeated measures timepoints (CONCEN = 0,1,2,3), dependent variable (OUT). I have previously discussed in the forum using xtmixed to create a multilevel mixed effects model ("Help with interpreting xtmixed output") with the following command:

    Code:
    xtmixed OUT i.GRP##c.CONCEN || id: CONCEN, mle variance
    I have belatedly realised that I presumably need to check assumptions for this model, and have been searching for assumptions of multilevel mixed effects models. I have come across several references to these, but no concise, understandable approach. And importantly, I have not been able to find what to do if your model violates one of the assumptions, e.g. if residuals are not normally distributed. This is what I have found so far:


    Regression
    from http://www.ats.ucla.edu/stat/stata/webbooks/reg/chapter2/statareg2.htm

    - linearity: relationships between predictors and outcome variable should be linear
    - Normality - the errors should be normally distributed
    - Homoscedasticity: homogeneity of variance
    - Independence: errors assoc with one observation are not correlated with erros of any other
    - Predictor variables are measured without error
    - Model specification: model should be properly specified, including all relevant variables


    Repeated measures ANOVA
    from https://statistics.laerd.com/statistical-guides/repeated-measures-anova-statistical-guide.php

    - NORMALITY: Simliar to ANOVA, each level of the independent variable needs to be approximately normally distributed
    - SPHERICITY: Similar to homogeneity of variances in between-subjects ANOVA, it is where the variances of the differences between all combinations of related groups (levels) are equal.Tested for with Mauchly's Test. Corrections can be employed in instances where the assumption of sphericity is violated: these are lower-bound estimate, Greenhouse-Geisser correction and Huynh-Feldt correction.


    Multilevel mixed effects models
    from Stata manual

    - The overall error distribution of the linear model is assumed to be Gaussian

    from http://www.griffith.edu.au/__data/assets/pdf_file/0011/439346/Stata_mixed_intro-1.pdf
    - recommends adding "cov(uns)" to the options to allow unstructured covariance matrix (this is explained further with syntax in ats.ucla.edu - Repeated Measures Analysis)



    My questions are:
    1. What are the assumptions underlying multilevel mixed effects models?
    2. Are these the same as for two way ANOVA with repeated measures?
    3. What is a realistic way to check these assumptions (i.e. to ensure compliance with essential criteria)?
    4. What to do if essential criteria are not met? (I say essential, as presumably there is some flexibility in adherence to criteria for most statistical tests)


    Would a realistic approach be:
    i. Check covariance matrix of data in wide format with

    Code:
    correlate OUT0-OUT5, cov
    If it does not appear to be compound symmetric/exchangeable, specify an unstructured covariance as follows:

    Code:
    xtmixed OUT i.GRP##c.CONCEN || id: CONCEN, cov(uns)
    ii. predict residuals with

    Code:
    predict resid, residuals
    iii. Check normality of residuals using

    Code:
    pnorm resid
    qnorm resid
    kdensity resid, normal
    Thanks

    Jem

  • #2
    Hi Jem,

    Perhaps the first thing to do is to check the distribution of residuals at upper and lower levels. Your code "predict res" will produce only the lower level residuals. Let's do it in a systematic way. You have a random slope model, therefore we need to check three types of residual distribution i) variance at lowest level 2) variance for random slopes iii) variance for random intercepts. You can estimate and see the histogram of them:

    Code:
    xtmixed OUT i.GRP##c.CONCEN || id: CONCEN, ml
    
    /*Predict level-1 residuals*/
    
    predict res,res
    
    /*Predict residuals for random slopes & intercepts:
    Note: As consistent with Stata outputs, intercepts/constants are estimated last*/
    
    Predict rslp rint,reffects
    
    /*Now you need to generate a variable to pick one individual for the histogram*/
    
    byso id: gen indv = _n == 1
    
    
    /*Create the histograms and combine the graphs for visualisation*/
    
    hist res, normal  title("Residuals") saving(x1,replace)
    
    hist rslp if indv==1, normal title("Residuals for random slopes") saving(x2,replace)
    
    hist rint if indv==1, normal title("Residuals for random intercepts") saving(x3,replace)
    
    graph combine x1.gph x2.gph x3.gph, col(2)
    If it does not appear to be compound symmetric/exchangeable, specify an unstructured covariance as follows:
    I don't think in mixed model you can do many things with covariance structures. The definable covariance structures in "mixed" are either "identity" or "unstructured". It is worth checking the covariance structures though. If you find that "unstructured" assumption is severely wrong, you can use Generalized Estimating Equation model to exclusively define the working correlation structure, which can be implemented via Stata's "xtgee". However, GEE will have its own limitations and assumptions different than "mixed". An often taken approach, against violation of assumptions or mis-specified model, is to use sandwich estimator for robust standard error which does not rely on incorrect specification of model. You can specify robust standard error using "vce(robust) option in the command line.
    Roman

    Comment


    • #3
      There is very interesting book on the matter (West, Welch, Gatecki: Linear Mixed Models: a practical guide using statistical software. 2nd edition).

      In page 41, we read a few comments on "checking model assumptions". The authors underlined that model diagnostics shall be employed after fitting a linear mixed model, so as "to check whether distributional assumptions for the residuals are satisfied and whether the fit of the model is sensitive to unusual observations. The process of carrying out model diagnostics involves several informal and formal techiques".

      What is more, it is said that, in spite of the well-stablished post-estimation methods concerning linear regression models, "diagnostics for linear mixed models are more difficult to perform and interpret".

      In short, they emphasized "residual diagnostics" based upon conditional residuals, standardized and studentized residuals, influence diagnostics, diagnostics for random effects (using eplubs - empirical best linear unbiased predictors), histograms, quantile-quantile plots and scatter plots.

      Indeed, most Stata commands are well presented by Roman in the message above.

      You can also make a scatter plot with the standardized residuals and the predicted values:

      . predict stand, rstandard
      . predict fit, fitted
      . twoway(scatter stand fit), yline(0)

      Check if there is constant variance. Otherwise, could we do something else? Yes, think about some sort of transformation of the dependent variable.

      Apart from "qnorm" for the standardized residuals, you can also employ margins and marginsplot for the "time" variable.

      Concerning the eblups, you can estimate them under the "reffects level()" option and, aftwards, have a graphical display with boxplots, now comparing the variation of the predicted values for the random effects in different levels.


      Best,

      Marcos
      Best regards,

      Marcos

      Comment


      • #4
        Hi Roman,

        Thanks for your response, and your code. Perhaps I've missed something, but from the Stata manual it says
        residuals calculates residuals, equal to the responses minus fitted values. By default, the fitted values take into account random effects from all levels in the model
        So shouldn't the following produce residuals based on all levels, not just level 1?
        Code:
        Predict resid, residuals
        Also, won't the command
        Code:
        Predict rslp rint, reffects
        just produce the BLUPs of the random effects at these two levels, not the residuals?

        And why do you generate a variable to pick one value for each id (sample) at a particular value of CONCEN? Why not produce a histogram of all the residuals together?

        With regard to covariance structures, as far as I can see, "Unstructured" seems to offer good flexibility, as in it does not necessitate compound symmetry/exchangeable covariance structure (as described at http://www.ats.ucla.edu/stat/stata/seminars/repeated_measures/repeated_measures_analysis_stata.htm). If I understand correctly, this means your covariance structure can have basically any values in all positions - is that right?
        Thanks for the "vice(robust)" command. Is this in simple terms a way of relaxing all assumptions, with the downside that it is more difficult to achieve significance?

        Jem


        Comment


        • #5
          Hi Marcos,

          thanks for your response - the info from the book is very helpful. Sounds like the checking of assumptions is more complicated than actually producing the model!
          I have read the Stata manual's entry regarding standardised residuals, though it does not explain why to use them rather than regular residuals (i.e. observed minus predicted values).
          Could I ask why you suggest a scatter plot of the standardised residuals vs fitted values? The way I had though to do it was to use
          Code:
          twoway (scatter resid CONCEN), by (GRP) yline(0)   /// where resid = observed - predicted values
          to allow visualisation of residuals at each level of CONCEN by GRP.

          Jem

          Comment


          • #6
            Jem and Roman are using terminology differently, and that is the source of some confusion. Some people refer to the random effects as residuals at those higher levels. The two things have the same meaning. So -predict resid, residuals- calculates the residuals at the observation level--which are what Jem is calling residuals. And -predict rslp rint, reffects- calculates the residuals at the id: level, which is what Jem is calling random effects.

            The reason for generating a variable to pick one value for each id (not each id at each level of CONCEN) is that -predict, reffects- places the value of the random effects for any id in all the observations associated with that id. When working with those, one generally doesn' want that same value counted as many times as there are observations. That would lead to id's with more observations being over-represented in whatever is calculated or graphed. And even if all id's have the same number of observations, it artificially inflates the sample size to include them all, which could be important if you are applying significance tests. (That said, I don't recommend applying significance tests to asses normality in this setting--I think graphical evaluation is better at capturing the degree of departure from non-normality that is most relevant in terms of the validity of model inferences. These models are not all that fragile when normality is mildly or even moderately violated.)

            It's vce(robust), not vice(robust) [and it's an option, not a command], and it produces a vce matrix whose estimates are robust to both violations of the usual assumptions about residual distribution and to misspecification of the model. Jem, read the manual, however, to get a good sense of what this really does. It does not necessarily make it more difficult to achieve significance, though it may do so. But the estimates it produces have a somewhat different meaning from what we ordinarily think of in a residual covariance matrix.

            Comment


            • #7
              Jem, I think you have the best person here to obtain a clear explanation of your queries and I can't think of anything better clar. Thanks Clyde for the explanations and correcting my typo (vice !!!) and wrong terminology. Jem, just to clarify to you a little more regarding picking one individual: Each individual measured 6 times, but each has one intercept (within) thus one variance parameter around the individual intercept. Now, after estimating that parameter, the estimated variance will be repeated for each individual 5 more times as rows for each individual are repeated and this will inflate the histogram. That is the idea behind picking one person. For mixed models, unstructured assumption is perhaps the best you option you can choose and Clyde's statement "These models are not all that fragile when normality is mildly or even moderately violated.)" rules.
              Roman

              Comment


              • #8
                @Roman: I wasn't correcting your terminology--I was just explaining that your terminology was different from Jem's, so there was a miscommunication there. But both your terminology and Jem's are widely used and I wouldn't characterize either as being correct or incorrect to the exclusion of the other.

                Comment


                • #9
                  Thanks Clyde, I was actually referring to my use of the word "command" instead "option" which you corrected and they are sensibly different. By the way, Santa is still around I guess, have you noticed my reply post #2, Santa corrected the whole thing leaving no trace and made portion of your comment irrelevant at #6 . Thanks Santa.
                  Roman

                  Comment


                  • #10
                    Hi Jem,

                    Concerning the scatterplot of standardised residuals versus fitted values, that was according to the textbook I mentioned.

                    I gather the interpretation may be somewhat related to the Bland-Altman graphs. Let's consider the 0 value in the y line as a measure of centrality of standardised values. Concerning the x line, the range of fitted values, from 0 to the highest value.

                    In other words, if variance in the residuals is not constant, for example, if dispersion of standardised values increases on a par with the fitted values, and provided we' ve inserted in the model the "right" set of covariates in the fixed effects portion, the graph suggests a transformation of the dependent variable (log, square root, etc.) might be interesting in terms of getting a better model.


                    Best,

                    Marcos


                    Best regards,

                    Marcos

                    Comment


                    • #11
                      Dear Roman, Clyde and Marcos,

                      Many thanks for your contributions - as always, very helpful and informative. By the way, I think the "vice(robust) command" may have been my error, the "vice" due to predictive text...
                      So, to try to summarise from your input above:

                      - for multilevel mixed effects models, the assumptions are essentially the same as for multivariate regression
                      - the key things to test post-estimation are:
                      • normality of the residuals distribution, at each level of the model (as per Roman's post)
                      • homoscedasticity, i.e. homogeneity of variance - e.g. with standardised residuals vs fitted values (as per Marcos)
                      The normality assumption is fairly flexible (Clyde's post), i.e. the model is still valid as long as the residuals show some degree of normality, graphically.
                      The homoscedasticity assumption can be relaxed by specifying an unstructured covariance matrix. An alternative approach would be to use the - vce(robust) - option, though interpretation is somewhat different.

                      Finally, just going back to my original post with regards to assumptions of multilevel models vs repeated measures ANOVA - for rmANOVA, am I correct in thinking the two main assumptions are:
                      • normality of the observations (rather than residuals) at each level of the independent variable
                      • sphericity?
                      Jem

                      Comment


                      • #12
                        One more thing about heteroscedasticity--at the level of observations (the lowest level residuals), if the heteroscedasticity is attributable to the effect of a variable in your data set, you can specify that with the -residuals()- option to -mixed-. That option also permits you to model the same types of residual covariance structure as are available for the higher-level random effects using the -covariance()- option, and a few others as well.. (At least in Stata 13 you can, I don't remember if they had this in Stata 12's -xtmixed-.)

                        As for your questions about rmANOVA, the sphericity assumption is indeed part of that, in theory. Sphericity is a very strong assumption and few real-world data sets will satisfy it. That is why the Huyhn-Feldt and Greenhouse-Geiser statistics are used.

                        And concerning normality of the observations at each level of the independent variable, when there is homoscedasticity (as rmANOVA assumes), that is equivalent to, not rather than, normality of the residuals.

                        Comment


                        • #13
                          Excellent, thanks Clyde. Just to clarify, one final question regarding your point about normality and homoscedasticity if I may - if I were to run a rmANOVA, in order to check normality of observations I should use -qnorm- or a histogram with normal overlay for the observed values at each level of the independent variable. And to check for homoscedasticity, I should check the variance of the observed values at each level. The key point being that for rmANOVA I use the observations/raw data, rather than residuals (c.f. multilevel models where residuals are used), although the principle of homogeneity of variance is the same for both types of value. Is that correct?

                          Jem

                          Comment


                          • #14
                            Hi Jem,

                            Thanks for approaching such an interesting subject. Among the differences between manova and multilevel models, I modestly wish to mention a few aspects.

                            First, when you have several repeated measures under a manova, sphericity is virtually unavoidable.
                            Second, when the measures are done in different lags of time, mixed models can better take them into consideration, I mean, spaces between measures are somewhat "weighted".
                            Third, missing data will have a trend to spoil much more a listwise manova than a mixed model.

                            Best,

                            Marcos
                            Best regards,

                            Marcos

                            Comment


                            • #15
                              Hi Jem,

                              Thanks for approaching such an interesting subject. Among the differences between manova and multilevel models, I modestly wish to mention a few aspects.

                              First, when you have several repeated measures under a manova, sphericity is virtually unavoidable.
                              Second, when the measures are done in different lags of time, mixed models can better take them into consideration, I mean, spaces between measures are somewhat "weighted".
                              Third, missing data will have a trend to spoil much more a listwise manova than a mixed model.

                              Best,

                              Marcos
                              Best regards,

                              Marcos

                              Comment

                              Working...
                              X