POWERLIB : SAS/IML Software for Computing Power in Multivariate Linear Models

The POWERLIB SAS/IML software provides convenient power calculations for a wide range of multivariate linear models with Gaussian errors. The software includes the Box, Geisser-Greenhouse, Huynh-Feldt, and uncorrected tests in the “univariate” approach to repeated measures (UNIREP), the Hotelling Lawley Trace, Pillai-Bartlett Trace, and Wilks Lambda tests in “multivariate” approach (MULTIREP), as well as a limited but useful range of mixed models. The familiar univariate linear model with Gaussian errors is an important special case. For estimated covariance, the software provides conﬁdence limits for the resulting estimated power. All power and conﬁdence limits values can be output to a SAS dataset, which can be used to easily produce plots and tables for manuscripts.


Description of POWERLIB
POWERLIB is a suite of SAS/IML (SAS Institute Inc. 2007b) modules which computes statistical power for hypothesis tests in a wide variety of univariate, multivariate, and repeated measures linear models with Gaussian errors and fixed predictors. This paper describes version 2.1 of POWERLIB. This version was developed using SAS version 9.1 (SAS Institute Inc. 2007a) and can be run on both Windows and UNIX systems. The code includes matrix names longer than 8 characters, a naming convention added in version 7; hence, the modules will not run in SAS versions 6.12 or earlier.

