[continuation of the previous post]
Now let's turn to your PROC POWER step where you obtained a sample size of 69 subjects per group. The corresponding power estimate is even more suitable for comparisons:
ods output output=ppow;
proc power;
twosamplemeans test=equiv_diff
dist=normal
lower = -0.223
upper = 0.223
meandiff = 0.05
stddev = 0.4
npergroup = 69
power = .;
run;
proc print data=ppow;
format power best16.;
run; /* 0.80179614325271 */
/* The same result (up to 14 decimals!) can be obtained by using a formula from section 3.2.3 of the book (p. 62): */
data _null_;
a=0.05; /* significance level alpha */
dU=0.223; /* bioequivalence limit (log(1.25)) */
d=0.05; /* mean difference delta */
s=0.4; /* common standard deviation in both groups */
n=69; /* sample size per group */
power=1-cdf('T',quantile('T',1-a,2*n-2),2*n-2,(dU-d)/(s*sqrt(2/n)))
-cdf('T',quantile('T',1-a,2*n-2),2*n-2,(dU+d)/(s*sqrt(2/n)));
put power= best16.;
run;
The power computed above is that of an equivalence test based on a much simpler model than model (10.2.1) in the book: Just two independent normal distributions (ideally) with equal variances (which are assumed to be unknown). No cross-over design, just overall mean plus fixed treatment effect plus normal random error ("two-sample parallel design").
/* Simulate data of 500,000 experiments in a two-sample parallel design */
%let n=69; /* sample size per group */
%let mu=1; /* overall mean, arbitrarily selected */
%let d=0.05; /* fixed treatment effect */
%let s=0.4; /* common standard deviation in both groups */
data sim2(drop=F:);
call streaminit('MT64',31415927);
retain sim k i y e;
array F[*] FT FR (%sysevalf(&d/2) -%sysevalf(&d/2));
do sim=1 to 500000;
do k=1 to 2;
do i=1 to &n;
e=rand('normal',0,&s);
y=&mu+F[k]+e;
output;
end;
end;
end;
run;
/* Perform equivalence test */
options nonotes;
ods select none;
ods noresults;
ods output EquivLimits=eqlim2(where=(method='Pooled'));
proc ttest data=sim2 tost(-0.223, 0.223);
by sim;
class k;
var y;
run;
options notes;
ods select all;
ods results;
proc freq data=eqlim2;
tables assessment / binomial;
run; /* proportion assessment="Equivalent", i.e., estimated power: 0.8011 (95%-CI: [0.8000, 0.8022]) */
Again, the result of the simulation matches the computed power 0.8017... obtained before.
... View more
