# HP Forums

Full Version: ANOVA
You're currently viewing a stripped down version of our content. View the full version with proper formatting.
hi all,
with the new firmware 7820 ANOVA one way is a part of Inference app (thanks to the HP Team to have put it!), but I would propose my version, that has more parameters and customization, for now...
Then I propose my version for ANOVA two Way and ANOVA for Regression Analysis too.
The data must be put into Statistic 2var App of the Prime.

Enjoy!
Salvo Micciché

ANOVA one way.

This version doesn't require the parameter alpha (significance level) at the prompt: it can be put inside of the program -> anova(), anova_synt() (and anova_help() )...
Data are taken by Statistic 2Var, and so in ANOVA_synt() the variables represent:
C1: means, C2: stddevs, C3: numerosity;
in ANOVA() every treatment are in columns C1, C2, ..., C0 (they can be 10 max), as usual...

This version calculates also a matrix of residues and export it.

Code:
``` export anova_means:= {0,0}; export anova_var_within:={0,0}; export anova_F_test:=0; export anova_m_s2:={0,0}; export anova_ssg_msg_sse_mse:={0,0,0,0}; export anova_pvalue:= 0; export anova_dfreedom:={0,0}; export anova_pooled_stddev:=0; export anova_LSD:=0; export anova_res:=[0,0]; ANOVA_help(); EXPORT ANOVA() // ANOVA 1 way a=α significance level (ex. 0.05) // use Statistic 2var to input groups (treatments) // by Salvo Micciché 2015 BEGIN local a, j, m, n, N, g, me:={}; local mm, msg, mse, sp, re:={}; local p, ssg, tn, xij, sst, sse, sw:={}; local xi2:={}, ss,  num, den, st, f, mesg; local DL:={eval('C1'),eval( 'C2'), eval('C3'), eval('C4'), eval('C5'), eval ('C6'), eval('C7'), eval('C8'), eval('C9'), eval('C0')}; INPUT(a, "ANOVA 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:= 0;  N:= 0; // m number of treqtments, N total number of items FOR j FROM 0 TO 9 DO IF size(DL(j) <> 0) THEN  m := m+1; N:= N+size(DL(j));   END; END; // count groups n:= N/m; // (max) numerosity of groups IF (m == 0 OR m == 1) THEN MSGBOX("First use Statistic 2var (C1 ... C0) to input groups (treatments), 2 groupr or more");  RETURN("Input data in Statistics 1var, 2 groups or more"); END; FOR j FROM 1 TO m DO me(j):= approx(mean(DL(j))); // list of means xi2(j):= (ΣLIST( DL(j) ) ^2) / SIZE(DL(j)) ;  //list of x_i^2/n_i  sw(j):= approx((stddevp(DL(j)))^2); // list of variances within groups re(j):=DL(j)-me(j); tn:= tn + ΣLIST(DL(j)); xij:= xij + ΣLIST(DL(j)^2); END; // calc  mean and variances anova_res:=list2mat(re); sto(me, anova_means); sto(sw, anova_var_within); ssg:= ΣLIST(xi2)-((tn^2)/N); // Sum of Square between sst:= xij-(tn^2/N); // Sum Square total sse:= sst - ssg; // Sum of square within mm:= ΣLIST(me)/m; // mean of means, sample mean ss:= approx(stddevp(me)^2); // variance between groups (σ2 of means) msg:= ssg/(m-1); // mean square between mse:= sse / (N-m); // mean square within (error) num:= msg; // numerator of test den:= mse; // denominator = mean of variances st:= num/den; // estimation F test f:= fisher_icdf(m-1, N-m, 1-a); // quantile F Fisher Snedecor p:= 1- fisher_cdf(m-1, N-m, st); // pvalue sp:= sqrt(mse);  // pooled std_dev sto(st, anova_F_test); sto({m-1, N-m} , anova_dfreedom); sto({mm, ss}, anova_m_s2); sto( {ssg, msg, sse, mse} , anova_ssg_msg_sse_mse); sto(p, anova_pvalue); sto(sp, anova_pooled_stddev); sto(sqrt((2*( ssg/(N-m))*fisher_icdf(1,N-m,1-a))/n) , anova_LSD); // Least Significant Difference IF st < f THEN MSGBOX("Test ST < F Fisher: Fail to reject H₀ at α=" + eval('a')); END; IF st >= f THEN MSGBOX("Test ST > F Fisher: Reject H₀ at α=" + eval('a')); END; mesg:= "Degrees of freedom: "+ eval( '(m-1)')+ " and " + eval( 'N-m'); MSGBOX(mesg); RETURN ("ST ="+ st + " F = "+ f + " p = " + p); END; EXPORT ANOVA_synt() BEGIN // ANOVA 1 factor with a=α significance level // use Statistics 2Var (C1, C2, C3) to input: C1  list m of means, C2 stddev (σ), C3 number of threatments IF ((size(D1) == 0) OR (size(D1) == 0) OR (size(D3) = 0)) THEN MSGBOX("ANOVA with synthetic data: use Statistics 2Var (C1, C2, C3) to input: C1: list m of means, C2: stddev (σ), C3: n numerosity of each threatment");  RETURN("Statistics 2var, C1 means, C2 stddev, C3 numerosity"); END; local a,  m, N, n, mm, ss, num, den; local me, s2, p, st, f, mesg; local ssg, sse, msg, mse, sp; INPUT(a, "ANOVA 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:= SIZE(C1); // num of treatments N:= ΣLIST(C3); // total number of n n:= mean(C3); // mean numerosity me:= C1; // means of tratments s2:= C2; // standard deviations of treatments sto(me, anova_means); sto((s2).^2, anova_var_within); mm:= mean(C1); // mean of means, sample mean ss:=  ΣLIST((me - mm).^2) /(m-1); // variance between groups ssg:= ΣLIST( (C1*C3).^2 / C3 ) - (ΣLIST( (C1*C3)).^2 /N); // sum squares for groups; sse:= ΣLIST( ((C2).^2)* (C3 -1)); // Sum square of errors msg:= ssg / (m-1); // mean square for groups mse:= sse / (N-m); // mean square error num:= msg; // numerator of test den:= mse; // denominator = mean of variances st:= num/den; // estimation f:= fisher_icdf(m-1, N-m, 1-a); p:= 1- fisher_cdf(m-1, N-m, st); sp:= sqrt(mse);  // pooled std_dev sto(st, anova_F_test); sto({m-1, N-m} , anova_dfreedom); sto({mm, ss}, anova_m_s2); sto( {ssg, msg, sse, mse} , anova_ssg_msg_sse_mse); sto(p, anova_pvalue); sto(sp, anova_pooled_stddev); sto(sqrt((2*( ssg/(N-m))*fisher_icdf(1,N-m,1-a))/n) , anova_LSD); // Least Significant Difference IF st < f THEN MSGBOX("Test ST < F Fisher: Fail to reject H₀ at α=" + eval('a')); END; IF st >= f THEN MSGBOX("Test ST > F Fisher: Reject H₀ at α=" + eval('a')); END; mesg:= "Degrees of freedom: "+ eval( '(m-1)')+ " and " + eval( 'N-m'); MSGBOX(mesg); RETURN ("ST ="+ st + " F = "+ f  +  " p = " + p); END; EXPORT ANOVA_help() BEGIN local mesg; mesg:="ANOVA 1 way (© Salvo Micciché 2015) ANOVA(): First use Statistic 2var (C1 ... C0) to input groups (treatments) Input a (α) = significance level (ex. 0.05) ANOVA_synt(): ANOVA with synthetic data: use Statistics 2Var (C1, C2, C3) to input: C1: list m of means, C2: stddev (σ), C3: n numerosity of each treatment Input a (α) = significance level (ex. 0.05)"; PRINT; PRINT(mesg); RETURN("Thanks for using"); END;```

