macro
   power column;
     sign alpha;
     graph.    

#Copyright © 11 December 2000 Mathews and Malnar 
#All rights reserved
#Rev. 1.0 for Minitab V13.3 (Paul and Bek)
#Rev. 1.1 (Paul)
#Rev. 1.2 Textbox in graph (Bek)
#Mathews, Malnar & Bailey
#217 Third Street, Fairport Harbor, OH 44077
#(440)350-0911
#pmathews@apk.net

#This macro calculates the power for one-way and multi-way fixed effects
#balanced ANOVA. The user specifies parameters that characterize the experiment
#and the size of the difference to be detected between two levels.
#The input specifications can be entered from a column of the worksheet
#or by responding to prompts in the session window. (The session window
#must have commands enabled. Do this from the Editor menu when the
#session window is active.) The default significance level is alpha=0.05
#but the user can specify other values with the sign subcommand. The
#graph subcommand plots the power vs. the number of replicates.

#The order of the input parameters is:
#Runs
#Replicates
#Levels
#dfmodel
#sigma
#delta  

############### Description of the Variables #################################

#Runs is the number of runs in one replicate of the experiment. For example:
#1) For a one-way analysis with 5 treatments Runs=5.
#2) For a 3x2x5 design Runs=30.

#Replicates is the number of complete replicates of the design to be used. The
#minimum number of replicates is 1.

#Levels is the number of levels of the variable to be used in the power
#calculation. If all variables are of interest then the largest number of
#levels should be specified. The exact power will be different for variables
#that have different numbers of levels. For example, a 3x8 design will have
#more observations per level of the 3 level variable than of the 8 level
#variable so the power to detect a difference between a pair of levels of the
#3 level variable will be higher than the power to detect a difference between
#a pair of levels of the 8 level variable. The user must decide which variable
#is most important and then provide the appropriate delta. If the user wants to
#determine the power to detect a difference delta for two-factor interactions
#then the number of Levels is equal to the product of the number of levels of
#those two variables.

#dfmodel is the number of model degrees of freedom. The user must make this
#calculation himself. For example:
#1) A one-way classification problem with 5 levels will have dfmodel=5-1=4.
#2) A 3x5 factorial design with interactions will have dfmodel=2+4+(2x4)=14.
#3) A 3x2x5 factorial design with two and three factor interactions will have
#   dfmodel=2+1+4+(2x1)+(2x4)+(1x4)+(2x1x4)=29.

#Sigma is the process standard deviation or precision error. The design space
#is assumed to be homoscedastic, that is, this one sigma characterizes the errors
#throughout the entire design space.

#Delta is the smallest difference that we want to detect between two levels of
#a variable. For example, if Levels is specified to be 5, then the macro reports
#the power or probability of detecting the difference delta between 2 of these
#5 levels. 

################################################################################

#This macro should be copied into the Minitab source macros directory for
#immediate use, otherwise a change directory (cd) command will be needed.

#Example call:
#	mtb>  %power c1;
#     	subc> sign 0.01;     
#	subc> graph.		

#Copy the following lines, excluding the leading # signs, to the command line:
#%power c1;
#sign 0.01;
#graph.    

###################### Examples ####################################################

#Example: (Montgomery, p. 217) 3x3 factorial design with 4 replicates. The standard
#deviation is 25 and we want to detect a difference of 40 between two levels of a
#variable.
#Solution: We have Runs=9, Replicates=4, Levels=3, dfmodel=2+2+(2x2)=8, sigma=25,
#and delta=40. The macro returns Power=0.9225.

#Example: (Mathews, p. 119) 5 level one-way classification with 7 replicates. The
#standard deviation is 4.2 and we want to detect a difference of 10 between two
#levels of the one variable.
#Solution: We have Runs=5, Replicates=7, Levels=5, dfmodel=5-1=4, sigma=4.2, and
#delta=10. The macro returns Power=0.9279.

#Example: 3x8 factorial design with 6 replicates. The standard deviation is 80 and
#we want the power to detect a difference of 60 between two levels of the first
#variable (3 levels) and a difference of 80 between two levels of the second variable
#(8 levels). Consider two-factor interactions in the model.
#Solution: This problem requires two runs of the macro, one run for each variable.
#The first run for the first variable with Runs=24, Replicates=6, Levels=3,
#dfmodel=2+7+(2x7)=23, sigma=80, and delta=60 returns Power=0.9126.
#The second run for the second variable with Runs=24, Replicates=6, Levels=8,
#dfmodel=2+7+(2x7)=23, sigma=80, and delta=80 returns Power=0.5347.

##################### Here's the code ###############################################

mcolumn column TermCol colRuns colrep colLev coldf colsig coldelta colal 
mcolumn replic colpow reps2 
mcolumn textcol labcol labels
mconstant one Quest Runs reps N Levels dfmodel sigma delta alpha 
mconstant Fcrit lambda dferror Notalpha df1 beta Power 
mconstant N2 lambda2 beta2 pow2 K1 dferr2
mconstant yruns ylevels ydfmodel yalpha ysigma ydelta dist 
mconstant txruns txlevel txdfmo txalp txsig txdel
Default alpha=0.05

Note
Note Do you want to specify your design from the terminal or from a column?
Note
Note If from terminal a column will be created in the column specified.
Note Otherwise the column specified will be the input.
Note
Note (terminal=1, column=2)
Set TermCol;
  FILE "Terminal";
  NOBS 1;
  END.
