//T-test online. Including odds-ratioos, risk-ratio's, and number 
//needed to treat (NNT)
//copyright: Daan G Uitenbroek PhD

function truncate(real,places)  // truncates 'real' numbers n 'places' behind the decimal dot
{
var multi=Math.pow(10,places);
return (Math.round(real*multi)/multi)
}

function cut(real,lenght)  // changes scientific into normal notation and makes to size
{
var str=""+real;
var r="0.";
var pt=str.indexOf("e");
if (pt>0 ) 
    {
     var a=str.substring(0,1)+str.substring(2,pt+1);
     for (var index=1;index<=(str.substring(pt+2,pt+5)*1-1);index++)
        r+=0;
     var st=r+a;
      for (var index=st.length;index<=lenght;index++)
         st+="0";
     return (st.substring(0,lenght));
     };
for (var index=str.length;index<=lenght;index++)
      str+="0";
pt=str.indexOf(".");
if (pt==0) str="0"+str           //solves that newer versions of Netscape give inconsistent numeric notations
return str.substring(0,lenght);
}

function format(integer,fieldlenght) // places leading spaces
{
var str=""+integer;
var r="";
for (var index=1;index<fieldlenght-str.length+1;index++)
 r+=" ";
str=r+str;
pt=str.indexOf(".");
if (pt==0) str="0"+str           //solves that newer versions of Netscape give inconsistent numeric notations
return str;
}

function odt(len1)   // to create a onedimensional table
{
  this.lenght=len1;
   for (var i=1;i<=len1;i++)
      this[i]=0;
  return this
}

function table(len1,len2)   // to create a two dimensional table
{
   this.lenght=len1;
      for (var j=1;j<=len1;j++)
           this[j]=new odt(len2);
    return this
}

function odd(number)    //returns zero if number is not odd, one if it is
{
var isodd=0*1;
if (Math.round(number/2+0.00001)==number/2)  isodd=1*1;
return isodd
}

function c2prob(chi2,df)
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

var pval=0;

if (chi2<1000)
{
var add=0;
var mult=1;
if (odd(df)<0.5)
 {
      add=1;
      mult=Math.sqrt(2/chi2/Math.PI);
   }
var denominator=1;
for (var index=df;index>1;index-=2)
         denominator=denominator*index;
var numerator=Math.pow(chi2,((df*1+add*1)/2))*Math.exp(chi2*-1/2);
sum=1*1;
m=1*1;

for (var count=df*1+2*1;m>0.00000000001;count+=2)
   {
      m=(m*chi2)/count;
      sum=sum*1+m*1;
   }
var pval=1-sum*mult*numerator/denominator
}
return pval
}

function rprob(corr,n)
{
return c2prob((corr*corr*n),1)
}

function probv(b,s,tval,ioe,df)  //support function for tprob
{
  var a=tval/Math.sqrt(df);
  if (ioe>0.5)
     {
       if (df==1) s=0;
       var hulp= 0.5*1+(a*b*s+Math.atan(a))/Math.PI;
       return hulp
     }
  else
   var hulp=0.5*1+0.5*a*Math.sqrt(b)*s;
   return hulp
}

function currentCI()
{
this.interval="95%"
this.intzval=1.96
this.intcval=3.8416
this.intup=0.975
this.intlow=0.025
}

function CI(length)
{
this.length=length
   for (var j=1;j<=length;j++)
      this[j]=new currentCI();

this[1].interval="80%"
this[1].intzval=1.28216
this[1].intcval=1.642
this[1].intup=0.90
this[1].intlow=0.10
this[2].interval="90%"
this[2].intzval=1.64485
this[2].intcval=2.706
this[2].intup=0.95
this[2].intlow=0.05
this[3].interval="95%"
this[3].intzval=1.96
this[3].intcval=3.8416
this[3].intup=0.975
this[3].intlow=0.025
this[4].interval="99%"
this[4].intzval=2.57583
this[4].intcval=6.635
this[4].intup=0.995
this[4].intlow=0.005
}

function changeCI(prop)
{
currentoption=prop.selectedIndex;
currentoption+=1*1;
CIcurrent=allCI[currentoption];
}

//following are statements to innitialize the confidence interval
currentoption=3  //3 is for 95% CI
CIcurrent=new currentCI()
allCI=new CI(4)  //because we allow for four different Confidence Intervals

function fprob(F,df1,df2)   //4    
{    
writeln="\n";    
if (navigator.appVersion.lastIndexOf('Win') != -1)  
writeln="\r\n";    
   
if (F < 0.2 || df1 < 1 || df2 < 1)   
 document.form1.outp.value +="input invallid" ;                
         //10   
else   
{   
var a = df1%2 ? 1 : 2;   
var b = df2%2 ? 1 : 2;   
var w = (F * df1) / df2;   
var z = 1.0 / (1.0*1 + w*1);   
   
if (a == 1)   
  {if (b == 1)   
     {               //20   
      var p = Math.sqrt (w);   
      var y = 1 / Math.PI   
      var d = y * z / p;   
      var p = 2.0 * y * Math.atan(p)   
     }   
  else   
     {   
      var p = Math.sqrt (w * z);   
      var d = 0.5 * p * z / w   
     }}   
else    
   {if (b == 1)   
     {   
      var p = Math.sqrt (z);   
      var d = 0.5 * z * p;   
      var p = 1.0 - p   
     }   
   else   
     {   
      var d = z * z;   
      var p = w * z   
     }}   
   
y = 2.0 * w / z;   
   
// if (a == 1)  Tolmans version   
//   for (j = b*1 + 2*1; j <= df2; j += 2)   
//      {   
//       d *= (1.0 + a / (j - 2.0)) * z;   
//       p += d * y / (j - 1.0);   
//      }   
//   else   
//   {   
//    zk = 1.0;   
//    for (j = (df2 - 1) / 2; j; j--)   
//         zk *= z;   
//    d *= zk * df2/b;   
//    p *= zk + w * z * (zk - 1.0)/(z-1.0);   
//   }   
   
for (j = b*1 + 2*1; j <= df2; j += 2)   
   {   
    d *= (1.0 + a / (j - 2.0)) * z;   
    p = (a == 1 ? p*1 + d * y / (j - 1.0) : (p*1 + w*1) * z);  
   } //Tolmans version replaces this loop   
   
y = w * z;   
z = 2.0 / z;   
b = df2 - 2;   
   
for (i = a*1 + 2*1; i <= df1; i += 2)   
   {   
    j = i*1 + b*1;   
    d *= y * j / (i - 2.0);   
    p -= z * d / j;   
   }   
   
if (p < 0.0) p = 0.0;   
   else    
   if (p > 1.0) p = 1.0;   
    
}    
if (p>0.5) p=1-p
if (p==undefined) p=0
if (p==NaN) p=0

return (p)    
}    