ANOVA two way (with blocks and treatments)

ANOVA with blocks and treatments (2 way).
Values are in a matrix (anovamat). The user can 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;```

ANOVA for the Regression Analysis.
The program takes data from C1 (XList) and C2 (List) (Statistics 2var), so the user can do the normal Regression analysis and tests via Statistics and Inference App (built in), and then get other parameters for ANOVA Table from the program itself:
test F0 (Fisher) and t (Student) for slope and intercept, SSR, SSE, SST, MSR, MSE, b0 (intercept), b1 (slope), sxx, sxy, p value, r and R^2 (correlation, determination)...

(ANOVA for Regression isn't ANOVA 1 way or ANOVA 2 way but a form of ANOVA related to the regression analysis)

Code:
``` export anova_F0_f:={0,0}; export anova_T0_T1_t:={0,0,0}; export anova_ssr_sse_sst:={0,0}; export anova_msr_mse:={0,0}; export anova_b0_b1_sxx_sxy:={0,0,0,0}; export anova_s2:=0; export anova_r_R2:=0; export anova_pvalue:=0; export anova_DF:={0,0}; ANOVA_help(); EXPORT ANOVA_reg() // ANOVA for REgression Analysis //use Statistic 2var to input C1 (XList) and C2 (YList) // by Salvo Micciché 2015 BEGIN local a, ym, xm, res, sxx, sxy; local SST, SSR, SSE, MSR, MSE; local CL:={eval('C1'),eval( 'C2')}; local n, s2, DF, f, F0, p, mesg; local b0, b1, T0, T1, t, r; INPUT(a, "ANOVA 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; n:=size(C1); // numerosity xm:=mean(C1); ym:= mean(C2); // MeanY res:= Resid(); // Residues SSE:=ΣLIST(res^2); // Square Sum of Errors SST:=ΣLIST((C2-ym)^2); //Square Sum Total SSR:=SST-SSE; // Square Sum Regression sxx:= ΣLIST((C1-xm)^2); // Sum of square deviation of x sxy:=ΣLIST( (C1-xm)*(C2-ym) );  // Sum of crossed product of deviations b1:=sxy/sxx; // β₁ (Slope) b0:= ym-b1*xm; // β₀  (Intercept) s2:=SSE/(n-2); // σ^² variance of the observations of response MSR:= SSR; MSE:=SSE/(n-2); R2:= 1- (SSE/SST); // Determination r:= b1*sqrt(sxx/SST); // correlation DF:={1, n-2}; F0:= MSR/MSE; // Test F₀ f:= fisher_icdf(1, n-2, 1-a); // quantile F Fisher Snedecor p:= 1- fisher_cdf(1, n-2, F0); // pvalue T0:= b1/(sqrt(s2/sxx)); // Test T₀ (Slope) T1:= b0/sqrt(s2*((1/n)+xm^2/sxx) ); // Test T₁ (Intercept) t:= student_icdf(n-2, 1-a/2); // Quantile Student t sto({"F₀=", F0,  "F=", f}, anova_F0_f); sto({"T₀=", T0,  "T₁=", T1, "t=", t}, anova_T0_T1_t); sto({SSR, SSE, SST}, anova_ssr_sse_sst); sto({MSR, MSE}, anova_msr_mse); sto({b0, b1, sxx, sxy}, anova_b0_b1_sxx_sxy); sto({r, R2}, anova_r_R2); sto(DF, anova_DF); sto(s2, anova_s2); sto(p, anova_pvalue); MSGBOX( "Test Hypothesis: H₀: β₁=0 vs H₁: β₁ ≠ 0"); IF F0 < f THEN MSGBOX("Test F₀ < F Fisher: Fail to reject H₀ at α=" + eval('a')); END; IF F0 >= f THEN MSGBOX("Test F₀ > F Fisher: Reject H₀ at α=" + eval('a')); END; mesg:= "Degrees of freedom: "+ eval( '(1)')+ " and " + eval( 'n-2'); MSGBOX(mesg); RETURN ("F₀ ="+ F0 + " F = "+ f + " p = " + p); END; EXPORT ANOVA_reg_help() BEGIN local mesg; mesg:="ANOVA for Regression Analysis (© Salvo Micciché 2015) ANOVA(): First use Statistic 2var to input C1 (XList) and C2 (YList), then input significance level (α) "; PRINT; PRINT(mesg); RETURN("Thanks for using"); END;```
Reference URL's
• HP Forums: https://www.hpmuseum.org/forum/index.php
• :