ANOVA

[salvomic]

ANOVA two way (blocks and treatments)

ANOVA with blocks and treatments (2 way).
Values in a matrix (anovamat). Edit it inside or outside the program.

Usage:
anova2(), edit matrix and input alpha significance level (i.e. 0.05)
anova2_help() give short help.

The program returns now also a matrix with residues and pooled stddev.

Code:


export anovamat:=[[0,0],[0,0]];
export anova_res:=[0,0];
export anova_test_row_col:={0,0};
export anova_F_crit:=0;
export anova_pvalue2:={0,0};
export anova_dfreedom2:={0,0,0,0};
export anova_pooled_stddev:=0;
export anova_ssr_ssc_sse:={0,0,0};
export anova_mst_msb_mse:={0,0,0};
export anova_Xio_Xoj_Xoo:={0,0,0};
export anova_var_rowcol:={0,0};
export anova_xio_xoj_xij:={0,0,0};

ANOVA_help();

EXPORT ANOVA2()
// ANOVA 2 way a=α significance level (ex. 0.05)
// by Salvo Micciché 2015
BEGIN
local a, j, k, n, m, N;
local tot, fr, fc, mesg, re;
local  xio:={}, xoj:={}, xij, p1, p2;
local Xio:={}, Xoj:={}, Xoo, sp;
local SSE, SSR,SSC, STR, STC, MST, MSB, MSE;
local varrow:={}, varcol:={};

EDITMAT(anovamat);
sto(anovamat, anovamat);
// a:= 0.05;
INPUT(a, "Significance Level", "Sig. Level α=", "Input α for the test (default value: 0.05)", 0.01, 0.05);
IF (a <=0 OR a>= 1) THEN RETURN("Significance level α must be > 0 AND <1"); END;

m:= rowDim(anovamat);
n:=colDim(anovamat);
N:= (n-1)*(m-1);
re:= makemat(0,2,2);
FOR j FROM 1 TO rowDim(anovamat) DO
xio(j):= ΣLIST(row(anovamat, j));
varrow(j):= (stddevp(row(anovamat, j))).^2; // Variance of rows
END; // for
FOR j FROM 1 TO colDim(anovamat) DO
xoj(j):= ΣLIST(col(anovamat, j));
varcol(j):= (stddevp(col(anovamat, j))).^2; // Variance of columns
END; // for
Xio:=xio/n;
Xoj:=xoj/m;
Xoo:= (ΣLIST(Xio))/m;
FOR j FROM 1 TO rowDim(anovamat) DO
FOR k FROM 1 TO colDim(anovamat) DO
SSE:= SSE + ((anovamat(j,k)-Xio(j)-Xoj(k)+Xoo)^2);  // Sum square of errors
re(j,k):=anovamat(j,k)-Xio(j)-Xoj(k)+Xoo; // matrix of residues
END; END; // for
SSE:= ΣLIST(SSE);
SSR:= n*(ΣLIST((Xio-Xoo)^2));  // Sum square of treatments
SSC:= m*(ΣLIST((Xoj-Xoo)^2)); //Sum square of blocks
STR:= (SSR/(m-1))/(SSE/N); // Test F for rows
STC:=(SSC/(n-1))/(SSE/N); // Test F for columns
MST:= (SSR)/(m-1); // Mean Square Treatments
MSB:= (SSC)/(n-1); // Mean Square Blocks
MSE:= (SSE)/((n-1)*(m-1)); // Mean Square Errors

fr:= fisher_icdf(m-1,N,1-a);
fc:= fisher_icdf(n-1, N, 1-a);
p1:= 1- fisher_cdf(m-1, N, STR);
p2:= 1- fisher_cdf(n-1, N, STC);
sp:= sqrt(MSE);  // pooled std_dev
sto(re, anova_res);
sto({STR, STC}, anova_test_row_col);
sto({fr, fc}, anova_F_crit);
sto({m-1, N, n-1, N} , anova_dfreedom2);
sto({SSR, SSC, SSE}, anova_ssr_ssc_sse);
sto({MST, MSB, MSE}, anova_mst_msb_mse);
sto({xio,xoj, ΣLIST(xio)},  anova_xio_xoj_xij);
sto({Xio, Xoj, Xoo}, anova_Xio_Xoj_Xoo);
sto({varrow, varcol}, anova_var_rowcol);
sto({p1, p2}, anova_pvalue2);
sto(sp, anova_pooled_stddev);

IF STR< fr THEN MSGBOX("Test ST_row < F Fisher: Fail to reject H₀ at α=" + eval('a')); END;
IF STR >= fr THEN MSGBOX("Test ST_row > F Fisher: Reject H₀ at α=" + eval('a')); END;
IF STC< fc THEN MSGBOX("Test ST_col < F Fisher: Fail to reject H₀ at α=" + eval('a')); END;
IF STC >= fc THEN MSGBOX("Test ST_col > F Fisher: Reject H₀ at α=" + eval('a')); END;
mesg:= "Degrees of freedom: "+ eval( '(m-1)')+ " and " + eval( 'N') + "; " + eval('(n-1)') + " and " + eval('N');
MSGBOX(mesg);

RETURN ("ST_row="+ STR+ " (F=" + fr + ") ;  ST_col="+ STC + " (F=" + fc + ")");
END;

EXPORT ANOVA2_help()
BEGIN
local mesg;
mesg:="ANOVA 2 way (© Salvo Micciché 2015)

ANOVA2(): Calc ANOVA with a (α) = significance level (ex. 0.05) and a 2 factors matrix

Program returns test’s results for rows and columns and exports anovamat matrix";
PRINT;
PRINT(mesg);
RETURN("Thanks for using");
END;