function tprob(tval,df)
{
if (df>3000) df=3000;  //change this if you feel its too limmited
var s=1/1;
var c=1/1;
var f=1/1;
var ioe=Math.round(Math.round(df/2+0.1*1)-df/2+0.1*1); //takes the modules
var b=df/(df*1+tval*tval);
var fk=2*1+ioe*1;
var pval=989;
if ((df-2)<2) 
   pval=probv(b,s,tval,ioe,df);
for (index=fk;index<(df-1);index+=2)
  {
    c=c*b*(fk-1)/fk;
    s=s*1+c*1;
    if (s==f) pval=probv(b,s,tval,ioe,df);
    f=s;
    fk=fk*1+2*1;
  }
if (pval==989) pval=probv(b,s,tval,ioe,df);

return pval;
}

function tt(x1,x2,n1,n2,varianc1,varianc2,EQBox,CIBox)    // does the t-test, this procedure is always executed
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

if (EQBox.checked)
{var se=Math.sqrt(((n1-1)*varianc1+(n2-1)*varianc2)/(n1*1+n2*1-2))
 var se=se*Math.sqrt(1/n1+1/n2)
}
else 
var se=Math.sqrt((varianc1/n1)+(varianc2/n2));

var t=(x1-x2)/se;

if (EQBox.checked) {df=n1*1+n2*1-2}
else
if (n1==1) {df=5000}
else
  df=((varianc1/n1)+(varianc2/n2))*((varianc1/n1)+(varianc2/n2))/(((varianc1/n1)*(varianc1/n1)/(n1-1))+((varianc2/n2)*(varianc2/n2)/(n2-1)));

var sd=se*Math.sqrt(df);
if (n1==1) sd=se
document.form1.outp.value+="mean1 eq: "+truncate(x1,3)
if ((x1>1)|(x1<-1)) document.form1.outp.value+=" (variance="+truncate(varianc1,3)+")"
else document.form1.outp.value+= " (sd="+truncate(Math.sqrt(varianc1),3)+")"
if (n1>1) document.form1.outp.value+=" (se="+truncate(Math.sqrt(varianc1/(n1-1)),4)+")"

document.form1.outp.value+=writeln+"mean2 eq: "+truncate(x2,3)
if ((x2<1)&(x2>-1)) document.form1.outp.value+= " (sd="+truncate(Math.sqrt(varianc2),3)+")"
else document.form1.outp.value+=" (variance="+truncate(varianc2,3)+")"
if (n2>1) document.form1.outp.value+=" (se="+truncate(Math.sqrt(varianc2/(n2-1)),4)+")"

if (!((x2<1)&(x2>-1)&(x1<1)&(x1>-1))&&(n1>2)&(n2>2))
{var ftest=fprob(varianc1/varianc2,n1-1,n2-1)
document.form1.outp.value+=writeln+writeln+"Single sided probability that the two variances"
                         +writeln+"are equal (ftest): "+truncate(ftest,5)}

document.form1.outp.value+=writeln+writeln+
                           "difference between means:"+writeln+truncate((x1-x2),4)+" (sd="+truncate(sd,4)+")"+" (se="+truncate(se,4)+")"+writeln+
                           CIcurrent.interval+" CI: "+truncate((x1-x2-(se*CIcurrent.intzval)),4)+"<diff<"+truncate((x1-x2+(se*CIcurrent.intzval)),4)+" (Wald)"+writeln

//fprob(F,df1,df2)

if (CIBox.checked)
{
  if ((x2<1)&(x2>-1)&(x1<1)&(x1>-1))
  {document.form1.outp.value+=CIcurrent.interval+" CI: "
                         +truncate((x1-x2-(se*CIcurrent.intzval+0.5*(1/n1+1/n2))),4)+"<diff<"
                         +truncate((x1-x2+(se*CIcurrent.intzval+0.5*(1/n1+1/n2))),4)+" (CC-Wald)"+writeln
   document.form1.outp.value+=CIcurrent.interval+" CI: "
                         +truncate((x1-x2-(Newcombe(x1,x2,n1,n2,CIcurrent.intzval,1)*CIcurrent.intzval)),4)+"<diff<"
                         +truncate((x1-x2+(Newcombe(x1,x2,n1,n2,CIcurrent.intzval,-1)*CIcurrent.intzval)),4)+" (N-W)"+writeln}}

var tprobval=tprob(t,df)
if (n1==1) document.form1.outp.value+=writeln+"z-value of difference: "
else document.form1.outp.value+=writeln+"t-value of difference: "

document.form1.outp.value+=truncate(t,3)
if (n1>1) document.form1.outp.value+="; df-t: "+Math.round(df-0.5)
document.form1.outp.value+=writeln+"probability: "+truncate(tprobval,6)
document.form1.outp.value+=" (left tail pr: "+truncate(1-tprobval,5)+")"
document.form1.outp.value+=writeln+"doublesided p-value: "
   if (tprobval>0.5) document.form1.outp.value+=truncate(2*(1-tprobval),4)+writeln
      else  document.form1.outp.value+=truncate(2*tprobval,4)+writeln
if (EQBox.checked) document.form1.outp.value+="Equal variance assumed"+writeln
return
}  // change the ,n in 'truncate' statement if you want a more precise output

function ttest(x1,x2,n1,n2,sd1,sd2,deff1,deff2,ORRBox,PARFBox,NNTBox,FBox,CBox,EQBox,CIBox)     // THE MAIN PROCEDURE THE MAIN PROCEDURE
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

