st: Re: Mixed models for repeated measures

["Joseph Coveney" <jcoveney@bigplanet.com>]

Avi Zakai wrote:

I have data for the comparison of a new treatment for eye disease versus two
control groups.

My response (Y) is an ordinal variable with 8 levels. I have the group variable
(X), and for every subject (S), I have two measurments, representing the two
eyes !

I need to run a model to compare the different treatments, while taking into
account the inner correlation within a subject, and I need to do it with STATA
12.

I was thinking about mixed models, however no matter what I tried, I got error
messages. I would like to ask, what is the syntax I need to write in order to
run a mixed model of Y=bX+e where e~N(0,R). The correlation between the eyes is
0.6, there is no structure since it's only 2 measurments. This model assumes
normallity, it's the best I could think of, is there a better model in STATA for
my ordinal scale ?

--------------------------------------------------------------------------------

In addition to what Jay and Joerg have mentioned, you can scan through the 
do-file below to glean some ideas on how to proceed.  (Start at the 
"Begin here" comment.  The first portion just sets up a fake dataset for 
illustration.)

If you have reason to believe that there will be a left-right difference (for 
example, different doses of a topical ophthalmic drug are instilled into each
eye), then you can use the treatment-by-side interaction term in the model.  
Otherwise, I would sum the left and right eyes' scores for a single response per
patient.  This has a couple of advantages.  First, sumscores are more likely to
be more normal-like (think: Rensis Likert), albeit, as Jay mentions, with eight
categories, you might not need to worry a whole lot about that.  (Do check your
residuals.)  Second, you can use ordinary least squares linear regression. 

If you're writing clinical protocol's "Statistical Considerations" section or 
the associated statistical analysis plan, then another consideration is the 
potential for a clinic-by-treatment interaction.  Referees (whether for 
scientific journals or in government health ministries) usually want this 
addressed.  For the same audience, you'll need to think carefully about the
delta, regardless of whether you're conducting so-called noninferiority 
comparisons between your different treatments.

Joseph Coveney

version 11.2

clear *
set more off
set seed `=date("2012-10-06", "YMD")'
quietly set obs 210

generate int pid = _n
generate byte group = mod(_n, 3)
generate double u = rnormal()
forvalues side = 0/1 {
    generate double latent`side' = ///
        0.5 * (group != 0) + ///
        0   * (`side' - 1/2) + ///
        1   * u + ///
        1   * rnormal()
}
quietly reshape long latent, i(pid) j(side)
generate double manifest = normal(latent)
generate double response = 0
quietly forvalues cut = 1/7 {
    replace response =`cut' if `cut' / 8 < manifest
}
label variable response "Therapeutic Response Score"
label define Groups 0 Placebo 1 Active 2 Experimental
label values group Groups
label variable group "Treatment Group"
label define Sides 0 L 1 R
label values side Sides
label variable side Eye
tabulate response group if side == "L":Sides
sort pid side
list group pid side response in 1/6, noobs sepby(pid)

*
* Begin here
*

// Linear model
xtreg response i.group##i.side, i(pid) mle nolog
estimates store Full
xtreg response i.group i.side, i(pid) mle nolog
lrtest Full
test 2.group // "Superiority" (to placebo control treatment)
lincom 2.group - 1.group, level(95) // "Noninferiority" (to reference drug)
lincom 2.group - 1.group, level(90) // "Therapeutic Equivalence"

// Generalized linear model
quietly xi i.group*i.side
gllamm response _I*, i(pid) family(binomial) link(ologit) adapt // nolog
estimates store Full
gllamm response _Igroup_? _Iside_1, i(pid) family(binomial) link(ologit) ///
    adapt // nolog
lrtest Full
test _Igroup_2
lincom _Igroup_2 - _Igroup_1, level(95)
lincom _Igroup_2 - _Igroup_1, level(90)

// Alternative (also see -findit gologit2-)
ologit response i.group##i.side, cluster(pid) nolog
test 1.group#1.side, notest
test 2.group#1.side, accumulate
ologit response i.group i.side, cluster(pid) nolog
test 2.group
lincom 2.group - 1.group, level(95)
lincom 2.group - 1.group, level(90)

exit


*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/faqs/resources/statalist-faq/
*   http://www.ats.ucla.edu/stat/stata/