Estimate factorial experimental design by using some methods with an application

Introduction

Experimental design plays a central role in applied statistics, providing a structured framework for planning experiments, controlling variability, and drawing reliable inferences. Among these designs, factorial experiments-originally introduced by Fisher (1926) [13] and further advanced by Yates (1937) enable the simultaneous investigation of multiple factors and their potential interactions. This multifactorial approach enhances efficiency by reducing the number of required experimental runs, thereby saving time, cost, and resources. In agricultural research, factorial designs are particularly valuable for exploring how various controllable factors-such as fertilizer types, irrigation levels, and environmental conditions-influence crop productivity and sustainability. However, traditional analytical methods like Analysis of Variance (ANOVA) often rely on assumptions such as normality and homoscedasticity, which may not hold in real-world data. When these assumptions are violated, the validity and reliability of the conclusions drawn from such methods are compromised. To address these limitations, this study evaluates the performance of six statistical approaches for analyzing data from factorial experiments in the presence of non-normal and heterogeneous errors: Ordinary Least Squares (OLS), Generalized Linear Models (GLM), Quantile Regression, the Least Absolute Shrinkage and Selection Operator (LASSO), Adaptive LASSO, and the MMestimator. The aim is to identify the most robust and accurate estimation technique under these challenging data conditions, thereby providing practical guidance for researchers conducting factorial experiments in applied fields such as agriculture.

Experiment Design

An experimental design identifies the independent and dependent variables and includes the method by which the statistical and randomization aspects of the experiment will be carried out. Its primary goal is to establish a causal relationship between the independent and dependent variables, while the secondary goal is to extract the maximum amount of information with minimal resource expenditure. (Kirk, 1995) [2] . Experiment design refers to the plan used to carry out the experiment. It involves a series of steps that include the allocation of treatments to experimental units for the purpose of collecting data. This process requires a comprehensive understanding of how to analyze the data after it has been gathered from the experiment. It is essential for the design to allow for the estimation of experimental error, as the inability to estimate this error renders the design invalid. There are various types of experimental designs, each tailored to a specific objective. The researcher's goal and the nature of the study determine the appropriate design. A well-chosen design ensures the ability to overcome challenges that may arise during statistical analysis. Thus, the optimal design for an experiment will vary based on the factors being studied.

Factorial Experimental Design

Factorial experiments are used to study the effect of two or more factors on a particular phenomenon or specific property. Instead of conducting a simple experiment for each factor separately, which requires a lot of effort, time, cost, and experimental materials, a single experiment that includes all factors at once can be performed, providing many advantages. These advantages include eliminating the complexities and problems of repeating the experiment for each factor individually, as well as reducing costs, effort, time, and materials used. Additionally, this method does not require testing the independence of factors, as the experimental results provide comprehensive information about the effect of each factor and its interactions with other factors. This experiment, which involves multiple factors, is known as a "factorial experiment" or "multifactorial experiment" (Factorial Experiment) (Kamal Alwan Al-Mashhadani, 2010) [1] . Factorial experiments can be conducted in a wide range of experimental designs, which are organized in a manner suitable for simple designs and are statistically analyzed. The factorial experiment provides valuable information about the interaction effects between factors, which are of great importance in interpreting the data.

Methods of estimating factorial experiments 1.3.1 Ordinary Least Squares (OLS)

It is considered one of the statistical methods used to estimate parameters in linear models by minimizing the sum of squared errors between observed and predicted values. The sample treatment mean and its observed value respectively represent the "least squares estimator" and "least squares estimate" of the treatment mean (Dean & Voss, 1999) [12] . In one-way analysis of variance models, least squares estimators are derived by optimizing the following sum of squared errors function: Treatment effect for The estimated model is given by:

Predicted value for treatment . Least squares estimate of the overall mean.

: Least squares estimate of the treatment effect .The calculation of least squares estimators is done by differentiating the sum of squared errors with respect to each parameter and setting the derivatives to zero, leading to the following system of normal equations:

(3) (4) i=1,……… Where,: Grand total of all observations.: Total of observations for treatment .: Total number of observations.

Generalized Linear Models (GLM)

