options ps=65 ls=65;
*This is powert2.pilot.sas.  Two files are read by this
program.  The first is read into correl and contains a
group indicator variable and scores on p variables.  Both
the vector of effect sizes and the correlation matrix are
calculated from these data.  The second is read into param and
contains the minimum sample size to be considered for the study,
the maximum the study, alpha, and the target power.;
data correl;
input group var1-var3;
datalines;
1 1  6 3
1 2  7 3
1 3  7 1
1 9  4 4
1 6  6 5
1 6  4 3
1 3  9 6
1 4  9 7
1 5  4 3
1 4 10 8
2 1  8 4
2 0  6 0
2 2  8 3
2 2  3 0
2 6  7 5
2 1 10 4
2 3  9 3
2 1 11 6
2 5  8 3
2 0  7 2
proc sort; by group;
proc corr cov data=correl outp=rs; by group;
  var var1-var3;
data param;
input minn maxn alpha tpower;
datalines;
2 200 .05 .80
proc iml;
edit rs;
read all into temp;
p=ncol(temp)-1;
cov1=temp[1:p,2:p+1];
cov2=temp[2*p+4:3*p+3,2:p+1];
n1=temp[p+3,2];
n2=temp[3*p+6,2];
cov=((n1-1)*cov1+(n2-1)*cov2)/(n1 + n2 - 2);
isd=(inv(sqrt(diag(cov))));
r=isd*cov*isd;
barx1=temp[p+1,2:p+1];
barx2=temp[3*p+4,2:p+1];
d=(barx1-barx2)*isd;
print barx1, barx2, cov1, cov2, cov, d, r;
edit param;
read all into temp1;
minn=temp1[1,1];
maxn=temp1[1,2];
prob=1-temp1[1,3];
target = temp1[1,4];
if minn <= p+1 then minn = p+2;
do n=minn to maxn;
  DF2 = N + N - P -1;
  nc=(n/2)*d*inv(r)*d`;
  cval=FINV(prob,P,DF2);
  POWER = 1 - PROBF(cval,P,DF2,NC);
  if n=maxn then do; if  power < target then do;
  print 'Sample too small, raise maximum size';
  print n nc cval power;
end;
end;
if power >  target then do;
  print 'Power calculations have been completed';
  print n nc cval power;
  print 'Note.  N is the reqired number of Ss in each treatment';
  n=maxn;
end;
end;
quit;