if ((deff1>0)&(deff1!=1)) {n1=n1/deff1;document.form1.outp.value+="Design effects in group 1. N1 reduced to: "+truncate(n1,1)+writeln}
if ((deff2>0)&(deff2!=1)) {n2=n2/deff2;document.form1.outp.value+="Design effects in group 2. N2 reduced to: "+truncate(n2,1)+writeln+writeln}

if ((x1>1)&(n1>1)&(sd1==0)) x1=x1/n1
if ((x2>1)&(n2>1)&(sd2==0)) x2=x2/n2

if ((x1<1)|(x2<1))   // prepares the t-test analysis
  if ((x2>1)|(x1>1))  
     document.form1.outp.value+=
         "WARNING one parameter could be a proportion and the other a number"+writeln;

if (n1<2) //work out the standardeviations in this block
{sd1=Math.sqrt(x1)
 n1=1
     document.form1.outp.value+=
         "\t*** POISSON ANALYSIS ***"+writeln+writeln};
if (n2<2) {sd2=Math.sqrt(x2)
     n2=1}
if (sd1==0) {variance1=x1*(1-x1)} else variance1=sd1*sd1;
if (sd2==0) {variance2=x2*(1-x2)} else variance2=sd2*sd2;
tt(x1,x2,n1,n2,variance1,variance2,EQBox,CIBox);

if ((x1<1)&(x2<1))   // prepares the table analysis
{
 obstable=new table(3,3);
 obstable[1][1]=x1*n1;
 obstable[1][2]=(1-x1)*n1;
 obstable[2][1]=x2*n2;
 obstable[2][2]=(1-x2)*n2;

obstable[1][3]=obstable[1][1]+obstable[1][2]*1;
obstable[2][3]=obstable[2][1]+obstable[2][2]*1;
obstable[3][2]=obstable[1][2]+obstable[2][2]*1;
obstable[3][1]=obstable[1][1]+obstable[2][1]*1;    // make marginals

if (ORRBox.checked) riskratio(obstable,CIBox)
if (ORRBox.checked) alpha(obstable,n1*1+n2*1,CIBox)
if (PARFBox.checked) atririsk(obstable,CIBox)
if (NNTBox.checked) NNT(x1,x2,variance1,variance2,n1,n2,CIBox,EQBox)

if ((FBox.checked)&(CBox.checked)) 
document.form1.outp.value+=writeln+"  FISHER EXACT & CHI SQUARE ANALYSIS"+writeln
else
{ if (FBox.checked) document.form1.outp.value+=writeln+"\tFISHER EXACT ANALYSIS"+writeln
if (CBox.checked) document.form1.outp.value+=writeln+"\tCHI SQUARE ANALYSIS"+writeln}

if ((FBox.checked)|(CBox.checked)) 
{document.form1.outp.value+=writeln+"Table used:"+writeln
document.form1.outp.value+="\t -----------------------"+writeln+
                                               "\t|"+format(Math.round(obstable[1][1]),5)+
                                               "\t|    "+format(Math.round(obstable[1][2]),5)+"\t|"+writeln+
                                               "\t|-----------------------"+writeln+
                                               "\t|"+format(Math.round(obstable[2][1]),5)+
                                               "\t|   "+format(Math.round(obstable[2][2]),5)+"\t|"+writeln+
                                               "\t -----------------------"+writeln}

if (FBox.checked) 
   exact(Math.round(obstable[1][1]),Math.round(obstable[1][2]),Math.round(obstable[2][1]),Math.round(obstable[2][2]))
if (CBox.checked) Chi2(obstable,n1*1+n2*1)
}
if (n1==1) rateratio2(x1,x2,FBox)
document.form1.outp.value+=writeln+"*******ready"+writeln+writeln
}

function rateratio2(Obs,Exp,FBox)
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

var smr=Obs/Exp

document.form1.outp.value+=writeln+"Rate Ratio (RR)= "+truncate(smr,4)+" ("+Obs+"/"+Exp+")"+writeln

if (FBox.checked)
{
var lowerbino=Ibeta(CIcurrent.intlow,Obs,Exp*1+1)
lowerbino=lowerbino/(1-lowerbino)
var upperbino=Ibeta(CIcurrent.intup,Obs*1+1,Exp)
upperbino=upperbino/(1-upperbino)
document.form1.outp.value+="The exact "+CIcurrent.interval+" confidence interval: "
                         +writeln+lowerbino+"  <RR<  "+upperbino+writeln
}

document.form1.outp.value+="The Fieller "+CIcurrent.interval+" confidence interval: "
                         +writeln
                          fieller(Obs,Exp,0)
}

function fieller(k1,k21,k22)
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

var k13=k1/k21
var k14=k21*(k21-CIcurrent.intcval)
var k15=-2*k1*(k21-(k22*CIcurrent.intcval))
var k16=k1*(k1-CIcurrent.intcval)

var k17=Math.sqrt((k15*k15)-(4*k14*k16))
document.form1.outp.value+=0.5*(-k15-k17)/k14+ " <RR< "+-0.5*(k15-k17)/k14 //+"; Chi2 "+CIcurrent.intcval
}

function Ibeta(x,p,q)
{
function invertbino(p,x,n) 
{
var minu=0.0*1 
var maxu=1.0*1 

var mean=x/n;   
for (var index=1;maxu-minu > 0.0000001;index+=1) 
  { 
   if (bino2(mean,x,n)>p)  
                 maxu=mean 
   else 
                  minu=mean 
   mean=(maxu*1 + minu*1)*0.5 
   } 
return mean
} 
return invertbino(x,p,(p*1+q*1-1))
}

function Wilson(p,n,z,ul)
{
var outcome=z*Math.sqrt(z*z+4*n*(p*(1-p)))*ul+2*n*p+z*z
outcome/=(2*(n*1+z*z))
return outcome
}

function Newcombe(p1,p2,n1,n2,z,sel)
{
var u1= Wilson(p1,n1,z,1)
var l1= Wilson(p1,n1,z,-1)
var u2= Wilson(p2,n2,z,1)
var l2= Wilson(p2,n2,z,-1)

if (sel>0) var outcome=Math.sqrt(Math.abs(l1*(1-l1)/n1)+(u2*(1-u2)/n2))
      else var outcome=Math.sqrt(Math.abs(l2*(1-l2)/n2)+(u1*(1-u1)/n1))

return outcome
}

