Section 2: Power Calculations for Two-Sample and Paired T-tests I. Two-Sample T-Test with Independent Groups II. Paired T-test: observations collected from paired data I. Two-Sample T-Test with Independent Groups Example 1: Computes the power tables found in Cohen's "Statistical Power Analysis for the Behavioral Sciences", Second Edition (1988) Results from this PROC POWER procedure match Tables 2.3.1 - 2.3.6 found on pages 28-39. The values of n that appear in the left hand column of the book's tables are group totals, so to match the book's results with what SAS provides, divide NTOTAL by 2. I. Two-sample T-test OPTIONS NOdate ps=65; TITLE1 'Chapter 2.3 Power Tables, Cohen'; TITLE2 'Tabled Values of Power * 100, rounded to 2 digits'; FOOTNOTE1 '* Power values are greater than .995'; PROC DATASETS NOlist; DELETE two_smpl_pwr; run; quit; %LET sides=2; * ENTER: sides= 1 or 2; %let alpha=.05; * ENTER: desired p-value (alpha = .01 .05 .10); ODS OUTPUT OUTPUT=two_smpl_pwr(drop=Analysis error index nulldiff Weight1 Weight2); ODS EXCLUDE OUTPUT; PROC POWER; TWOSAMPLEMEANS TEST = diff ALPHA = &alpha. SIDES = &sides. MeanDiff = .1 .2 .3 .4 .5 .6 .7 .8 1. 1.2 1.4 /* muA - muB */ STDDEV = 1 /* sigma (standard deviation) */ NTotal = 16 to 100 by 2 /* Total Number in both groups */ POWER = . /* Solve for power */ ; RUN; * adjust the data to match Cohen's results; DATA two_smpl_pwr; SET two_smpl_pwr; n_per_group = nTotal/2; * divide total sample size by 2; pwr=100*ROUND(power,.01); * round power to the nearest two-digit integer; IF power > .995 THEN pwr = .; IF info eq ' ' THEN OUTPUT ; RUN; PROC TABULATE DATA=two_smpl_pwr NOseps; CLASS Alpha MeanDiff n_per_group Sides ;* StdDev; VAR Power pwr; TABLE n_per_group, sides="Sides"*alpha="Alpha: Type I Error"* meandiff='D: effect size'*pwr=' '*sum=' '*f=3.0 / rts=12 Box='Power of T-test for HO: muA=muB' misstext='*'; RUN; PROC print DATA=two_smpl_pwr noobs; where meandiff=.5 and n_per_group=20; var Sides Alpha MeanDiff StdDev NTotal Power n_per_group ; run; This file contains: Mean Std n_per_ Sides Alpha Diff Dev NTotal Power group 2 0.05 0.5 1 40 0.338 20 Example 2 Enter two group means with the GROUPMeans = option ODS OUTPUT OUTPUT=two_smpl_pwr; TITLE1 'Two Independent Samples Power Analysis for Means'; proc power; TWOSAMPLEMEANS TEST = diff ALPHA = .05 SIDES = 2 GROUPMeans = 5 | 10 STDDEV = 9 NTotal = 20 to 60 by 10 POWER = . ; RUN; PROC PRINT DATA=two_smpl_pwr NOobs ; var sides stddev NTotal Power; run; PROC TABULATE DATA=two_smpl_pwr NOseps; CLASS Sides StdDev NTotal; VAR Power; TABLE sides, StdDev*NTotal*Power=' '*sum=' '*f=6.4 / rts=12 Box='Power'; run; Here is an example how to compute power to detect the difference in two means with a t-test running PROC MIXED Assume medium effect size of MeanDiff/sigma = 2 / SQRT(16) = .5 PROC FORMAT; VALUE $gp '1'='H0' '2' = 'H1' ; RUN; %LET nl=0; * enter null value; %LET mn=2; * enter mean diff; %LET vr=16; * enter common variance; * First, make exemplary dataset with an equal number of subjects in each group; %LET n_grp = 86; * enter sample size for each group; %LET mn_A = 1; * enter group 1 mean; %LET mn_B = 2; * enter group 2 mean; DATA _mns; group='1'; DO subj = 1 TO &n_grp.; * generate records for n_grp subjects; y_mn= &mn_A. ; OUTPUT ; END; group='2'; DO subj = 1 TO &n_grp.; * generate records for n_grp subjects; y_mn= &mn_B. ; OUTPUT ; END; run; Since power calculations are based on population parameters, enter the population _variance_ in the PARMS statement and the NoProfile option on the PROC and NOIter on the PARMS statements, as indicated below: ODS OUTPUT tests3=ts3; ODS LISTING CLOSE; PROC MIXED data=_mns NOitprint NOPROFILE ; CLASS group ; MODEL y_mn = group; PARMS ( &vr. ) / NOITER ; * enter the population variance (NOT Std Dev); run; ODS LISTING; * compute power with Noncentrality value; DATA pwr; SET ts3; alpha = .05; Noncen = NumDF*Fvalue; Fcrit = FINV(1-alpha, numdf, dendf, 0); power = 1 - PROBF(Fcrit, NumDF, DenDF, Noncen); n_total = 2*&n_grp. ; RUN; proc print DATA=pwr NOobs; run; DATA _null_; set pwr; CALL SYMPUT("pwr",LEFT(ROUND(power,.001))); RUN; %PUT &pwr ; result: mn=2, alpha=.05, n_grp=20 ---> power=.338 Note this calculation of Power with PROC MIXED is entered for illustration of its capabilities. The primary purpose is to introduce how it works with a simple design. As will be demonstrated (eventually) on another page of this web-site, PROC MIXED has the capability to compute power for a variety of very complex experimental designs, including repeated measures. II. Paired T-test: observations collected from paired data Test the mean difference of a collection of paired observations Assume you have one group of randomly chosen subjects that you want to test the efficacy of a new treatment. (You should have a control group as well but for ethical reasons you are not able to withhold the treatment to anyone in this experiment.) A process for computing power with paired measurements from multiple groups are illustrated in a section yet to appear. Define the response variables: Measure response before and after therapy y1: Value before treatment y2: Value after treatment Hypotheses: Let diff = y1 - y2 The researcher is interested if was there is a change in average value of the response over time either positive or negative, so you have a 2-sided alternative hypothesis, HA. H0: mu_D = mu_1-mu_2 = 0 vs. HA: mu_D <> 0 To compute a power analysis, evaluate: [mu1, mu2] == before and after mean response values (to compute difference) [sigma1, sigma2] == standard deviations before and after (they may be equal) corr(Y1, Y2) == correlation of the before and after measurements Suppose [mu_1, mu_2] = [150, 145], so that mu_Diff = mu_2 - mu_1 = -5. The researcher believes that the standard deviation (sigma) will be greater in week 1 prior to the treatment rather than at the end of the study in week 2. Suppose [sigma1, sigma2] = [25, 20]. The correlation between the two weights (y1 and y2) across individuals is suspected to be fairly strong, that is, %rho = 0.8 A 2x2 Variance/Covariance matrix is formed where the covariance term is: Cov(y1,y2)=(25)*(20)*(0.8) = 400 The matrix contains the variances on the diagonal: _ _ Var/Cov = | 625 400 | | 400 400 | - - The standard deviation of a difference is computed as: sigma_Diff = SQRT[Var(y1) + Var(y2) - 2*Cov(y1,y2)] sigma_Diff = SQRT[625 + 400 - 2*400] = 15 Apply a one-sample T-test to assess H0: mu_Diff = 0 assuming Diff is Normal with mu_Diff = -5 and sigma_Diff = 15 What is the power to detect this difference is not equal to 0 for sample sizes of N = 10, 11, 12, 13, 14, 15 subjects? What if the correlation is different, say rho= 0.90? 0.60? What if the standard deviations are 20% larger? TITLE1 "Beginning/Ending Measurements"; TITLE2 "Power to detect difference in paired means with T-test"; PROC DATASETS NOlist; DELETE prd_smpl_pwr; run; quit; ODS OUTPUT OUTPUT=prd_smpl_pwr; ODS EXCLUDE OUTPUT; PROC POWER; PAIREDMEANS TEST = diff SIDES = 2 ALPHA = .05 PairedMeans = 150 | 145 pairedSTDDEVs = (25 20) CORR = .8 NPAIRS = 10 11 12 13 14 15 POWER = . ; RUN; PROC PRINT DATA=prd_smpl_pwr NOobs ; var sides Alpha Corr Mean1 stddev1 Mean2 stddev2 NPairs NullDiff Power ; run; results in this output Std Std Null Sides Alpha Corr Mean1 Dev1 Mean2 Dev2 NPairs Diff Power 2 0.05 0.8 150 25 145 20 10 0 0.157 2 0.05 0.8 150 25 145 20 11 0 0.171 2 0.05 0.8 150 25 145 20 12 0 0.184 2 0.05 0.8 150 25 145 20 13 0 0.198 2 0.05 0.8 150 25 145 20 14 0 0.212 ** 2 0.05 0.8 150 25 145 20 15 0 0.225 ** see this calculation with PROC MIXED below DATA prd_smpl_pwr; SET prd_smpl_pwr; std=CATX('_|_',PUT(stddev1,4.2),PUT(stddev2,4.2)); PROC TABULATE DATA=prd_smpl_pwr NOseps; CLASS sides std alpha corr NPairs; VAR power; TABLE alpha*Npairs, corr*std*Power=' '*sum=' '*f=12.3 / rts=20 Box='Power to test paired means'; TITLE3 "alpha = &alpha., Sides = &sides. "; LABEL corr='Correlation'; run; Compute power of a paired t-test with PROC MIXED * Make exemplary dataset with the total number of paired observations ; %LET n_pr = 14; * enter number of subjects; DATA _mns; DO pair = 1 TO &n_pr.; * generate records for the number of pairs, n_pr; time = 1; y_mn= 150; OUTPUT; * assign the expected mean at time 1; time = 2; y_mn= 145; OUTPUT; * assign the expected mean at time 2; END; proc print data=_mns NOobs; run; * Power calculations are based on population parameters, so enter the population _variances_ in PARMS statement and the NoProfile and NOIter, as indicated below. Also, since the variances are different at times 1 and 2, enter the arh(1) correlation structure; ODS OUTPUT tests3=ts3; ODS LISTING CLOSE; PROC MIXED data=_mns NOitprint NOPROFILE ; CLASS time pair; MODEL y_mn = time; *enter these two statements for unequal variance model; REPEATED / subject=pair type=arh(1) rcorr r ; PARMS ( 625) ( 400 ) ( 0.8) / NOiter; * these are variances/correlations Var(y1) Var(y2) corr from the problem description; *enter these two statements to compute power for the model when variances are assumed equal; *REPEATED / subject=pair type=ar(1) rcorr; *PARMS ( 625) ( 0.8) / NOiter; run; ODS LISTING; * compute power with Noncentrality value; DATA pwr; SET ts3; alpha = .05; drop probF noncen Fcrit alpha; Noncen = NumDF*Fvalue; Fcrit = FINV(1-alpha, numdf, dendf, 0); power = 1 - PROBF(Fcrit, NumDF, DenDF, Noncen); n_pairs = &n_pr. ; RUN; proc print DATA=pwr NOobs; var effect n_pairs power; run; result: Effect n_pairs power time 14 0.21175