Post Reply 
ANOVA
02-28-2015, 05:47 PM (This post was last modified: 03-02-2015 10:43 PM by salvomic.)
Post: #1
ANOVA
hi all,
I'm reading this thread, and I'd like to calculate ANOVA with Prime.
Any hints, help, tricks to operate with the commands and apps that Prime has is welcome; also hints to make a simple program.

It would be very interesting can use Statistic (one or two var) and Inference App for that purpose also, otherwise to have a dedicated app to analyze variance...

Help is much appreciated, thanks!

Salvo

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
03-04-2015, 06:57 PM (This post was last modified: 03-05-2015 08:04 PM by salvomic.)
Post: #2
RE: ANOVA
please, advice to me some tips...
I'm starting to focalize the problem for ANOVA...
I would get two results:
1. input G groups of data with items, calculating "ANOVA table" (mean, Sum of squares, degree of freedom, SSA, SSW, SST, F and so on)
2. giving the minimal useful data (summary result), obtain the same table.

About 1st purpose, the user must insert items, and a good location could be "Statistics 2Var" and it's columns, so we can use C1, C2,... and treat also C1(1), C1(size(C1) and so on...
Only we cannot use more than 10 groups (from C1 to C0).
I would like to have a sufficient numbers of groups and don't have limit for their size (items) and make that they could be possibly of different size...
However, then we should perform a cycle to count how much column are been used and so on...

I wonder if there is a better way to do that, as Spreadsheet perhaps is not suited...
Any opinion?

After that we must created formulae (we could start from here or here; second link is in Italian, sorry), but I'm a bit confused for now and I need help...

Last but not least thing: to make an inference test confronting F with Fisher (already in Prime) to validate or not H0 hypothesis...

Eventually we need a nice and well ordered routine to present the results.

I know that it's a not easy job, almost for me. I sincerely hope someone would collaborate, help, advice...

Thanks in advance :-)

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
03-08-2015, 04:16 PM
Post: #3
RE: ANOVA
ok,
a first draft...

Input data items via Statistic 1var in groups (columns) >= 2
Then use anova(a) where a is alpha, confidence level (i.e. 0.05) with a>0 and <1
Program give a message to say if the test for H0 is to accept (estimation ST > F Fisher) or to refuse, another with degree of freedom, then the value for ST (estimation) and F (with m-1 and m(n-1) d. of fr.)

Code:

EXPORT ANOVA(a)
// ANOVA 1 factor a=α confidence level (es. 0.05)
// use Statistic 1var to input groups
// by Salvo Micciché 2015
BEGIN
local j, k, m, n, g, me:={}, mm;
local s2:={}, ss,  num, den, st, f, msg;
local DL:={eval('D1'),eval( 'D2'), eval('D3'), eval('D4'), eval('D5'), eval
('D6'), eval('D7'), eval('D8'), eval('D9'), eval('D0')};
IF (a <=0 OR a> 1) THEN RETURN("Confidence level must be > 0 AND <1"); END;
m:= 0; 
FOR j FROM 0 TO 9 DO
IF size(DL(j) <> 0) THEN  m := m+1;  END; 
END; // count groups
n:= size(D1); // numerosity

FOR j FROM 1 TO m DO
me(j):= ΣLIST(DL(j))/n; // list of means
s2(j):= ΣLIST((DL(j)-me(j))^2) / (n-1); // list of variances within groups
END; // cacl mean

mm:= ΣLIST(me)/m; // mean of means, sample mean
ss:=  ΣLIST((me - mm)^2) /(m-1); // variance between groups
num:= n*ss; // numerator of test
den:= ΣLIST(s2/m); // denominator = mean of variances
st:= num/den; // estimation
f:= fisher_icdf(m-1,m*(n-1),1-a);
IF st < f THEN MSGBOX("Estimate ST < F Fisher: H0 to accept"); END;
IF st >= f THEN MSGBOX("Estimate ST > F Fisher: H0 to refuse"); END;
msg:= "Degrees of freedom: "+ eval( '(m-1)')+ " and " + eval( 'm*(n-1)');
MSGBOX(msg);
RETURN ("ST ="+ st + " F = "+ f);
END;

I need help for beta testing and a big tip to get RETURN() with multiple values: for now program give only a string, I would like it'd give numeric values with descriptive strings for ST, F, and also others (list of means, variances...)
Thanks to the user DrD for the tip with the columns!