function riskratio(table,CIBox)
{
var rratio=  (table[1][1]/(table[1][1]+table[1][2]))/(table[2][1]/(table[2][1]+table[2][2]));
document.form1.outp.value+=writeln+"Risk Ratio (p1/p2): "+truncate(rratio,2) +writeln

var varrisk=Math.sqrt((1/table[1][1])-(1/(table[1][1]+table[1][2]))+(1/table[2][1])-(1/(table[2][1]+table[2][2])))
// document.form1.outp.value+=writeln+"(1/table[1][1])= "+(1/table[1][1])   for debug
// +writeln+"(1/(table[1][1]+table[1][2]))= "+(1/(table[1][1]+table[1][2]))+writeln for debug
document.form1.outp.value+=CIcurrent.interval+" CI: "+truncate(Math.exp(Math.log(rratio)-(varrisk*CIcurrent.intzval)),4)+
                           "<R.R.<"+truncate(Math.exp(Math.log(rratio)+(varrisk*CIcurrent.intzval)),4)+writeln
+"Invers R.R (p2/p1): "  +truncate((1/rratio),2)+writeln+CIcurrent.interval+" CI: "
+truncate(Math.exp(Math.log(1/rratio)-(varrisk*CIcurrent.intzval)),4)+                          "<R.R.<"+truncate(Math.exp(Math.log(1/rratio)+(varrisk*CIcurrent.intzval)),4)+writeln 

//if (CIBox.checked)
//{document.form1.outp.value+=writeln+"Alternative C.I. for Risk Ratio:"
//+writeln
//                           +CIcurrent.interval+" CI: "+truncate(rratio-//(varrisk*CIcurrent.intzval),4)}
//                           //+"<"+truncate(rratio,4)+"<"+truncate(rratio+(varrisk*CIcurrent.intzval),4)+" //(Normal Apr)"+writeln}

document.form1.outp.value+=writeln
return
};

function atririsk(table,CIBox)
{
var rratio= (table[1][1]/(table[1][1]+table[1][2]))/(table[2][1]/(table[2][1]+table[2][2]));
rratio=1/rratio
var pc=table[2][1]/(table[2][1]+table[1][1])
var atrisk=pc*(rratio-1)/rratio

document.form1.outp.value+=writeln+"Population Attributable Risk Fraction (PARF)"+writeln+"\t((pTOT-p1)/pTOT): "+truncate(atrisk,4)+writeln
if (atrisk<0) document.form1.outp.value+="PARF is negative and meaningless"+writeln+"Please remember:"+writeln+"Mean 1 are the unexposed"+writeln+"Mean 2 are the exposed"+writeln+"and p2>p1."+writeln
   else 
   {var A=table[2][1]
   var B=table[2][2]
   var C=table[1][1]
   var D=table[1][2]
   var N=A+B+C+D
   var u=CIcurrent.intzval*(A*1+C*1)*(C*1+D*1)/((A*D)-(B*C))*Math.sqrt((A*D*(N-C)+(C*C*B))/(N*C*(A*1+C*1)*(C*1+D*1))) 
   var rrlow=(A*D-B*C)*Math.exp(u*-1)/(N*C+(A*D-B*C)*Math.exp(u*-1))
   var rrhigh=(A*D-B*C)*Math.exp(u)/(N*C+(A*D-B*C)*Math.exp(u))
   document.form1.outp.value+=CIcurrent.interval+" CI: "+truncate(rrlow,4)+"<PARF<"+truncate(rrhigh,4)+writeln
   if (CIBox.checked)
   {var varrisk=Math.sqrt((1/table[1][1])-(1/(table[1][1]+table[1][2]))+(1/table[2][1])-(1/(table[2][1]+table[2][2])))
   rrlow=Math.exp(Math.log(rratio)-(varrisk*CIcurrent.intzval))
   rrhigh=Math.exp(Math.log(rratio)+(varrisk*CIcurrent.intzval))
   var atriskl=pc*(rrlow-1)/rrlow
   var atriskh=pc*(rrhigh-1)/rrhigh
   document.form1.outp.value+=CIcurrent.interval+" CI: "+truncate(atriskl,4)+"<PARF<"+truncate(atriskh,4)+" (RRCI)"+writeln
   }}

document.form1.outp.value+=writeln
return
};

function alpha(table,n,CIBox)
{
var oratio=(table[1][1]/table[1][2])/(table[2][1]/table[2][2]); 
document.form1.outp.value+="Odds Ratio p1<->p2: "+truncate(oratio,4)  

var  varodd=0;
for (index=1;index<3;index++)    
    for (jndex=1;jndex<3;jndex++) 
        varodd=varodd+(1/table[index][jndex]);   //standard deviation of odds-ratio
// document.form1.outp.value+=writeln+index+","+jndex+": "+table[index][jndex]} for debug

var varodd1=Math.sqrt(varodd*oratio*oratio); //(old not so good method to estimate standard error) 
var varodd2=Math.sqrt(varodd*(1/oratio)*(1/oratio)); //((this one for inverse (came all from Reynolds))

varodd=Math.sqrt(varodd)
document.form1.outp.value+=" (~Std.Err. "+truncate(varodd1,4)+")"+writeln+
                           CIcurrent.interval+" CI: "+truncate(Math.exp(Math.log(oratio)-(varodd*CIcurrent.intzval)),4)+
                           "<O.R.<"+truncate(Math.exp(Math.log(oratio)+(varodd*CIcurrent.intzval)),4)+" (Wald)"+writeln

document.form1.outp.value+="inverse OR p2<->p1: "+truncate((1/oratio),4)+" (~SE="+truncate(varodd2,2)+")" +writeln+
                           CIcurrent.interval+" CI: "+truncate(Math.exp(Math.log(1/oratio)-(varodd*CIcurrent.intzval)),4)+
                           "<O.R.<"+truncate(Math.exp(Math.log(1/oratio)+(varodd*CIcurrent.intzval)),4)+" (Wald)"+writeln

if (CIBox.checked)
{
document.form1.outp.value+=writeln+"Alternative C.I. for Odds Ratio:"+writeln
                           +CIcurrent.interval+" CI: "+truncate(oratio-(varodd1*CIcurrent.intzval),4)
                           +"<"+truncate(oratio,4)+"<"+truncate(oratio+(varodd1*CIcurrent.intzval),4)+" (Normal Apr)"+writeln

document.form1.outp.value+= CIcurrent.interval+" CI: "+truncate(1/oratio-(varodd2*CIcurrent.intzval),4)
                           +"<"+truncate(1/oratio,4)+"<"+truncate(1/oratio+(varodd2*CIcurrent.intzval),4)+" (Normal Apr)"+writeln
var chi2=(table[1][1]*table[2][2]-table[1][2]*table[2][1])*(table[1][1]*table[2][2]-table[1][2]*table[2][1])
chi2=chi2*(table[3][1]+table[3][2])
chi2=chi2/table[3][1]/table[3][2]/table[1][3]/table[2][3]
chi2=Math.sqrt(chi2)
if (oratio>1) {varodd=1} else varodd=-1
document.form1.outp.value+=CIcurrent.interval+" CI: "+truncate(Math.pow(oratio,(1*1-varodd*CIcurrent.intzval/chi2)),4)
                           +"<"+truncate(oratio,4)+"<"+truncate(Math.pow(oratio,(1*1+varodd*CIcurrent.intzval/chi2)),4)+" (Chi-sq based)"+writeln
document.form1.outp.value+=CIcurrent.interval+" CI: "+truncate(Math.pow(1/oratio,(1*1+varodd*CIcurrent.intzval/chi2)),4)
                           +"<"+truncate(1/oratio,4)+"<"+truncate(Math.pow(1/oratio,(1*1-varodd*CIcurrent.intzval/chi2)),4)+" (Chi-sq based)"+writeln
}
};