Available models and hypothesis tests
POWERLIB computes power for the four commonly used tests for the "univariate" approach to repeated measures: Box, Geisser-Greenhouse, Huynh-Feldt, and uncorrected. These four tests and associated methods will be collectively referred to as the UNIREP approach.
The program also computes power for the three most popular multivariate test statistics: Hotelling Lawley Trace, Pillai-Bartlett Trace, and Wilks Lambda. In this manual, these three tests and related methods are collectively referred to as the MULTIREP approach.
The UNIREP approach can be shown to be equivalent to a restricted class of linear mixed models with Gaussian errors which meet the following restrictions (Gurka, Coffey, and Muller 2007 The increasing popularity of mixed models for data analysis naturally motivates the need for reliable power analysis based on mixed models. Unfortunately, as Verbeke and Molenbergs (2000) noted, very little is known about non-null distributions in mixed models. Consequently, we know of no software that dependably calculates power for general mixed models. Use of the methods proposed by Gurka et al. (2007) within POWERLIB allows for power analysis for a restricted class of linear mixed models. Even though often in practice the required restrictions are not met, a power analysis based on this restricted class of mixed models may provide reasonable guidelines for sample size, even in more general scenarios (e.g, missing data). Thus, POWERLIB in this setting becomes a valuable tool for researchers who need to plan studies in which mixed models will be fit to the collected data.
The familiar univariate linear model with Gaussian errors is a special case of the multivariate and mixed model formulation. The program's output and syntax simplifies for computation of power for this case.
and E (N × p) containing random errors. Rows of Y , X and E correspond to independent sampling units (such as subjects), columns of Y , B and E to levels of the multivariate response (often time), and columns of X and rows of B to predictors. With r = rank(X), the methods used in this program assume N > r, row i (E) ∼ N p (0, Σ), i ∈ {1, 2, . . . , N }, all rows independent, and no missing data in Y or X.
In the context of data analysis, values in B and Σ are estimated. In contrast, in the context of power analysis B and Σ are assumed known constants, although in practice, estimates of B and Σ from a previous study are often used. The software can compute confidence limits for power to reflect the uncertainty in the estimation of Σ, and in some cases, also in B.
The corresponding general linear hypothesis (GLH) involves Θ = CBU , with C an a × q matrix of known constants defining "between subject" contrasts, and U a p×b matrix of known constants defining "within subject" contrasts. The power program requires rank(C) = a ≤ q and rank(U ) = b ≤ p. The GLH is Most often Θ 0 = 0.

Statistical theory
Detailed knowledge of the statistical theory behind power computations in multivariate linear models is not required to use this software; however, for more sophisticated users, POWERLIB provides options to choose distributional approximations. Defaults have been chosen to reflect methods the authors believe to be the best available, and should not be altered without explicit rationale. Muller, Lavange, Ramey, and Ramey (1992) gave a review of the theory behind the power methods implemented. Muller and Benignus (1992) provided a brief introduction to the most basic ideas of power in the context of toxicology, while O'Brien and Muller (1993) provided a lengthy tutorial in linear models power. Additionally, Sections 2.7 and 2.9 give several references for various distributional approximations for UNIREP and MULTIREP tests and their confidence limits.  (2007), Hedeker, Gibbons, and Waternaux (1999) and Spybrook, Raudenbush, Liu, Congdon, and Martinez (2008) also provide additional useful free software. POWERLIB currently has the following advantages over these and other products: 1. POWERLIB implements the power approximations described in Muller, Edwards, Simpson, and Taylor (2007), which has considerably better test size accuracy in some cases.
(a) PROC GLMPOWER from SAS does not compute power for multivariate Gaussian models, only univariate.
(b) PASS implements older multivariate power approximations from Muller and Barton (1989).
2. POWERLIB is easy to embed in other SAS/IML code for use in simulations.
3. POWERLIB is the only software that computes confidence limits for power in Gaussian linear models.
2. How to use POWERLIB

Execution
The POWERLIB21.IML file included in the distribution contains several modules and is the only file required to run POWERLIB. The POWER module performs all power calculations. All other modules included in the POWERLIB21.IML file, except for some independent utility modules discussed later in this paper, are called by the POWER module and are transparent to the user.
With basic knowledge of SAS/IML, POWERLIB is simple to use via the general program skeleton given in Table 1. The first two program statements are always required. They initialize IML, ask for extra symbol space, and make the power modules available for use. The RUN POWER statement executes the POWER module and is also required. Note that the directory listed in the %INCLUDE statement is the directory where POWERLIB21.IML has been copied; this most likely must be modified from this skeleton.
The software can also be run from within another module by using a CALL statement. This feature is particularly helpful when running simulations. Table 2 illustrates how to define a user module that calls the power software.

Inputs overview
User inputs to the POWERLIB modules are made through 25 global matrices listed in Table 3, grouped by application. Sections 2.4-2.11 discuss how to use these matrices in the order the groups are presented in Table 3. Because these matrices are global, use of these matrix names for reasons other than their intended purpose in POWERLIB will result in an error. When PROC IML SYMSIZE=2000; %INCLUDE "..\IML\POWERLIB21.IML"/NOSOURCE2; ...programming statements to assign global matrices that describe model and choose options...

FINISH user_defined_module;
In previous versions of POWERLIB, the current input matrix TOLERANCE was spelled as TOLERANC. We retained a matrix with this spelling in the CALL statement so that version 2.1 is compatible with programs written with previous versions of POWERLIB.
Unlike the input matrices listed in Table 3, these are not globally defined. As such, these four matrices must be listed if POWERLIB is executed with a CALL statement. As before, use of these four matrix names for reasons other than their intended purpose in POWERLIB may result in an error.
All power computations are saved to the matrix _HOLDPOWER with labels given in the vector _HOLDPOWERLBL. Section 2.8 discusses options for choices of power computations to perform as well as elements to include in _HOLDPOWER.
By default, _HOLDPOWER is sent to the output window. Section 2.8 shows how to save the output matrix to a SAS dataset, which can be used, among other things, for graphing or simulation purposes. By default, several model matrices and warnings are sent to the output window prior to printing the output matrix. Section 2.8 discusses options that control matrix printing and warning notification.
Information about any numerical difficulties and approximation accuracies is stored in the global output matrix _POWERWARN with labels given in the vector _POWERWARNLBL. These

Required matrices
Computing power for a GLH with fixed predictors requires knowing only seven variables: Σ, X, B, C, U , α and Θ 0 . The user must always specify Σ, X, B, C, via the matrices SIGMA, ESSENCEX, BETA, and C, respectively. The power program assumes default values of U = I p (which corresponds to a MANOVA hypothesis with multivariate responses), Θ 0 = 0, and α = 0.05, specified with the matrices U, THETA0, and ALPHA, respectively.
Note that for univariate models, SIGMA is the variance, not the standard deviation, because Σ = σ 2 is 1 × 1. POWERLIB gains substantial advantages from treating the univariate case as a special case of the multivariate, which requires the same variable name (SIGMA) for both settings, causing this potential source of confusion.

Input checks
If a required input has not been given by the user or if any input matrix does not conform to the dimensions specified in Table 3, POWERLIB will stop with an error. POWERLIB also checks for many other possible errors in the user input matrices to assure they contain plausible values.
Specifying the design matrix X The essence matrix (Helms 1988) contains one and only one copy of each unique row of the original matrix. We use this notation for X to ease choice of sample sizes for which to compute power. For example, the following shows the X matrix for an independent groups t test with 10 observations, denoted as X 1 , in terms of its essence matrix: X 1 = 1 10 0 10 0 10 1 10 = 1 0 0 1 ⊗ 1 10 = Es (X 1 ) ⊗ 1 10 .
In the above, 1 10 denotes a 10 × 1 vector of 1's and 0 10 a 10 × 1 vector of 0's. The ⊗ symbol denotes the Kronecker multiplication operator, where A ⊗ B = {a ij · B} for any matrices A and B. See Muller and Fetterman (2002) for an introduction to ANOVA coding in one and two way designs.
The program assumes that X is specified in terms of its essence matrix and row replication factor. These are specified with the matrices ESSENCEX and REPN, respectively. In terms of these matrices, X 1 above is specified with: The X matrix can also always be given as an essence matrix equal to the entire matrix and replication of one. In this way, X 1 can be specified with: ESSENCEX = I(2) @ J(10,1,1); REPN = 1; In the above, X 1 had equal cell sizes. The program can tolerate unequal cell sizes through no special coding. As an example, consider: One way to specify this is: Alternately, as mentioned previously, X can be specified with ESSENCEX as the entire X matrix and REPN = 1. In this representation, X 2 is coded as: Finally, POWERLIB allows X to have fractional cell sizes, e.g, REPN = DO(1 TO 5 BY .5) , by specifying the FRACREPN option. Section 2.8 discusses how to designate this and other options available with POWERLIB.
2.5. A simple power program -One power value from a two sample t test Table 4 gives the complete code needed to compute power for an independent groups t test with ten observations per group and cell mean cell coding. This program computes one power value and is an example of the simplest program that can be written to call POWERLIB. Note that for univariate models, Σ = σ 2 , a variance, not a standard deviation.

Producing power for a range of scenarios
The program makes it easy to compute power for ranges of values of sample size, type I error, mean difference, variance, and correlation level among response variables, through the global matrices REPN, ALPHA, BETASCAL, SIGSCAL, and RHOSCAL, respectively.
REPN and ALPHA contain simple lists of desired sample size and significance levels. BETASCAL and SIGSCAL contain multipliers for the user-specified values of BETA and SIGMA, respectively. BETASCAL could be called THETASCAL because it also multiplies Θ by the same amount.   more than one of these matrices has been changed from the default, power is computed for all factorial combinations of inputs.
The matrices described in this section are especially useful for producing power curves. Such graphical displays of power over a wide range of study design parameters are far more helpful than computation of point estimates. Figure 1 gives the output of Example 1 in the Examples folder of the software distribution. This program computes power over several values of mean difference (linear combination of elements of B).

Computing confidence limits for power
In the context of power analysis, B and Σ are assumed known constants. In practice, however, estimates of Σ only or of both B and Σ from a previous study are often used. The randomness of these estimates create randomness in resulting power values. Confidence intervals about power values, and confidence regions about power curves, greatly help in planning by accounting for the uncertainty due to using estimates. POWERLIB is the only current software which allows for computation of confidence intervals for power in Gaussian linear models. Muller (1995, 1996) give theory for calculation of exact confidence intervals for univariate linear models with Gaussian errors. Gribbin (2007) and Park (2007) develop approximate confidence intervals for UNIREP and MULTIREP tests, respectively. These methods are exact for univariate models and for the MULTIREP tests when s = 1. For confidence CLTYPE = 1; N_EST = 21; *# Obs for variance estimate; RANK_EST = 1; *# model df for study giving variance estimate; ALPHA_CL = 0.025; *Lower confidence limit tail size; ALPHA_CU = 0.025; *Upper confidence limit tail size; Table 6: Additional code to create confidence limits for power.
limits for multivariate tests, we assume the user has given the (unbiased) REML estimate for Σ. Note that the current version of the program does not allow compensating for the bias of truncation discussed in Taylor and Muller (1996) and Muller and Pasour (1997 (Lenth 2001).
If CLTYPE ≥ 1, RANK_EST and N_EST are required inputs. RANK_EST describes the rank of the design matrix from the previous study in which B and Σ were obtained. N_EST gives the total sample size of the previous study.
ALPHA_CL, and ALPHA_CU specify the desired lower and upper tail probabilities, respectively, for confidence limits computations. Both have default values equal to 0.025. Table 6 gives an example of the additional code needed to create confidence limits for power. For the entire context of this code, see Example 4 in the Examples folder of the software distribution. Figure 2, also produced by code in Example 4, shows an example of a useful plot for determining the impact of estimation of power parameters.

Overview
The next five sections describe the two global input matrices OPT_ON and OPT_OFF. Entries of these matrices are options which the user may specify to modify output produced by POWERLIB. These options allow the user to specify: Print power program warnings NOPRINT Suppress all printed output DS Write _HOLDPOWER to a SAS file with variable names _HOLDPOWERLBL When b = 1 (and a ≥ 1), the UNIREP and MULTIREP tests coincide to give the familiar univariate linear model test. If the COLLAPSE option is on (as is the default), then if b = 1 (and a ≥ 1) all MULTIREP and UNIREP powers coincide, and so are combined into one output column labeled POWER. If the collapse option is on with b > 1 and a = 1 then all MULTIREP variable powers are combined into one column labeled POWER_MULT.

Additional model specification options
By default: 1. X must be full rank.
The LTFR option allows choice of less than full rank X. The FRACREPN option allows choice of fractional values in REPN. To tolerate non-orthonormal U , specify the ORTHU option. Note that user utility modules UPOLY1, UPOLY2, and UPOLY3 (documented in section 2.13) may also be used to create orthonormal contrast matrices for 1-3 within or between factors. Table 8 lists all possible columns in the output matrix _HOLDPOWER. Most are not included by default. To add or remove columns, include the option of the same name from the third section of Table 7 in OPT_ON or OPT_OFF. Note that asking the program to compute confidence limits for power values adds the columns with suffixes _L and _U to the output matrix.

Options that control printing and warnings
POWERLIB prints several matrices to the screen by default after power computations are completed. The fourth section of Table 7 lists possible choices and defaults. User inspection of these matrices is important to verify that the model matrices specified to POWERLIB are those intended.
Specifying the NOPRINT option in OPT_ON ensures that no matrices are printed, including the final _HOLDPOWER output matrix. This option is especially useful when power values are output to a dataset and when many power values are computed in a simulation study. By default, POWERLIB writes helpful warnings to the screen even when no fatal syntax error is present. To suppress printing of these, specify the WARN option in OPT_OFF.

Output to a dataset
To create a dataset, the user must specify the name of the dataset in the input matrix DSNAME as well as include the DS option in the matrix OPT_ON. The user can name the data file by defining DSNAME as follows:  If DSNAME is not defined and the DS option is also selected or if "membername" already exists in the library specified by "libref", a default file name is used. The default file names are numbered and of the form WORK.PWRDT### (where ### is a number). The program scans the library for the largest numbered data file and assigns the next number to the new data file. The maximum ### is 999. If PWRDT999 exists no more data files can be created. Note that the program uses the name _PWRDTMP as an intermediate file. If this file already exists in the specified library no files can be created. To use a library other than WORK as the default, define:

DSNAME = {libref membername defaultlib};
The software will not write over existing files. To continually write to the same file with multiple runs of the power software, the user must consciously delete the existing file. To delete a file IN1.MYDATA, for example, execute the statement: prior to executing the software.

Choosing power approximations
Overview All powers are approximated by noncentral F probabilities. Approximations are used for all tests whenever s = min (a, b) > 1, and for the UNIREP tests whenever b > 1 (whether or not a = 1). All MULTIREP tests' powers (and test sizes) coincide whenever s = 1, while MULTIREP and UNIREP powers (and test sizes) all coincide if b = 1.
Naturally, the default approximation methods have been chosen, given the current state of knowledge, to be the best available; hence, most users will never need to change the approximation methods from the default. The information here is included to allow comparing to earlier versions of POWERLIB and other program output, and to allow for future developments.
The next two sections discuss the input matrices UMETHOD, UCDF, and MMETHOD used to specify the distributional approximations used in computing power for the UNIREP and MULTIREP tests, respectively. Choices for these matrices are summarized in Table 9.

Choosing UNIREP approximations
The software allows choosing one of four approximations for the distribution of the UNIREP test statistic under the alternative via the values of UCDF given in Table 9.
The default is UCDF = {2, 2, 2, 2}, with all four UNIREP tests using the Muller et al. (2007) approximation for the distribution of the test statistic under the alternative. Exact results may be achieved, at the cost of computing time, for the uncorrected (UN) and Box tests by specifying UCDF = {3, 2, 2, 3}.
UMETHOD specifies whether to use the Muller and Barton (1989) or Muller et al. (2007) approximations for E (˜ ) and E ( ). In turn, this option implies the approximate critical value used for the GG and HF tests. The default is UMETHOD = {2,2}.

Choosing MULTIREP approximations
For the MULTIREP tests, the Muller and Peterson (1984) approach requires specifying approximate degrees of freedom, which implies a critical value via the FINV() function, and an   Table 9: UNIREP and MULTIREP distributional approximation methods. approximate noncentrality. By default, the program uses two moment approximations (Rao 1951;McKeon 1974;Muller 1998) for the null distributions, which imply degrees of freedom and critical values. Optionally, older and less accurate one moment methods are also available (Pillai 1954(Pillai , 1955Pillai and Samson 1959). By default, the MULTIREP tests use the Muller and Peterson (1984) noncentrality approximations. Optionally, each MULTIREP noncentrality may be multiplied by N/[N − rank(X)] as recommended by O'Brien and Shieh (1992).
Using the O'Brien and Shieh (1992) multiplier gives slightly larger approximate powers. Especially for the Hotelling-Lawley test, the original Muller and Peterson (1984) noncentralities can be somewhat conservative in small samples.
Power approximations used for MULTIREP tests can be specified in the matrix MMETHOD. MMETHOD is a 3×1 vector whose elements correspond to one of four choices for the method used for the Hotelling Lawley Trace, Pillai-Bartlett Trace, and Wilks Lambda tests, respectively, as given in Table 9. The duplication of settings for Wilks' test is merely for programming convenience.
The default setting is MMETHOD = {4,2,2}. Use of the O'Brien and Shieh (1992) multiplier for all three multivariate tests may be chosen by setting MMETHOD = {4,4,4} or MMETHOD = 4, since MMETHOD can be specified as a 1 × 1 matrix if all the entries are the same.

F distribution probability calculations
In all but the most extreme cases, POWERLIB computes the probabilities from an F distribution necessary for power computations using the SAS supplied CDF function. If a test evaluating extremeness of conditions indicates likelihood of the CDF function failing, then POWERLIB computes F probabilities using the Tiku approximation (Kotz, Balakrishnan, and Johnson 2000). In situations where the Tiku approximations fails or will be inaccurate, POWERLIB uses a Gaussian approximation via the CDF function.
The user can see which method has been used by specifying the FMETHOD option in the OPT_ON matrix. This includes the FMETHOD columns listed in Table 8 into the output matrix _HOLDPOWER. Note that these columns merely describe which methods have been used; the user cannot modify the method chosen.
Difficulties with power calculations occur almost entirely when power approaches 0 or 1.

Numerical accuracy
The input matrices ROUND and TOLERANCE control the rounding of output values and the threshold for judging whether a numerical value is judged to be zero, respectively. The TOLERANCE matrix is included to provide a user with sophistication in computing methods with some flexibility when working with large models.

Element Description 1
A Tiku approximation was used in calculating power. 2 A Z approximation was used in calculating power. 3 A Z approximation was used in calculating power and |Z| > 6, so that the power returned is exactly 0 or 1. 4 Power is missing because the FINV function returned a missing value. 5 The lower confidence limit on power is conservative. 6 A Tiku approximation was used in calculating the lower confidence limit on power. 7 A Z approximation was used in calculating the lower confidence limit on power. 8 A Z approximation was used in calculating the lower confidence limit on power and |Z| > 6, so that the power returned is exactly 0 or 1. 9 The lower confidence limit on power is missing because the FINV function returned a missing value. 10 The upper confidence limit on power is conservative. 11 A Tiku approximation was used in calculating the upper confidence limit on power. 12 A Z approximation was used in calculating the upper confidence limit on power. 13 A Z approximation was used in calculating the upper confidence limit on power and |Z|> 6, so that the power returned was exactly 0 or 1. 14 The upper confidence limit on power is missing because the FINV function returned a missing value. 15 Power is missing because because the noncentrality could not be computed. 16 Confidence limits are missing because power is missing. 17 The approximate expected value of estimated epsilon was truncated up to 1/b. 18 The approximate expected value of estimated epsilon was truncated down to 1. 19 Power missing due to failure of Davies' algorithm. 20 Inputs give off-diagonal correlation = 1 in RHO.

21
(N -R) <= 5, so power approximations may be inaccurate, especially  Power values were rounded to 1 using the value contained in ROUND and should not be be reported as Power = 1. For example, if ROUND = 3 then report Power > 0.999.

Error checking
Counts of certain numerical difficulties are stored in entries of the output vector _POWERWARN. Table 10 describes the elements of this matrix.

User utilities
Four modules are included which are useful for creating U contrast matrices. These are: UMEAN is a function module which generates a p × 1 U matrix which computes the average response, p −1 1 1 . . . 1 . For example, it could create the matrix: U = 1/5 1/5 1/5 1/5 1/5 .
It has one input, p, which tells the size of the matrix to create. For any arbitrary value of p, the user may execute: For example, the previous matrix would be obtained with U = UMEAN(5).
UPOLY1, UPOLY2, and UPOLY3 each generate U contrast matrices with orthogonal polynomial coding for one, two, or three repeated factors, respectively, via the SAS ORPOL function.
UPOLY1 takes two inputs, VALUES and NAME. Here VALUES is a k × 1 or 1 × k vector which gives the k levels of the single repeated factor. NAME is a 1 × 1 character matrix describing the repeated factor. The module outputs two matrices U and ULBL. Columns of U contain up to level k − 1 polynomial contrasts; c contains labels for the order of the polynomial contrast each column represents. The UPOLY1 module may be called using the following syntax: CALL UPOLY1 (VALUES, NAME, U, ULBL); As an example, the following code creates a U matrix with orthogonal polynomial contrasts for four levels, 1, 10, 100, and 1000, of the factor DOSE: The columns of U (as described by ULBL) are the linear, quadratic, and cubic polynomial trends for DOSE. This produces the matrix: UPOLY2 and UPOLY3 generate U contrast matrices with orthogonal polynomial coding for two and three repeated factors, respectively. Modules UPOLY2 and UPOLY3 work the same way as the UPOLY1 module except they require 4 or 6 input matrices, respectively, and produce 6 or 12 matrices, respectively, due to adding one or two more factors. These modules are called with the following syntax: CALL UPOLY2 (VALUES1, NAME1, VALUES2, NAME2, U1, U1LBL, U2, U2LBL, U12, U12LBL); CALL UPOLY3 (VALUES1, NAME1, VALUES2, NAME2, VALUES3, NAME3, U1, U1LBL, U2, U2LBL, U3, U3LBL, U12, U12LBL, U13, U12LBL, U23, U23LBL, U123,U123LBL); Here VALUES1, VALUES2, and VALUES3 give the levels for factors 1, 2, and 3, respectively, and NAME1, NAME2, NAME3 describe factors 1, 2, and 3, respectively.
Also, U1, U2, and U3 give the main effect contrasts for factors 1, 2, and 3, respectively and U12, U13, and U23 give the two way interaction contrasts for factor 1 with 2, 1 with 3, and 2 with 3, respectively. U123 gives the three way interaction contrasts for factor 1 with 2 and 3.
U1LBL, U2LBL, U3LBL, U12LBL, U13LBL, U23LBL, and U123LBL give column labels for matrices U1, U2, U3, U12, U13, U23, and U123, respectively. Again, numbers in the column labels describe the degree of the polynomial trend corresponding to that column. Labels with one variable indicate a main effect; labels with two variable names indicate a two-way interaction; labels with three variables indicate a three-way interaction. The following statements create a U orthogonal polynomial trends matrix including all trends for main effects, two-way interactions, and the three-way interaction:

Additional examples
The code, log, output listings, and required pre-existing datasets for all the following power programs are found in the ZIP file available for download with this paper. They can be run in an interactive or a batch environment. One change is needed for the user to run these programs: the folder where the POWERLIB files and folders have been copied must be specified in the macro variable ROOT with: %LET ROOT = Your location here..; This variable is used when bringing in the POWERLIB21.IML code as in the statement: %INCLUDE "&ROOT.\Iml\POWERLIB21.IML"/NOSOURCE2; Additionally, in order that the statement will run, the POWERLIB21.IML program is assumed to reside in a sub-folder of the ROOT directory named IML. Similarly, many programs assume that necessary pre-existing datasets reside in a sub-folder of the ROOT directory named DATA. The programs that produce plots will likely require changing the FILENAME statement, as well as a few GOPTIONS, such as DEVICE, to tailor the output to the particular computer. The DISPLAY (in GOPTIONS) may also need to be controlled.

Power for a t test with overlay plot
Example 1 calculates power for a two-sample t test. The hypothesis tested is whether the two group means are equal. Power is computed for three values of σ 2 and several values of mean difference (B). Powers for these values are then plotted on a power curve.

Power for a paired t test
Example 2 performs power calculations for a simple paired t test using a general linear hypothesis in a multivariate setting. The second section of code produces results equivalent to those produced by the first section; however, it uses difference scores to test the null hypothesis of no difference between group means.

Power for a t test with 3 dimensional plot
Example 3 produces a 3-dimensional graph that illustrates power trade offs among total sample size and the hypothesized difference between two group means in an independent groups t test.

Confidence limits for a univariate model test
Example 4 produces three graphs showing two-sided or one sided lower or upper confidence limits for power, which reflect uncertainty in power calculations due to use of estimated variance parameters. The program utilizes the power confidence limit calculations available in POWERLIB and replicates the figures seen in Taylor and Muller (1995).

Power for a test of an interaction term in a multivariate model
Example 5 performs a more complicated set of power calculations for a test of the hypothesis of no time by treatment interaction in a multivariate model.

Confidence limits for a UNIREP test in a multivariate model
Example 6 utilizes a dataset that contains cerebral vessel tortuosity measures for subjects in four regions of the brain. A set of power calculations is performed for the test of the hypothesis of no gender by region interaction with five age groups. Confidence limits are computed for these power values.
Four graphical displays are produced: 1. Three 3 dimensional plots of power by sample size by mean difference, displayed on different axes 2. A plot of the hypothesized gender by region interaction with a sample size of 100 and an approximate Geisser-Greenhouse power of 0.90 3. A plot of Geisser-Greenhouse power curves for sample sizes of 20, 40, and 80 for the gender by region interaction 4. A plot with confidence limits for power with N=40 from the third plot.

Illustrate use of the UPOLY1 module
Example 7 demonstrates the use of the UPOLY1 module when performing power calculations for a time by treatment interaction.

Illustrate use of the UPOLY3 module
This first part of Example 8 demonstrates the direct creation of three way contrast matrices of two types: orthonormal polynomials and pair-wise differences to a reference level. As long as cell mean coding is used for a factorial design (including the special case of a oneway design), the approach taken to create U matrices may be applied to C matrices, and vice versa, with the obvious change of transposing matrices. Although we recommend using UPOLY3, this example is intended to provide a basis for creating contrasts for more unusual designs. The second part of Example 8 demonstrates creation of the same contrast matrices using the UPOLY3 module.

Test in a multivariate model with two within factors
Example 9 illustrates use of UPOLY2 for a design with two within and no between subject factors. It uses SIGSCAL and combines results from multiple runs of the power module. The results reproduce, except for some rounding differences, the predicted GG and HF powers in Table III in Coffey and Muller (2003), which used version 1 of POWERLIB, based on Muller and Barton (1989) methods. The example program also produces predicted Huynh-Feldt and Geisser-Greenhouse powers using the methods from Muller et al. (2007), which are the default in POWERLIB version 2.1. The new methods are far more accurate, especially for very small (near 1/b) or very large (near 1) values of .

Concluding remarks
For those familiar with SAS, POWERLIB substantially increases the power analysis capabilities related to univariate and multivariate linear models. With some simple manipulations of code in SAS/IML, one can calculate power for a wide range of tests and for a variety of approximation methods associated with them. At the same time, users who wish to use the best approximation methods available can simply use the default options in POWERLIB. POWERLIB offers cutting-edge power analysis capabilities, such as confidence intervals for power and power analysis for a special class of linear mixed models. For these reasons among others, POWERLIB is a capable tool for any SAS user who needs comprehensive sample size calculations when planning a study. Additionally, the ability to easily implement POWERLIB in simulation studies makes it a useful tool for those researching statistical power in general.