Copy TermCol Quest
Note

If Quest=1
  Note How many runs are in one replicate?
  Set colRuns;
    File "Terminal";
    NOBS 1.
  Copy colRuns Runs
  Note

  Note How many replicates?
  Set colrep;
    File "Terminal";
    Nobs 1.
  Copy Colrep reps
  Note

  Note How many levels does the variable have?
  Set colLev;
    File "Terminal";
    Nobs 1.
  Copy colLev Levels
  Note

  Note How many model degrees of freedom are there?
  Set coldf;
    File "Terminal";
    Nobs 1.
  Copy coldf dfmodel
  Note

  Note What is the standard deviation?
  Set colsig;
    File "Terminal";
    Nobs 1.
  Copy colsig sigma

  Note
  Note What is the smallest difference that you
  Note want to detect between two levels?
  Set coldelta;
    File "Terminal";
    Nobs 1. 
  Copy coldelta delta

  
Stack colRuns colrep colLev coldf colsig coldelta column.

else

  Let Runs=column(1)		#identifies variables from column 
  Let reps=column(2)
  Let Levels=column(3)
  Let dfmodel=column(4)
  Let sigma=column(5) 
  Let delta=column(6) 

endif

################### Calculations

Let N=Runs*reps			#Total number of runs in the experiment
Let Notalpha=1-alpha		#Complement of alpha 
Let df1=Levels-1		#Numerator degrees of freedom for F distribution
Let dferror=N-1-dfmodel		#Error degrees of freedom

invcdf Notalpha Fcrit;		#Critical F value for accept/reject bound
  f df1 dferror.

Let lambda=(N/2/Levels)*(delta/sigma)**2	#F distribution noncentrality parameter

#lambda is the noncentral F distribution noncentrality parameter. It is related to
#the expected F value when Ho is false by lambda=(Levels-1)F. Note that Levels-1 is
#the numerator degrees of freedom for the F distribution. Be careful referencing
#lambda from the literature - sometimes its square or square root is called the
#noncentrality parameter.

cdf Fcrit beta;			#Finds beta value to left of Fcrit for noncentral F
  f df1 dferror lambda.

Let Power=1-beta		#Here's the power

Print N Runs reps Levels dfmodel dferror Fcrit lambda sigma delta Power


################# Graphing
If Graph
  Note
  Note Please wait.
  Note Graph in progress.

    set reps2		#sets replicate values
        (2:10 15:100/5)
    end

     Do K1=1:27    
       Let N2=Runs*reps2(K1)
       Let dferr2=N2-1-dfmodel

       invcdf Notalpha Fcrit;
         f df1 dferr2.
 
       Let lambda2=(N2/2/Levels)*(delta/sigma)**2

       cdf Fcrit beta2;
        f df1 dferr2 lambda2.

       Let Pow2=1-beta2
       Let colpow(K1)=pow2
     enddo
############################## Textbox setup!
Let labcol(1)="Runs       = "      #creates column with first portion of text line
Let labcol(2)="Levels    = "
Let labcol(3)="dfmodel = "
Let labcol(4)="Alpha     = "
Let labcol(5)="Sigma    = "
Let labcol(6)="Delta      = "

Let textcol(1)=runs		#creates column for number in textline
Let textcol(2)=levels
Let textcol(3)=dfmodel
Let textcol(4)=alpha
Let textcol(5)=sigma
Let textcol(6)=delta

Text textcol textcol; 		#changes numbers to text
  Maxwidth 80.

concatenate labcol-textcol labels	#merges two columns together

Let txruns=labels(1)		#plucks cell into variable
Let txlevel=labels(2)
Let txdfmo=labels(3)
Let txalp=labels(4)
Let txsig=labels(5)
Let txdel=labels(6)


Let ydelta=colpow(1)		#finds range of y scale
let dist=(1.0-ydelta)*0.06	#finds distance correction factor for spacing of text
Let ysigma=ydelta+dist		#y coordinates for text
Let yalpha=ysigma+dist
Let ydfmodel=yalpha+dist
Let ylevels=ydfmodel+dist
Let yruns=ylevels+dist
##############################

   Plot colpow*reps2;
   LogScale 1;
   Connect;
      Size 2;
   Title "Power vs. Replicates";
   Footnote "Mathews and Malnar December 2000";
     TSize 0.7;
     Right;
   Text 40.0 yruns txruns;
   Text 40.5 ylevels txlevel;
   Text 40.5 ydfmodel txdfmo;
   Text 40.5 yalpha txalp;
   Text 40.5 ysigma txsig;
   Text 40.5 ydelta txdel;
   Figure 0.0 1.00 0 1.00;
   Data 0.15 0.85 0.15 0.85;
   ScFrame;
   ScAnnotation;
   Tick 1 1 2 3 4 5 6 7 8 9 10 20 30 40 50 60 70 80 90 100;
     Labels "1" "2" "3" "4" "5" '' "7" '' '' "10" "20" "30" "40" "50" '' "70" '' '' "100";
     TSize 0.8;
   Tick 2;
     NMinor 1;
   Axis 1;
     Label "Number of Replicates";
     TSize .9;
   Axis 2;
     Label "Power";
     TSize .9;
   Grid 1;
     Type 1 1;
     Color 15 15;
   Grid 2;
    Type 1 1;
     Color 15 15.

endif

endmacro