function NNT(x1,x2,varianc1,varianc2,n1,n2,CIBox,EQBox)
{
var nntr=1/(x1-x2)
var sd=nntr*Math.sqrt(1-(1/nntr))

if (EQBox.checked)
{var se_difference=Math.sqrt(((n1-1)*varianc1+(n2-1)*varianc2)/(n1*1+n2*1-2))
 var se_difference=se_difference*Math.sqrt(1/n1+1/n2)  //
}
else 
var se_difference=Math.sqrt((varianc1/n1)+(varianc2/n2));

var se_schulzer=1/((x1-x2)*(x1-x2))*se_difference

document.form1.outp.value+=writeln+"Number Needed to Treat (NNT): "
                +writeln+"NNT: "+truncate(nntr,3)+"; sd: "+truncate(Math.abs(sd),3)+"; se: "+truncate(se_schulzer,3)+writeln
if (CIBox.checked)
document.form1.outp.value+=CIcurrent.interval+" CI (method Schulzer): "+truncate(nntr-(CIcurrent.intzval*se_schulzer),3)+"<NNT<"+truncate(nntr+(CIcurrent.intzval*se_schulzer),3)+writeln
document.form1.outp.value+=CIcurrent.interval+" CI (method Cook): "+truncate(1/(x1-x2+(CIcurrent.intzval*se_difference)),3)+"<NNT<"+truncate(1/(x1-x2-(CIcurrent.intzval*se_difference)),3)+" (Wald)"+writeln
if (CIBox.checked)
{
document.form1.outp.value+=CIcurrent.interval+" CI (method Cook): "
                           +truncate(1/(x1-x2+(se_difference*CIcurrent.intzval+0.5*(1/n1+1/n2))),3)+"<NNT<"
                           +truncate(1/(x1-x2-(se_difference*CIcurrent.intzval+0.5*(1/n1+1/n2))),3)+" (CC-Wald)"+writeln
document.form1.outp.value+=CIcurrent.interval+" CI (method Cook): "
                           +truncate(1/(x1-x2+(Newcombe(x1,x2,n1,n2,CIcurrent.intzval,-1)*CIcurrent.intzval)),3)+"<NNT<"
                           +truncate(1/(x1-x2-(Newcombe(x1,x2,n1,n2,CIcurrent.intzval,1)*CIcurrent.intzval)),3)+" (N-W)"+writeln
}
if (EQBox.checked) document.form1.outp.value+="Equal variance assumed"+writeln
document.form1.outp.value+=writeln;
return
}

function probval(data,n)
 {
// for (var index=1;index<4;index++) 
//  for (var jndex=1;jndex<4;jndex++)    
//    document.form1.outp.value+="i= "+index+" j= "+jndex+" data "+data[index][jndex]+writeln
// return Math.exp(factors[data[1][1]+1])
return Math.exp(factors[data[1][3]+1*1]-factors[data[1][1]+1*1]+factors[data[2][3]+1*1]-factors[data[1][2]+1*1]
     +factors[data[3][1]+1*1]-factors[data[2][2]+1*1]+factors[data[3][2]+1*1]-factors[data[2][1]+1*1]
      -factors[n+1*1]);  // +1'S are because as index [0] contains len-array the whole thing has to be shifted one place
 };

function Chi2(data,n)
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

document.form1.outp.value+="The total number of cases is: "+n+writeln+writeln;

exptable=new table(2,2)
exptable[1][1]=
   (data[1][1]+data[1][2])*(data[1][1]+data[2][1])/(n);
exptable[1][2]=
   (data[1][2]+data[1][1])*(data[1][2]+data[2][2])/(n);
exptable[2][1]=
   (data[2][1]+data[1][1])*(data[2][1]+data[2][2])/(n);
exptable[2][2]=
   (data[2][2]+data[1][2])*(data[2][2]+data[2][1])/(n);

var chi2=1*0;
var lrx=1*0;

for (var index=1;index<3;index++) 
 for (var jndex=1;jndex<3;jndex++)    
{
  chi2=chi2+(((data[index][jndex]*data[index][jndex])
              -(exptable[index][jndex]*exptable[index][jndex]))
              /exptable[index][jndex]);
 var lrx=lrx+(data[index][jndex]*(Math.log(data[index][jndex])
            -Math.log(exptable[index][jndex])));
}

lrx=lrx*2;

var Yates=((data[1][1]*data[2][2])-(data[2][1]*data[1][2]))
if (Yates<0) (Yates=Yates*-1)
Yates=Yates-0.5*n
Yates=Yates*Yates
Yates=Yates/(data[3][1]*data[2][3]*data[3][2]*data[1][3])*n;

