macro
grrsim nparts kops lrepli msmt part op

mcolumn opmeans prtmeans msmt part op
mconstant nparts kops lrepli nl nomsmts
mconstant index thispart thatop

#By PGMathews, 217 Third Street, Fairport Harbor, OH 44077 440-350-0911
#Copyright © February 2000 PGMathews
#Written for Version 12 of Minitab
#Rev. 2.0, PGM, Added input of nparts, kops, and lrepli to calling statement
#and validated for V14.


#This macro is a simulator for a gage error long study. It creates a
#data set with specified numbers of parts, operators, and replicates.
#The user can modify the values of the repeatability, reproducibility, and
#part variation by editing the macro. The user must specify the number of parts,
#the number of operators, the number of trials, and three output columns for
#the corresponding measurement values, the part ID, and the
#operator ID.

#If there are n parts, k operators, and l replicates there are a total
#of nomsmts=n*k*l observations. If the analysis is by two way ANOVA the
#part variation can be recovered from: 
#     sigma(parts)=sqrt(MS(parts)/k/l)
#The reproducibility can be recovered from:
#     sigma(operators)=sqrt(MS(operators)/n/l)
#The repeatability can be recovered from:
#     sigma(error)=sqrt(MS(error))

#Example call:
#   mtb> %grrsim 10 3 2 c1 c2 c3

#Follow this macro with a two-way random effects ANOVA or GRR study

name msmt 'Msmt'
name part 'Part'
name op 'Op'

let nl=nparts*lrepli			#special variable required to use set command below
let nomsmts=nparts*kops*lrepli		#total number of measurements

set op					#ID numbers for the operators
   (1:kops)nl
end
set part				#ID numbers for the parts
   kops(1:nparts)lrepli
end

rand nparts prtmeans;			#prtmeans is column of part means
   normal 1440 60.			#part mean and variability
rand kops opmeans;			#opmeans is column of operator means
   normal 0 20.				#reproducibility
rand nomsmts msmt;
   normal 0 6.				#repeatability

do index=1:nomsmts
   let thispart=part(index)		#the part ID
   let thatop=op(index)			#the operator ID
   let msmt(index)=prtmeans(thispart)+opmeans(thatop)+msmt(index)	#part+operator+error terms
enddo

endmacro