## 2 Evolutionary Polynomial Regression (EPR)  
EPR is a smart regression method that automatically searches for the best model for the relationship among the input and output variables. The main advantage of EPR is that it does not need to consider a nonlinear regression model. The basis of this algorithm is that it combines the genetic algorithm (GA) and the regression analysis method [55]. EPR works based on the creation of several candidate relations between input and output data using an evolutionary process that utilizes the GA. The developed relations depend on the number of data points, the kind of relations among input and output variables, the suggested limit for constant exponents, and the number of suggested terms in the final relation. Equation (2) shows the general scheme of EPR [33]:  
\[y=\sum_{j=1}^{m}F(X,f(X),a_{i})+a_{0} \tag{2}\]  
where, \(y\) is the output vector, \(a_{i}\) is a constant value, \(X\) is the input variables matrix, and \(m\) is the number of the terms in suggested relation. Furthermore, \(F\) is a relation created by the process and \(f(X)\) is a function defined by user. To develop Equation (2), the GA is used. However, \(f(X)\) and m are determined by the user based on their understanding of the intended physical phenomenon or a trial and error method. The first step in finding the relation between the input and output data is rewriting Equation (2) as a vector, as shown in Equation (3) [53]:  
\[Y_{N\times 1}(\theta,Z)=\Big{[}I_{N\times 1}\ Z_{N\times M}^{j}\Big{]}[a_{0} \ a_{1}\ \ldots\ a_{m}]^{T}=Z_{N\times d}\times\theta_{d\times 1}^{T} \tag{3}\]  
where, \(Y_{N\times 1}(\theta,Z)\) is a vector for the estimation of the least squares for \(N\) target values, and \(\theta_{1\times d}\) is a vector consisting of \(a_{j}\) and \(a_{0}\). Finally, \(Z_{N\times d}\) is a matrix that is composed of the identity matrix for \(a_{0}\) and m number of \(Z^{j}\) variables. For a fixed value of j, \(Z^{j}\) is the product of vectors of independent variables [53]. Equation (4) presents the input data of \(X\)[53]:  
\[X= \begin{bmatrix}x_{11}&\cdots&x_{1k}\\ \vdots&\ddots&\vdots\\ x_{N1}&\cdots&x_{Nk}\end{bmatrix} =[X_{1}\ X_{2}\ X_{3}\ldots\ X_{k}] \tag{4}\]where the \(k\)th column of \(X\) represents the candidate variables for the \(j\)th term of Equation (3). Thus, \(Z_{N\times M}^{j}\) in Equation (3) can be written as Equation (5) [56]:  
\[Z_{N\times 1}=\left\lfloor(X_{1})^{ES(j1)}.(X_{2})^{ES(j2)}.\ \ldots.\ (X_{k})^{ES(jk)}\right\rfloor \tag{5}\]  
where \(Z^{j}\) is the vector of the \(j\)th column whose elements are a product of candidate independent inputs, ES is a matrix of exponents, and k is the number of independent variables. To develop the final equation, assuming that the vector of exponent constants defined by the user is \(EX=[0.5,\ 1,\ 2]\), the number of terms specified by m (without bias) is 4, the number of independent variables used in the analysis (\(k\)) is 3, the number of columns is 3, and the number of the lines of the ES matrix is 4. The coefficients of the ES exponent, for example, will be as follows:  
\[ES= \begin{bmatrix}0.5&1&2\\ 0.5&2&2\\ 1&2&0.5\\ 1&1&0.5\end{bmatrix}  \tag{6}\]  
Applying the matrix presented in Equation (6) to Equation (5), the four mathematical equations (Equations (7)-(10)) are:  
\[Z_{1}=(X_{1})^{0.5}.(X_{2})^{1}.(X_{3})^{2}=X_{1}^{0.5}.X_{2}.X_{3}^{2} \tag{7}\]  
\[Z_{2}=(X_{1})^{0.5}.(X_{2})^{2}.(X_{3})^{2}=X_{1}^{0.5}.X_{2}^{2}.X_{3}^{2} \tag{8}\]  
\[Z_{3}=(X_{1})^{1}.(X_{2})^{2}.(X_{3})^{0.5}=X_{1}.X_{2}^{2}.X_{3}^{0.5} \tag{9}\]  
\[Z_{4}=(X_{1})^{1}.(X_{2})^{1}.(X_{3})^{0.5}=X_{1}.X_{2}.X_{3}^{0.5} \tag{10}\]  
Therefore, Equation (3) becomes:  
\[Y=a_{0}+a_{1}.Z_{1}+a_{2}.Z_{2}+a_{3}.Z_{3}+a_{4}.Z_{4}\] \[=a_{0}+a_{1}.X_{1}^{0.5}.X_{2}.X_{3}^{2}+a_{2}.X_{1}^{0.5}.X_{2}^ {2}X_{3}^{2}+a_{3}.X_{1}.X_{2}^{2}.X_{3}^{0.5} \tag{11}\] \[+a_{4}.X_{1}.X_{2}.X_{3}^{0.5}\]  
In order to determine the constant coefficients of Equation (11), the least-squares method is used [53].