var Mantel=chi2/n*(n-1);

document.form1.outp.value+=
 "Chi squares (all with 1 (one) degree of freedom):"+writeln+
 "Pearson's= "+truncate(chi2,3)+" (p= "+cut(c2prob(chi2,1),6)+") (GFX)"+writeln+
 "Likelihood Ratio= "+truncate(lrx,3)+" (p= "+cut(c2prob(lrx,1),6)+") (LRX)"+writeln+
 "Yate's= "+truncate(Yates,3)+" (p= "+cut(c2prob(Yates,1),6)+") (cont corr GFX)"+writeln+
 "Mantel Haenszel= "+truncate(Mantel,3)+
 " (p= "+cut(c2prob(Mantel,1),6)+")"+writeln+writeln
 +"Measures based on Chi-square:"
 +writeln+"\tPhI-sq: "
 +cut(chi2/n,7)+writeln
 +"\tPearson's correlation: "+cut(Math.sqrt(chi2/n),7)+" (p="
 +cut(rprob(Math.sqrt(chi2/n),n),6)+")"+writeln
}
//<605
function exact(n11,n12,n21,n22)
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

data=new table(3,3);
var n=0;
data[1][1]=(n11/1);
n=n*1+n11*1;
data[1][2]=(n12/1);
n=n*1+n12*1;
data[2][1]=(n21/1);
n=n*1+n21*1;
data[2][2]=(n22/1);
n=n*1+n22*1;
var smallest=Math.min(n11,n12);
smallest=Math.min(n21,smallest);
smallest=Math.min(n22,smallest);

data2=new table(3,3);
data2[1][1]=n11/1;
data2[1][2]=n12/1;
data2[2][1]=n21/1;
data2[2][2]=n22/1;

data2[1][3]=data2[1][1]+data2[1][2]*1;
data2[2][3]=data2[2][1]+data2[2][2]*1;
data2[3][2]=data2[1][2]+data2[2][2]*1;
data2[3][1]=data2[1][1]+data2[2][1]*1;    // make marginals

document.form1.outp.value+=writeln+"The total number of cases is: "+n+writeln;

document.form1.outp.value+="The smallest value: "+smallest+writeln;

for (var index=1;smallest<data[1][1];index++)
  {
    var number=data[1][1];
    data[1][1]=data[1][2];
    data[1][2]=data[2][2];
    data[2][2]=data[2][1];
    data[2][1]=number;
  };    // rotate table to put smallest number in upper left corner

data[1][3]=data[1][1]+data[1][2]*1;
data[2][3]=data[2][1]+data[2][2]*1;
data[3][2]=data[1][2]+data[2][2]*1;
data[3][1]=data[1][1]+data[2][1]*1;    // make marginals

var r=data[3][1];
if (data[3][1]>data[1][3]) r=data[1][3];

document.form1.outp.value+="The smallest marginal: "+r+writeln+writeln;

factors=new odt(n*1+1*1);   
factors[1]=0;
factors[2]=0;
for (var i=1;i<n+2;i++)
  factors[i+1]=factors[i]+Math.log(i); //make an array of log factorials

pexact=probval(data,n);
plarger=0;
pdblside=0;  
pdbls2=0//initialize

//document.form1.outp.value+="Calculated pexact and innitialized "+writeln;

for (var index=0;index<smallest+1;index++) 
  {
    data[1][1]=index;
    data[1][2]=data[1][3]-index;
    data[2][1]=data[3][1]-index;
    data[2][2]=data[3][2]-data[1][2];
    plarger=plarger*1+probval(data,n)*1;
    if (probval(data,n)<(pexact+1.0e-15*1))  
      pdblside=pdblside*1+probval(data,n)*1;
    if (probval(data,n)<(pexact))  
      pdbls2=pdbls2*1+probval(data,n)*1;
  };

//document.form1.outp.value+="Did first loop "+writeln;

psmaller=pexact;

for (var index=smallest+1;index<r+1;index++)
  {
    data[1][1]=index;
    data[1][2]=data[1][3]-index;
    data[2][1]=data[3][1]-index;
    data[2][2]=data[3][2]-data[1][2];
    psmaller=psmaller+probval(data,n)*1;
    if ((probval(data,n)<(pexact+1.0e-15*1)))  
      pdblside=pdblside+probval(data,n)*1;
    if (probval(data,n)<(pexact))  
      pdbls2=pdbls2*1+probval(data,n)*1;
  };

//document.form1.outp.value+="Did second loop "+writeln;

document.form1.outp.value+="The p for exactly this table: "+cut(pexact,8)+writeln+writeln

document.form1.outp.value+="The p-value for the same or a stronger"+writeln+" association= "+cut(plarger,8)+writeln
document.form1.outp.value+="The p-value for a stronger association= "+cut(plarger-pexact,8)+writeln
document.form1.outp.value+="The mid p-value "+cut(plarger-pexact/2,8)+writeln+writeln

document.form1.outp.value+="The p-value for the same or the reverse"+writeln+" association=  "+cut(psmaller,8)+writeln+writeln

document.form1.outp.value+="Two sided p-values for p(O>=E|O<=E) "+writeln+"  p-value= "
                         +cut(pdblside,12)+" (the sum of small p's)"+writeln

if (psmaller<plarger )
             document.form1.outp.value+="  p-value= "+cut((2*psmaller),8)
    else  document.form1.outp.value+="  p-value= "+cut((2*plarger),8)
document.form1.outp.value+=" (double the single sided p)"+writeln
document.form1.outp.value+="Two sided p-value for p(O>E|O<E) "+writeln+"  p-value= "
                         +cut(pdbls2,12)+" (the sum of small p's)"+writeln
document.form1.outp.value+="Two sided p-value mid-p"+writeln+"  p-value= "
                         +cut(pdbls2/2+pdblside/2,12)+" (the sum of small p's)"+writeln+writeln
} //end exact procedure
  