Salvo

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
03-08-2015, 04:20 PM
Post: #4
RE: ANOVA
Salvomic, are you a purely HP man?

The Casio Algebra FX 2.0Plus has a great ANOVA & statistics environment.

It's useful to have a comparison anyway.
Find all posts by this user
Quote this message in a reply
03-08-2015, 04:27 PM
Post: #5
RE: ANOVA
(03-08-2015 04:20 PM)Gerald H Wrote:  Salvomic, are you a purely HP man?

The Casio Algebra FX 2.0Plus has a great ANOVA & statistics environment.

It's useful to have a comparison anyway.

hi Gerald,
I'm using also TI (89 and Voyager), but mostly HP (I've HP 12C, 15C, 50g, Prime), they are good for Statistics (bu haven't got ANOVA), as my old Casio FC 200...

I don't know FX 2.0plus, I'm searching for it on Internet (thanks) Smile

Also StatMate on iOS (iPhone, iPad) is a good tool...

but I like much Prime, so I want also make ANOVA with it Smile

Salvo

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
03-08-2015, 10:07 PM (This post was last modified: 03-14-2015 09:42 AM by salvomic.)
Post: #6
RE: ANOVA
Second (and 3rd, and 4th) try.
Now the program works also with groups of different numerosity.

Added a new function to calc ANOVA(a) when (in Statistic 1var, D1, D2, D3) are given a list of means, a list of variances and a list of numerosity: ANOVA_list(a)
a is alpha (es. 0.05).
Added the calc of p-value, "post hoc" analysis test -balanced- LSD (least significant difference) and "Pooled standard deviation" (sp).

Reviewed and simplified formulae.
Now the program exports some useful variables: list of means of the groups, p_value, list of variances within groups, test's result, mean and variance between groups, SSG, MSG, SSE, MSE (sum square between, within and total), mean square between and within)...
It gives also degrees of freedom in a msgbox and as variable.

Added help inline: ANOVA_help()

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;

ANOVA_help();

EXPORT ANOVA(a)
// ANOVA 1 way a=α significance level (ex. 0.05)
// use Statistic 1var to input groups (treatments)
// by Salvo Micciché 2015
BEGIN
CASE
IF (SIZE(a)==0) THEN a:=0.05; END;
IF (TYPE(a)<>0) THEN ANOVA_help(); RETURN("Usage: ANOVA(a)");  END;
IF (a <=0 OR a>= 1) THEN RETURN("Significance level must be > 0 AND <1"); END;
DEFAULT
END; // case

local j, m, n, N, g, me:={};
local mm, msg, mse, sp;
local p, ssg, tn, xij, sst, sse, sw:={};
local xi2:={}, ss,  num, den, st, f, mesg;
local DL:={eval('D1'),eval( 'D2'), eval('D3'), eval('D4'), eval('D5'), eval ('D6'), eval('D7'), eval('D8'), eval('D9'), eval('D0')};
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 1var (D1 ... D0) 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
tn:= tn + ΣLIST(DL(j));
xij:= xij + ΣLIST(DL(j)^2);
END; // calc  mean and variances
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 H0 at α=" + eval('a')); END;
IF st >= f THEN MSGBOX("Test ST > F Fisher: Reject H0 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_list(a)
BEGIN
// ANOVA 1 factor with a=α significance level
// use Statistics 1Var (D1, D2, D3) to input: D1  list m of means, D2 stddev (σ), D3 number of threatments
CASE
IF (SIZE(a)==0) THEN a:=0.05; END;
IF (TYPE(a)<>0) THEN ANOVA_help(); RETURN("Usage: ANOVA(a)");  END;
IF (a <=0 OR a>= 1) THEN RETURN("Significance level must be > 0 AND <1"); END;
IF ((size(D1) == 0) OR (size(D1) == 0) OR (size(D3) = 0)) THEN MSGBOX("ANOVA with synthetic data: use Statistics 1Var (D1, D2, D3) to input: D1: list m of means, D2: stddev (σ), D3: n numerosity of each threatment"); 
RETURN("Statistics 1var, D1 means, D2 stddev, D3 numerosity"); END;
DEFAULT
END; // case

local  m, N, n, mm, ss, num, den;
local me, s2, p, st, f, mesg;
local ssg, sse, msg, mse, sp;
m:= SIZE(D1); // num of treatments
N:= ΣLIST(D3); // total number of n
n:= mean(D3); // mean numerosity
me:= D1; // means of tratments
s2:= D2; // standard deviations of treatments
sto(me, anova_means);
sto((s2).^2, anova_var_within);

mm:= mean(D1); // mean of means, sample mean
ss:=  ΣLIST((me - mm).^2) /(m-1); // variance between groups
ssg:= ΣLIST( (D1*D3).^2 / D3 ) - (ΣLIST( (D1*D3)).^2 /N); // sum squares for groups;
sse:= ΣLIST( ((D2).^2)* (D3 -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 H0 at α=" + eval('a')); END;
IF st >= f THEN MSGBOX("Test ST > F Fisher: Reject H0 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(a): Input a (α) = significance level (ex. 0.05)
First use Statistic 1var (D1 ... D0) to input groups (treatments)

ANOVA_list(a): Input a (α) = significance level (ex. 0.05)
ANOVA with synthetic data: use Statistics 1Var (D1, D2, D3) to input:
D1: list m of means, D2: stddev (σ), D3: n numerosity of each threatment";
PRINT;
PRINT(mesg);
RETURN("Thanks for using");
END;

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
03-12-2015, 06:12 PM (This post was last modified: 03-12-2015 10:10 PM by salvomic.)
Post: #7
RE: ANOVA
Tim, Parisse, Han and everybody, definitely this is not the best possible program for the calculation of ANOVA (I would like only to give "my two cents", but I'd like someone in the Forum could improve it or make another more powerful program to be able to put (in a next FW update) into XCas and Inference application of the Prime.

However I don't know if it is possible to do it, so, please, let us know your opinion, if you want Smile
Thanks in advance.

Salvo

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
03-13-2015, 12:44 AM
Post: #8
RE: ANOVA
(03-12-2015 06:12 PM)salvomic Wrote:  Tim, Parisse, Han and everybody, definitely this is not the best possible program for the calculation of ANOVA (I would like only to give "my two cents", but I'd like someone in the Forum could improve it or make another more powerful program to be able to put (in a next FW update) into XCas and Inference application of the Prime.

However I don't know if it is possible to do it, so, please, let us know your opinion, if you want Smile
Thanks in advance.

Salvo

I have thought about various kinds of popular ANOVA calculations. My approach is to have a function that does the calculations and another that reports the results back to the caller (or client).
Find all posts by this user
Quote this message in a reply
03-13-2015, 09:47 AM (This post was last modified: 03-13-2015 06:09 PM by salvomic.)
Post: #9
RE: ANOVA
(03-13-2015 12:44 AM)Namir Wrote:  I have thought about various kinds of popular ANOVA calculations. My approach is to have a function that does the calculations and another that reports the results back to the caller (or client).

good idea, Namir.
My program does the calculation, then put the most important result also in a few variables, so the user can use them for other calculation or for a "post hoc" analysis.
A program for some kind of post hoc analysis would be good also; my program gives also two variables for test: LSD and pooled standard deviation, but no tests, for now...

I thank you also for the good advice (and short programs) on this thread!

have a nice day

Salvo

***
EDIT: then I need also hints to try writing a program for ANOVA "two factors": how better is to input data? a matrix? Statistics? ...
In this case we have a table in witch the two factors are in the rows and in the columns, and we must calculate data (sum, squares, means) for rows, for columns and total, but a matrix perhaps is not so easy for this purpose (fancy a matrix with 8 rows and 10 columns...)

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
03-17-2015, 08:01 PM (This post was last modified: 03-22-2015 11:02 PM by salvomic.)
Post: #10
RE: ANOVA
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;

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
03-18-2015, 06:09 PM (This post was last modified: 03-25-2015 10:04 AM by salvomic.)
Post: #11
RE: ANOVA
ANOVA one way revisited.

This version doesn't require the parameter alpha (significance level) at the prompt: input it inside of the program -> anova(), anova_synt() (and anova_help() )...
Data now use Statistic 2Var, and so with ANOVA_synt() C1: means, C2: stddevs, C3: numerosity; with ANOVA() every treatment in columns C1, C2, ..., C0 (10 max), al always...

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

I corrected also a few typos.

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;

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
01-07-2021, 05:22 PM
Post: #12
RE: ANOVA
You can try some free tool to calculate ANOVA such as https://statsjournal.com/anova-test

Let me know is it work.
Find all posts by this user
Quote this message in a reply
Post Reply 




User(s) browsing this thread: 1 Guest(s)