Analysis of variance (ANOVA), orthogonal polynomials, and analysis of means (AOM)#
Literature#
See also: Howell, 2010
ANOM: Hothorn_AnalysisOfMeans, Elamir_2016
Nonparametric: Govindarajulu_2007
Genz_Bretz_2003
General: CRC - standard probability and statistic tables
Kirk, page 415
Lübschen
Power of Ancova: https://cran.r-project.org/web/packages/Superpower/vignettes/ANCOVAs.html (Shieh, 2020)
See also: https://en.wikipedia.org/wiki/Multiple_comparisons_problem
Anova: overview#
See also: https://en.wikipedia.org/wiki/One-way_analysis_of_variance
See also: https://en.wikipedia.org/wiki/Repeated_measures_design#Repeated_measures_ANOVA
We consider two different models:
The completely randomized model (CR): \(k\) independent samples of size \(n_i\) with mean \(\bar{x}_{i}\) and standard deviation \(s_i\). We define \(\displaystyle N = \sum_{i=1}^k n_i\) and \(\displaystyle \bar{\bar{x}} = \frac{1}{N} \sum_{i=1}^k n_i \bar{x}_{i}\). For \(k=2\) this is equivalent to a two-sided t-test for 2 independent samples.
The randomized blocks model (RB): \(k\) correlated samples of common size \(n_i=n\) and common correlation \(\rho\), with \(\bar{x}_{i}\), \(s_i\), \(N\) and \(\bar{\bar{x}}\) as above. For \(k=2\) this is equivalent to a two-sided t-test for 2 correlated samples.
TODO: Generalized randomized block design (Kirk, page 293).
We assume normally distributed data, homogeneity of variances for CR and the validity of the circularity assumption for RB. When referring to populations, we will use \(\bar{\mu}\), \(\mu_{i}\) and \(\sigma_i\) for means and standard deviations instead of \(\bar{\bar{x}}\), \(\bar{x}_{i}\) and \(s_i\). In the following sections, we will refer to the following definitions:
Anova (completely randomized and randomized blocks): p-value#
- ctx.anova_test(mean, sd, n, rho=none)#
where
ctxisdec,mpm, orgmp.Returns the results of tests for a CR or RB Anova.
We consider two different models:
With this definition of \(\mathit{MS}_{\text{Err}}\) and \(\mathit{df}_{\text{Err}}\), we can state that the test-criterion \(\displaystyle F = \frac{\mathit{MS}_{\text{Treat}}}{\mathit{MS}_{\text{Err}}}\) follows (under \(H_0\)) a central F-distribution with \(\mathit{df}_{\!1} = \mathit{df}_{\text{\!Treat}}\) and \(\mathit{df}_{\!2} = \mathit{df}_{\text{\!Err}}\) degrees of freedom, for both models.
A call to the function, requesting a CR Anova with the critical value for a two-sided test, the p-value for \(H_{03}\) , for 2 independent samples of size 10 and standard deviation 1 each, with means 2.3 and 4.5, and a type I error \(\alpha=0.05\) would be
>>> from mpfunlab import mpm >>> xreal.AnovaTest(means:=[5.24, 4.05, 7.01], sd:=1.5, n:=[22,11,16]) Input Variable Variable 1 Common N 22 Mean. Group 1 5.24 Mean. Group 2 4.05 StDev. Group 1 1.5 StDev. Group 2 3.5 Anova. k independent groups Parameter Result SSTreatment 15.5771 SSError 304.5 SSTotal 320.0771 dfTreatment 1 dfError 42 dfTotal 43 MSTreatment 15.5771 MSError 7.25 F-Value 2.148565517 p-value 0.150150344 F1-alpha 4.072653759
A call to the function, requesting a RB Anova with common correlation 0.5, the critical value for a two-sided test, the p-value for \(H_{03}\), for 3 independent samples of size 10 and standard deviation 1 each, with means 2.3 and 4.5, and a type I error \(\alpha=0.05\) would be
>>> from mpfunlab import mpm >>> xreal.AnovaTest(means:=[5.24, 4.05, 7.01], sd:=1.5, n:=22, rho=0.5) Input Variable Variable 1 Common N 22 Mean. Group 1 5.24 Mean. Group 2 4.05 StDev. Group 1 1.5 StDev. Group 2 3.5 Anova. k independent groups Parameter Result SSTreatment 15.5771 SSError 304.5 SSTotal 320.0771 dfTreatment 1 dfError 42 dfTotal 43 MSTreatment 15.5771 MSError 7.25 F-Value 2.148565517 p-value 0.150150344 F1-alpha 4.072653759
Anova (completely randomized and randomized blocks): power#
- ctx.anova_power(mean, sd, alpha=0.05, beta=0.1, rho=none, eta=0)#
where
ctxisdec,mpm, orgmp.Returns the results of power calculations for a CR or RB Anova.
See also: https://tjmurphy.github.io/jabstb/posthoc.html
>>> from mpfunlab import mpm >>> xreal.AnovaTest(means:=[5.24, 4.05, 7.01], sd:=1.5, n:=[22,11,16]) Input Variable Variable 1 Common N 22 Mean. Group 1 5.24 Mean. Group 2 4.05 StDev. Group 1 1.5 StDev. Group 2 3.5 Type 1 Error 0.05 Anova. k independent groups Parameter Result SSTreatment 15.5771 SSError 304.5 SSTotal 320.0771 dfTreatment 1 dfError 42 dfTotal 43 MSTreatment 15.5771 MSError 7.25 F-Value 2.148565517 F1-alpha 4.072653759 NC parameter 2.148565517 Power 0.299230817
Anova (completely randomized and randomized blocks): sample size#
- ctx.anova_samplesize(mean, sd, alpha=0.05, beta=0.1, rho=none, eta=0)#
where
ctxisdec,mpm, orgmp.Returns the results of sample size calculations for a CR or RB Anova.
Anova, trend tests using orthogonal polynomials: p-value#
- ctx.orthogonal_poly_test(mean, sd, n, rho=none)#
where
ctxisdec,mpm, orgmp.Returns the results of orthogonal polynomials for a CR or RB Anova.
See Kirk, page 152 and page 773.
\(\mathit{SS}_{\text{Treat}}\) can be partitioned into \(i = 1, \ldots, p-1\) sums of squares \(\displaystyle SS_{\psi_i} = \left(\sum_{j=1}^k n_j c_{ij} \bar{x}_{i} \right)^2 \bigg/ \sum_{j=1}^k c_{ij}^2\), reflecting orthogonal trend (linear, quadratic, cubic etc.) contrasts, such that \(\mathit{SS}_{\text{Treat}} = SS_{\psi_1} + SS_{\psi_2} + \ldots SS_{\psi_{k-1}}\), where the coefficients \(c_{ij}\) can be computed according to the procedure given in Narula.
>>> from mpfunlab import mpm >>> xreal.AnovaTest(means:=[5.24, 4.05, 7.01], sd:=1.5, n:=[22,11,16]) Input Variable Variable 1 Common N 22 Mean. Group 1 5.24 Mean. Group 2 4.05 Mean. Group 3 4.05 StDev. Group 1 1.5 StDev. Group 2 3.5 StDev. Group 3 3.5 Anova. k independent groups. orthogonale polynomials Parameter Polynomial 1 Polynomial 2 Treatment Error Total Variable Variable Variable Variable Variable Variable df 2 2 2 63 65 SSl 15.5771 5.192366667 20.76946667 561.75 582.5194667 MS 15.5771 5.192366667 10.38473333 8.916666667 F-Value 1.746964486 0.582321495 1.164642991 p-value. 2-sided (without adjustment) 0.19103836 0.448251923 0.318660397 p-value. 1-sided (SMM) 0.181523004 0.397610496 p-value. 2-sided (SMM) 0.344015295 0.693938137
Anova, trend tests using orthogonal polynomials: confidence interval#
- ctx.orthogonal_poly_ci(mean, sd, n, rho=none)#
where
ctxisdec,mpm, orgmp.Returns the results of orthogonal polynomials for a CR or RB Anova.
Anova, trend tests using orthogonal polynomials: power#
- ctx.orthogonal_poly_power(x, k, n, method='default')#
where
ctxisdec,mpm, orgmp.Returns the results of orthogonal polynomials for a CR or RB Anova.
Anova, trend tests using orthogonal polynomials: sample size#
- ctx.orthogonal_poly_samplesize(x, k, n, method='default')#
where
ctxisdec,mpm, orgmp.Returns the results of orthogonal polynomials for a CR or RB Anova.
Analysis of means (ANOM): p-value#
- ctx.anom_test(mean, sd, alpha=0.05, beta=0.1, rho=none)#
where
ctxisdec,mpm, orgmp.Returns the results of the ANOM-test for a CR Anova.
See also: Hothorn_AnalysisOfMeans, Elamir_2016
Distribution fitting#
- ctx.dist_fit(mean, sd, alpha=0.05, beta=0.1, rho=none)#
where
ctxisdec,mpm, orgmp.Returns the results of a distribtion fit
See also: https://en.wikipedia.org/wiki/Probability_distribution_fitting
See also: https://en.wikipedia.org/wiki/Method_of_moments_(statistics)
See also: https://en.wikipedia.org/wiki/Maximum_likelihood_estimation
See also: https://towardsdatascience.com/maximum-likelihood-estimation-in-r-b21f68f1eba4
See also: https://stat.ethz.ch/R-manual/R-devel/library/stats4/html/mle.html
Additional info: Narula’s algorithm#
The problem of approximating the relation of a response variable to a predictor variable consists of (i) the determination of the degree of the polynomial that “best” describes any trends in the data and (ii) the evaluation of the coefficients of the approximating polynomial. One may approach the problem in one of two ways: (a) the usual polynomial regression, or (b) the orthogonal polynomial approach.
The usual polynomial regression approach gives rise to the following difficulties: (i) the tests of significance for the various parameters are not independent, (ii) the estimates of the parameters of a polynomial of degree \((k-1)\) cannot be used for estimating coefficients of a polynomial of degree \(k\), and (iii) for degrees of polynomial greater than 4, the estimates may be inaccurate due to round-off errors.
The orthogonal polynomial approach does not suffer from the shortcomings associated with the usual polynomial regression.
The following algorithm has been given by Narula [445]: Let \(x_i\) denote the \(i^{\text{th}}\) value \((i=1,\ldots,k)\) of the predictor variable observed \(n_i\) times. Then
represents a set of orthogonal polynomials if
The orthogonal polynomial model can be written as
where \(Y_{ij}\) denotes the \(j^{\text{th}}\) values \((j=1,\ldots,n_i)\) of the response variable corresponding to \(x_i\), the \(i^{\text{th}}\) value of the predictor variable, and \(P_j\) denote the orthogonal polynomials as given above.
The least quares estimators \(b_j^*\) of \(\beta_j^*\) \((j=0,1,\ldots,m)\) are given by
and the sum of squares due to the \(j^{\text{th}}\) degree polynomial, \(P_j\) , is