function invertChi2(p,df) 
{   
var minchisq=0.0*1 
var maxchisq=2000*1 

var chisqval= (df/2) / Math.sqrt (p); 
for (var index=1;maxchisq-minchisq > 0.0000001;index+=1) 
  { 
   if (c2prob(chisqval,df)<p)  
                 maxchisq=chisqval 
   else 
                  minchisq=chisqval 
   chisqval=(maxchisq + minchisq)*0.5 
   } 
return chisqval
} 

function rateratio3(Obs,Exp,FBox,CIBox)
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

var normdev=Math.abs((Obs-Exp)/Math.sqrt(Exp))

var smr=Obs/Exp
var varsmr=Math.sqrt(1/Obs)

var logrr=Math.log(smr)
var selogrr=Math.sqrt((1/Obs)+(1/Exp))

document.form1.outp.value+=writeln+"Obs= "+Obs+"; Exp= "+Exp+writeln+"Rate Ratio (RR)= "+truncate(smr,4)+" ("+Obs+"/"+Exp+")"+writeln+"Std Dev: "+truncate(varsmr,4)+writeln
var pplow=Obs*Math.pow((1-(1/(9*Obs))-CIcurrent.intzval/Math.sqrt(9*Obs)),3)/Exp
var pphigh=(Obs*1+1*1)*Math.pow((1-(1/(9*(Obs*1+1*1)))+(CIcurrent.intzval/Math.sqrt(9*(Obs*1+1*1)))),3)/Exp
document.form1.outp.value+=CIcurrent.interval+" Poisson process CI: "+truncate(pplow,4)
      +" <RR< "+truncate(pphigh,4)+writeln

if (FBox.checked)
{document.form1.outp.value+=writeln+"Exact "+CIcurrent.interval+" Poisson confidence interval: "+writeln
                        +(0.5*invertChi2(CIcurrent.intup,(Obs*2))/Exp)+" <RR< "
                        +(0.5*invertChi2(CIcurrent.intlow,(Obs*2+2))/Exp)+writeln
if (Obs>120) document.form1.outp.value+="Too many cases for reliable estimate"+writeln}

if (CIBox.checked)
{
document.form1.outp.value+=writeln+"Often used approximations:"+writeln
      +CIcurrent.interval+" Walds C.I. "+truncate(Math.exp(Math.log(smr)-(varsmr*CIcurrent.intzval)),4)
      +" <RR< "+truncate(Math.exp(Math.log(smr)+(varsmr*CIcurrent.intzval)),4)+writeln


      +CIcurrent.interval+" Logarithmic Transformation: "+truncate(Math.exp(logrr-(selogrr*CIcurrent.intzval)),4)
      +" <RR< "+truncate(Math.exp(logrr+(selogrr*CIcurrent.intzval)),4)+writeln


      +CIcurrent.interval+" Square-root Transformation: "
      +truncate((Math.sqrt(Obs)-0.5*CIcurrent.intzval)*(Math.sqrt(Obs)-0.5*CIcurrent.intzval)/Exp,4)
      +" <RR< "
      +truncate((Math.sqrt(Obs)+0.5*CIcurrent.intzval)*(Math.sqrt(Obs)+0.5*CIcurrent.intzval)/Exp,4)
      +writeln

      +CIcurrent.interval+" Normal Deviate: "
      +truncate((Obs-Math.sqrt(Obs)*CIcurrent.intzval)/Exp,4)+" <RR< "
      +truncate((Obs*1+Math.sqrt(Obs)*CIcurrent.intzval)/Exp,4)+writeln

 //     +"Normal Dev (meth Rothman)= "+truncate(normdev,2)
 //     +"; probv: "+truncate(tprob(normdev,4000),5)+writeln
 //     +"log ND (conservative)= "+truncate((logrr/selogrr),4)
 //     +"; pr. LND: "+truncate(tprob((logrr/selogrr),4000),5)+writeln
}
}