Are a method for modeling the relationship between a dependent variable and several independent variables. This type of model is characterized by the dependent variable following one of the exponential family distributions, making generalized linear models less restrictive than traditional regression models. These models provide a solution for modeling data that does not follow a normal distribution. Generalized linear models were first introduced by Neider and Wedderburn in 1972 in an attempt to expand the assumptions of traditional regression models to fit practical applications (Barnett & Dobson, 2008) [3] . Generalized linear models are less restrictive for several reasons, including [Evelien, 2014]:  They do not require constant variance, allowing for different variances.

 The relationship between the dependent and independent variables does not have to be linear; rather, there must be a linear relationship between the link function and the set of independent variables, which allows for the application of nonlinear models using generalized linear models.  Model parameters are estimated using Maximum Likelihood Estimation (MLE) or traditional Ordinary Least Squares (OLS).  In these models, the expected value of the response variable is replaced by the link function, making the error variance more stable. The link function is a linear combination of the explanatory variables. Components of the GLM Model.

Quantile Regression (QR)

For over two centuries, regression techniques in statistical modeling have been developed to determine the relationship between a dependent variable (outcome) and independent variables (covariates). Classical regression has been one of the most commonly used statistical methods for estimating effects at the mean, assuming that regression coefficients or covariate effects remain constant across the population. However, these average effects are not always of interest in many fields and may often be highly heterogeneous. Quantile Regression (QR), introduced by Koenker and Bassett in 1978, has gained increasing attention and has been widely used in exploring relationships such as those between foreign direct investment and economic growth, as well as in precision health and personalized medicine. QR is capable of capturing heterogeneous effects through conditional quantiles of the outcome variable, offering a more comprehensive understanding of its distribution across different levels.

In cases of heteroscedasticity or data with asymmetries and heavy tails, the sample median (the 50th percentile) often provides a more accurate measure of central tendency than the mean. Therefore, QR is more flexible and robust compared to classical mean regression models, as it does not require strict distributional assumptions. While traditional regression methods still dominate statistical modeling, QR serves as a valuable extension and complement, particularly when classical assumptions are violated. QR has attracted considerable research interest and has been widely applied to independent data and time-to-event (survival) data. More recently, its applications have expanded to include longitudinal data analysis (Huang et al., 2023) [9] . In the context of improving product reliability, recent studies have shown that QR is an effective tool for analyzing lifetime data, particularly when dealing with issues such as heteroscedasticity or censored data. Unlike traditional methods such as Maximum Likelihood Estimation (MLE), which require strict distributional assumptions, QR offers a more flexible and robust approach to estimating quantiles without relying on such assumptions. To further enhance the performance of QR, the Weighted Quantile Regression (WQR) model has been proposed, introducing an adaptive weight parameter to improve stability and accuracy. Both theoretical and practical results have shown that WQR outperforms traditional QR, especially in estimating lower quantiles (e .g., 0.1), which are particularly crucial in reliability analysis (Lu & Zhao, 2023) [8] .

1.3.4

The Lasso Lasso (Least Absolute Shrinkage and Selection Operator) is a linear model estimation method proposed by Tibshirani (1996) [14] . This method shares with ridge regression the property of minimizing the sum of squared errors while imposing a constraint on the estimated coefficients. More specifically, while ridge regression uses a quadratic constraint, Lasso applies an upper bound o n the sum of absolute values of the coefficients. Tibshirani presents two main advantages of Lasso compared to ordinary least squar es (OLS) regression: First: OLS models that are over specified often exhibit high variance despite low bias, which reduces prediction accuracy. Lasso addresses this issue by shrinking some coefficients or setting them to zero (removing them from the model), enabling a trade-off between bias and variance. Second: OLS models typically contain many small coefficients that add little value to the model and increase its complexity without providing significant explanatory power. Lasso shrinks some coefficients and sets others to zero, automatically performing variable selection by including only non-zero effects in the final model. Let (x i, y i ), i=1,2…,N represent a traditional regression dataset with x i =(x 1, …,x ip ) as independent variables and y i as the response variable. Assuming standardized independent variables:

And assuming centered y i where ӯ =0, The Lasso estimator is defined as the set =( 1,….., p ) that minimizes:

Where t is the tuning parameter that controls the amount of shrinkage applied to the coefficient estimates. If are the full OLS estimates and t 0 = , then values of t<t 0 lead to shrinking coefficients toward zero, effectively excluding variables whose coefficients become zero from the model. Tibshirani recommends three different methods for estimating the bound t: K-fold cross-validation, generalized cross-validation, and a method based on Stein's unbiased risk estimator. K-fold cross-validation is used to estimate t. Suppose (7) Where E( =0 and var ( where ME is the mean squared error. To estimate (x) as ( 8)

And the prediction error PE as follows:

