Re: st: ttest or xtmelogit?

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

I wrote:

   generate double mu = 0.5 + trt * `delta' + ///
     `s' * _pi / sqrt(3) * invnormal(uniform())

It would be better to have:

   generate double mu = invlogit(trt * `delta' + ///
     sqrt(`rho' / (1 - `rho') * _pi^2 / 3) * invnormal(uniform()))

Conclusions regarding the t test with sum scores versus transformed mean 
scores are the same.

Joseph Coveney

clear *
set more off
set seed `=date("2008-03-11", "YMD")'
set seed0 `=date("2008-03-11", "YMD")'
*
capture program drop simem
program define simem, rclass
  version 10
  syntax [, Delta(real 0) Rho(real 0)]
   drop _all
   tempname p_t p_l
   set obs 12
   generate int pid = _n
   generate byte trt = mod(_n, 2)
   generate double mu = invlogit(trt * `delta' + ///
     sqrt(`rho' / (1 - `rho') * _pi^2 / 3) * invnormal(uniform()))
   assert !missing(mu) // if rho is misspecified
   generate byte den = 50
   rndbinx mu den
*
   generate double pos_prim = bnlx
   replace pos_prim = 50 - 1 / 2 / 50 if pos_prim == 50
   replace pos_prim = 1 / 2 / 50 if !pos_prim
   generate double logit_pi = logit(pos_prim / den)
   ttest bnlx, by(trt)
   scalar define `p_t' = r(p)
   ttest logit_pi, by(trt)
   scalar define `p_l' = r(p)
*
   rename bnlx count1
   generate byte count0 = den - count1
   reshape long count, i(pid) j(positive)
   drop if !count
   expand count
   xtlogit positive trt, i(pid) intpoints(30)
   return scalar rho = e(rho)
   return scalar p_r = chi2tail(1, e(chi2))
   return scalar delta = _b[trt]
   return scalar p_t = `p_t'
   return scalar p_l = `p_l'
end
*
forvalues delta = 0/2 {
   forvalues rho = 0(0.35)0.7 {
       display in smcl as text _newline(1) ///
         "delta = " %03.1f `delta' " rho = " %03.1f `rho'
       quietly simulate delta = r(delta) rho = r(rho) p_r = r(p_r) ///
         p_t = r(p_t) p_l = r(p_l), reps(300): ///
         simem , rho(`rho') delta(`delta')
       summarize delta rho
       foreach var of varlist p_* {
           generate byte pos_`var' = `var' < 0.05
       }
       summarize pos_*
   }
}
exit

It takes time to run the do-file, because of the addition of -xtlogit- in 
order to verify delta and rho.  For anyone interested but time-pressed, 
results are listed below:

Null hypothesis; independence

delta (mean SD):  0.01    0.16
 rho (mean SD):  0.00    0.01

Type I error rate (nominal 5%)

 random effects logistic regresion:  3%
 Student's t test on untransformed:  4%
t test on logistic-transformed mean:  4%


Null hypothesis; rho = 30%

delta (mean SD): -0.02    0.77
 rho (mean SD):  0.30    0.10

Type I error rate (nominal 5%)

 random effects logistic regresion:  9%
 Student's t test on untransformed:  2%
t test on logistic-transformed mean:  2%


Null hypothesis; rho = 70%

delta (mean SD):  0.08    1.71
 rho (mean SD):  0.64    0.12

Type I error rate (nominal 5%)

 random effects logistic regresion: 10%
 Student's t test on untransformed:  4%
t test on logistic-transformed mean:  5%


Delta = 1; independence

delta (mean SD):  1.0     0.18
 rho (mean SD):  0.00    0.00

Power (12 mice)

 random effects logistic regresion: 100% (300/300)
 Student's t test on untransformed: 100%
t test on logistic-transformed mean: 100%


Delta = 1; rho = 30%

delta (mean SD):  0.99    0.83
 rho (mean SD):  0.29    0.11

Power (12 mice)

 Student's t test on untransformed:  24%
t test on logistic-transformed mean:  21%


Delta = 1; rho = 70%

delta (mean SD):  1.1     1.6
 rho (mean SD):  0.64    0.12

Power (12 mice)

 Student's t test on untransformed:   9%
t test on logistic-transformed mean:   6%


Delta = 2; independence

delta (mean SD):  2.0     0.21
 rho (mean SD):  0.00    0.01

Power (12 mice)

As for delta = 1


Delta = 2; rho = 30%

delta (mean SD):  1.9     0.82
 rho (mean SD):  0.28    0.11

Power (12 mice)

 Student's t test on untransformed:  58%
t test on logistic-transformed mean:  51%


Delta = 2; rho = 70%

delta (mean SD):  2.1     1.7
 rho (mean SD):  0.61    0.13

Power (12 mice)

 Student's t test on untransformed:  23%
t test on logistic-transformed mean:  23%


*
*   For searches and help try:
*   http://www.stata.com/support/faqs/res/findit.html
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/