function popul(x1,x2,n1,n2,sd1,sd2,deff1,deff2,ORRBox,FBox,CIBox)     // population analysis, complete procedure
{
writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n";

var type=3
if ((sd1==0)&(sd2==0)) type=2
if ((n1==0)&(n2==0)) type=1


if ((x1>1)&(n1>1)&(sd1==0)) x1=x1/n1
if ((x2>1)&(n1>1)&(sd1==0)) x2=x2/n1

if ((x1<1)|(x2<1))   // some error trapping
  if ((x2>1)|(x1>1))  
     document.form1.outp.value+=
         "WARNING one parameter could be a proportion the other a number"+writeln;
if ((n2>1)&(n1>1))  
     document.form1.outp.value+=
         "WARNING two numbers of cases given. Top number ("+n1+") used."+writeln;

if (type==1) //work out the standardeviations in this block
{
document.form1.outp.value+="\t*** POISSON ANALYSIS ***"+writeln+writeln;
n1=1
sd1=Math.sqrt(x1)
}

if (type==2) //work out the standardeviations in this block
{
document.form1.outp.value+="\t*** BINOMIAL ANALYSIS ***"+writeln+writeln;
}

if ((sd1==0)&(x1<1)) {sd1=x1*(1-x1)} else sd1=sd1*sd1;
sd2=sd2*sd2;
if (sd2>0) sd1=0

document.form1.outp.value+="mean1 eq: "+truncate(x1,5)+"; "+"mean2 eq: "+truncate(x2,5)+writeln

if (n1==0) (n1=n2)

if ((deff1>0)&(deff1!=1)) {n1=n1/deff1;document.form1.outp.value+="Considering design effects. N reduced to: "+truncate(n1,1)+writeln}

if (sd1!=0)
{
var dif=x1-x2
var sedif=Math.sqrt(sd1/n1)
var zval=dif/sedif
  document.form1.outp.value+="difference equals: "+truncate(dif,5)+writeln
          +"standard deviation of difference: "+truncate(Math.sqrt(sd1),3)+writeln
          +"standard error of difference: "+truncate(sedif,3)+writeln
          +CIcurrent.interval+" CI: "+truncate(dif-(CIcurrent.intzval*sedif),3)+" <dif> "+truncate(dif+(CIcurrent.intzval*sedif),3)+writeln
          +writeln+"normal deviate (z-value): "+truncate(zval,4)+writeln
tprobval=tprob(zval,3000)
document.form1.outp.value+="prob-z: "+truncate(tprobval,5)

document.form1.outp.value+=" (left tail pr: "+truncate(1-tprobval,5)+")"
document.form1.outp.value+=writeln+"doublesided p-value: "
   if (tprobval>0.5) document.form1.outp.value+=truncate(2*(1-tprobval),5)+writeln
      else  document.form1.outp.value+=truncate(2*tprobval,5)+writeln
}

if ((type==3)&(sd1==0))
{var dif=x1-x2
sedif=Math.sqrt(sd2/n1)
var tval=dif/sedif
  document.form1.outp.value+=writeln+"Population standard deviation estimated using sample"+writeln
             +"T-distribution used"+writeln+writeln+"difference equals: "+truncate(dif,3)+writeln
             +"standard deviation of difference: "+truncate(Math.sqrt(sd2),3)+writeln
             +"standard error of difference: "+truncate(sedif,3)+writeln
             +CIcurrent.interval+" CI: "+truncate(dif-(CIcurrent.intzval*sedif),3)+" <dif> "+truncate(dif+(CIcurrent.intzval*sedif),3)+writeln+writeln
             +"t-value equals: "+truncate(tval,4)+writeln
probva=tprob(tval,n1-1)
document.form1.outp.value+="df-t: "+Math.round(n1-1)+"; probability: "+truncate(probva,5)+writeln

document.form1.outp.value+=" (left tail pr: "+truncate(1-probva,5)+")"+writeln
document.form1.outp.value+="doublesided p-value: "
   if (probva>0.5) document.form1.outp.value+=truncate(2*(1-probva),5)+writeln
      else  document.form1.outp.value+=truncate(2*probva,5)+writeln
}
if (type==1)
 {if (FBox.checked)
     {var probva=pois2(x1,x2-1)
      if (x2<0.5) probva=1
      var probva2=pois2(x1,x2)
      var pp=pois3(x1,x2)
      document.form1.outp.value+=writeln+"Exact Poisson probabilities"
      document.form1.outp.value+=writeln+"Pointprob: "+cut(pp,10)
      document.form1.outp.value+=writeln+"Single sided (>=): "+cut(probva,8)+"; left (<): "+cut(1-probva,8)
      document.form1.outp.value+=writeln+"Single sided (>): "+cut(probva2,8)+"; left (=<): "+cut(1-probva2,8)
      document.form1.outp.value+=writeln+"Mid-p: "+cut((probva+probva2)/2,8)+"; left: "+cut(1-(probva+probva2)/2,7)+"; double: "
      if ((probva+probva2)/2<0.5) {document.form1.outp.value+=cut(2*((probva+probva2)/2),7)+writeln}
        else  document.form1.outp.value+=truncate(2*(1-(probva+probva2)/2),6)+writeln
     }
   x2=Math.round(x2)
   if (ORRBox.checked) rateratio3(x2,x1,FBox,CIBox)
  }

if ((type==2)&(FBox.checked))
 {var probva=bino2(x1,x2*n1,n1)
  if (x2==0) probva=1
  var probva2=bino2(x1,x2*n1+1,n1)
  var pp=bino3(x1,x2*n1,n1)
  document.form1.outp.value+=writeln+"Exact Binomial probabilities"
  document.form1.outp.value+=writeln+"Pointprob: "+cut(pp,10)
  document.form1.outp.value+=writeln+"Single sided (>=): "+cut(probva,8)+"; left (<): "+cut(1-probva,8)
  document.form1.outp.value+=writeln+"Single sided (>): "+cut(probva2,8)+"; left (<=): "+cut(1-probva2,8)
  document.form1.outp.value+=writeln+"Mid-p: "+cut((probva+probva2)/2,8)+"; left: "+cut(1-(probva+probva2)/2,7)+"; double: "
  if ((probva+probva2)/2<0.5) {document.form1.outp.value+=cut(2*((probva+probva2)/2),7)+writeln}
      else  document.form1.outp.value+=truncate(2*(1-(probva+probva2)/2),6)+writeln

  var le=bino5(x1,x2*n1,n1,pp)
  document.form1.outp.value+=writeln+"Double sided (>=): "+cut(le,10)
  var lt=bino4(x1,x2*n1,n1,pp)
  document.form1.outp.value+=writeln+"Double sided (>): "+cut(lt,10)
  document.form1.outp.value+=writeln+"Mid-p: "+cut((lt+le)/2,10)
  document.form1.outp.value+=writeln+"Double sided according to the sum of small p's"+writeln
 }
}

function pois2(u,x)
{
p0=Math.exp(-u)
var tot=1-p0;
for (var index=1;index<=x;index++) 
{    p0*=u/index; 
     tot-=p0; 
}
return tot
}

function pois3(u,x)
{
p0=Math.exp(-u)
for (var index=1;index<=x;index++) 
   p0*=u/index; 
return p0
}

function bino2(m,k,nk)
{
var p0=Math.exp(Math.log(1-m)*nk)    //(1/m)^nk
var tot=1-p0;
for (var index=0;index<(k-1);index++) 
     {p0*=((nk-index)/(index*1+1*1))*(m/(1-m));
     tot-=p0 
     }
return tot;
}

function bino3(m,k,nk)
{
var p0=Math.exp(Math.log(1-m)*nk)    //(1/m)^nk
for (var index=0;index<k;index++) 
     p0*=((nk-index)/(index*1+1*1))*(m/(1-m));
return p0;
}

function bino4(m,k,nk,pointprob)
{
var p0=Math.exp(Math.log(1-m)*nk)    //(1/m)^nk
tot=0
for (var index=0;index<=nk;index++) 
     {if (p0<pointprob) tot+=p0 
      p0*=((nk-index)/(index*1+1*1))*(m/(1-m));
     }
return tot;
}

function bino5(m,k,nk,pointprob)
{
var p0=Math.exp(Math.log(1-m)*nk)    //(1/m)^nk
tot=0
for (var index=0;index<=nk;index++) 
     {if (p0<=pointprob) tot+=p0 
      p0*=((nk-index)/(index*1+1*1))*(m/(1-m));
     }
return tot;
}

function ClearOutp()
{
var writeln="\n";
if (navigator.appVersion.lastIndexOf('Win') != -1) writeln="\r\n"; 
document.form1.outp.value="";
}