This prediction error is then estimated using K-fold cross-validation. If (10) Leave-one-out cross-validation (where K=n) is used to determine the optimal value of t The Lasso problem represents a quadratic programming problem that can be solved using numerical analysis algorithms (Efron et al., 2004, p. 37). However, the modified least angle regression (LARS) method can also be used, which provides the complete solution path for all values of t (from 0 to ∞).

Adaptive lasso

Adaptive LASSO is considered an advanced statistical technique used to enhance variable selection and parameter estimation in statistical models, particularly in cases involving a large number of explanatory variables, such as in factorial designs. This method was first introduced by Zou (2006), who developed an improved version of the traditional LASSO by incorporating adaptive weights calculated from initial estimates of the parameters, allowing for more precise and selective shrinkage. This method is distinguished by its ability to handle multicollinearity and improve computational efficiency, making it suitable for complex analytical scenarios (Ren & Zhang, 2013) [5] . The strength of Adaptive LASSO lies in its capability to simultaneously perform estimation and variable selection by integrating adaptive weights into the penalty function, which results in shrinking the coefficients of irrelevant variables to zero. This improves both model interpretability and predictive accuracy (Qian & Yang, 2013) [6] . The quantile regression model at a specific quantile level is represented by the following equation:

Where, :is the conditional quantile function of the response variable y,: is the i-th row vector of the design matrix X,: is the parameter vector at quantile . the linear regression model can be represented as follows (12) To estimate the regression coefficients while performing variable selection, the Adaptive LASSO method solves the following penalized optimization problem:

:The initial estimate of the coefficient , usually obtained from an Ordinary Least Squares (OLS) or Ridge regression fit.

These estimates are used to compute the adaptive weights.

:A positive tuning parameter (commonly set to 1) that controls the influence of the initial estimates on the weights . It is often selected by the user or through cross-validation. The following equation represents the Lasso-penalized quantile regression model (Quantile Regression with Lasso Penalty), which combines the advantages of quantile regression with the regularization capability for selecting the most important variables:

The quantile loss function (\(\rho_{T}\)) used in the quantile regression model is defined by the following relation:

Quantile Factorial Design

Quantile regression models are increasingly employed in factorial designs, particularly when three factors are involved, to overcome the limitations of traditional models-especially in the presence of heteroscedasticity and outliers. Quantile regression estimates conditional quantiles of the response, allowing researchers to understand response behavior across the distribution's tails and gain insights into the effects of different factors at various levels of the response. The general form of the quantile regression model is as follows:

The conditional quantile of the response variableY at quantile level , given factorsA,B,C.

: The quantile of interest (e.g., 0.25, 0.50, 0.75).

: Intercept term; the baseline value of the -th quantile when all factors are at their reference level (e.g., zero).

: Main effect of factor A at quantile .

Simulation Example

We simulated data from a full factorial experimental design with 2 3 factors (A, B, and C), each taking two levels (-1 for low and +1 for high). The true model included all main effects, two-way interactions, and the three-way interaction:

Where the coefficients were fixed at known values and the error term followed a non-normal, heavy-tailed distribution to mimic real-world deviations from Gaussian noise. The experiment was simulated under three sample sizes (N=50,100 and 200) by replicating the factorial runs proportionally.

For each dataset, six estimation methods were applied: 1. Ordinary Least Squares (OLS) 2. Generalized Linear Model (GLM, Gaussian family) 3. Quantile Regression 4. Least Absolute Shrinkage and Selection Operator (LASSO) 5. Adaptive LASSO 6. MM Estimator (robust regression)

Model performance was using mean squared error (MSE), coefficient bias, and coefficient of determination (R 2 ). Each scenario was replicated 500 times for the The simulation results show that Quantile Regression consistently outperformed other methods in terms of bias and explained variance (R 2 ) across all sample sizes, demonstrating its robustness to the non-normal error structure. While OLS showed higher mean squared error (MSE) and lower R 2 , reflecting sensitivity to skewed noise, the Generalized Linear Model (GLM) and MM Estimator provided more stable and accurate estimates. Regularization methods behaved differently: Adaptive LASSO improved performance over standard LASSO, particularly at the largest sample size (N=200), where it achieved the lowest MSE and competitive R 2 . The MM estimator, with its robustness to outliers and heavy-tailed errors, offered consistent performance across all sample sizes, balancing bias and variance effectively.

Real data

