 |
[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]
Re: st: ttest or xtmelogit?
.
I'm running now at lower values of rho, which are more typical of my
data (.05 < rho < .15).
I agree simulating at larger rho's is smart because estimating rho
accurately is not easy with small samples, so I should look at the
upper bound of the estimated CI of rho.
This is fascinating stuff! I wonder if reviewers will buy this kind of
supporting documentation for analysis choice down the road.
Thanks Joseph.
-Dave
On Mar 11, 2008, at 9:25 AM, Joseph Coveney wrote:
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/
--
David C. Airey, Ph.D.
Pharmacology Research Assistant Professor
Center for Human Genetics Research Member
Department of Pharmacology
School of Medicine
Vanderbilt University
Rm 8158A Bldg MR3
465 21st Avenue South
Nashville, TN 37232-8548
TEL (615) 936-1510
FAX (615) 936-3747
EMAIL [email protected]
URL http://people.vanderbilt.edu/~david.c.airey/dca_cv.pdf
URL http://www.vanderbilt.edu/pharmacology
*
* 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/