The dataset employed in this analysis is the classic Iris dataset, originally introduced by R.A. Fisher (1936) [13] , which contains morphological measurements of 150 iris flowers from three species: Iris setosa, Iris versicolor, and Iris virginica. For each specimen, four continuous variables are recorded: sepal length (cm), sepal width (cm), petal length (cm), and petal width (cm). The present study utilizes sepal length as the response variable, with petal length and species as predictors, including thei r interaction. This configuration enables the evaluation of both main effects and interaction effects between a continuous morphological trait and species identity in predicting sepal length, within a well-balanced, real-world dataset. Figure 4 demonstrate that the Sepal length grows with petal length in all species and the Versicolor and virginica start with lower baseline sepal lengths but have steeper slopes than setosa. Whereas, the graph will show these three lines starting at different intercepts but all increasing positively, with virginica being the steepest. This scatter plot in Figure 4, shows the relationship between Petal Length (x-axis) and Sepal Length (y-axis) for the three iris species; setosa (red), versicolor (green), and virginica (blue) based on the classic iris dataset. A quantile regression line for the median (τ=0.5) is fitted separately for each species. Overall, the plot illustrates that as petal length increases, sepal len gth also tends to increase, with the strength of this relationship varying by species. Quantile regression provides a robust estimate of the median relationship, making it less sensitive to outliers than ordinary least squares regression.

Conclusion

These findings suggest that for factorial designs with non-normal errors, robust and regularized methods like Quantile Regression and Adaptive LASSO are preferable to classical least squares. Simulation results demonstrate that Quantile Regression consistently yields the lowest bias and highest explained variance (R 2 ) across all sample sizes under non-normal errors. Classical OLS exhibited higher mean squared error and lower R 2 , reflecting sensitivity to skewness. Adaptive LASSO outperformed standard LASSO, especially at larger samples, achieving the lowest MSE. The MM estimator provided robust and stable estimates throughout. Overall, robust and regularized approaches like Quantile Regression and Adaptive LASSO are recommended for factorial designs with non-normal noise, outperforming classical methods.

References

1. Al-Mashhadani AK. Design and analysis of experiments using computer. 1 st Ed. Baghdad: Al-Jazeera Printing and Publishing Office, Deposit no. 864, House of Books and Documents, Baghdad.; 2010.
2. Kirk RE. Experimental design: Procedures for the behavioral sciences. 3rd ed. Pacific Grove (CA): Brooks/Cole; 1995.
3. Dobson AJ, Barnett AG. An introduction to generalized linear models. 3rd ed. London: Chapman & Hall; 2008.
4. Brisard E. Pricing of car insurance with generalized linear models. Manchester: University of Manchester; 2014.
5. Ren Y, Zhang X. Model selection for vector autoregressive processes via adaptive lasso. Commun Stat Theory Methods. 2013;42(13):2423-36. doi:10.1080/03610926.2011.611317
6. Qian W, Yang Y. Model selection via standard error adjusted adaptive lasso. Ann Inst Stat Math. 2013;65(2):295-318. doi:10.1007/s10463-012-0370-0
7. Pepler PT. Detecting active effects in unreplicated 2k factorial experiments: A simulation study. Stellenbosch: University of Stellenbosch; 2006.
8. Lu A, Zhao S. Robust quantile regression analysis on experimental design. Hebei: School of Mathematical Sciences, Hebei Normal University; 2023.
9. Huang Q, Zhang H, Chen J, He M. Quantile regression models and their applications: A review; 2023.
10. Xu M. Quantile regression model and its application research. Acad J Sci Technol; 2023. DOI: 10.54097/vt1qpm59
11. Silva PV da. Uso da regressão não linear quantílica na descrição de dados de crescimento. São Paulo: Universidade de São Paulo; 2022. DOI: 10.11606/t.11.2022.tde-10052022-155257
12. Dean A, Voss D. Design and analysis of experiments. New York: Springer-Verlag; 1999.
13. Fisher RA. The use of multiple measurements in taxonomic problems. Ann Eugen. 1936;7(2):179-88.
14. Tibshirani R. Regression shrinkage and selection via the lasso. J R Stat Soc Series B Stat Methodol. 1996;58(1):267-88.
15. Hinkelmann K, Kempthorne O. Design and analysis of experiments. Vol. 1: Introduction to experimental design. 2 nd Ed. Hoboken (NJ): John Wiley & Sons; 2008.
16. Maronna RA, Martin RD, Yohai VJ. Robust statistics: Theory and methods. Hoboken (NJ): John Wiley & Sons; 2006.