Comparison of data-driven prediction methods for comprehensive coke ratio of blast furnace

2.1 MLR

MLR attempts to establish the relationship between a response variable and two or more explanatory variables by fitting a linear equation to the observed data [20]. The value of the response variable y has a relationship with every value of the independent variable x. The population regression line for p explanatory variables x 1 , x 2 , … , x p is defined to be μ y = β 0 + β 1 x 1 + β 2 x 2 + ⋯ + β p x p . This line describes how the mean response μ y changes with the explanatory variables. The observed values for y vary about their means μ y and are assumed to have the same standard deviation σ. Formally, the MLR model for n observations can be represented as follows:

where ε is the notation for the model deviations.

2.2 CART

CART [21,22] is a supervised learning technique that is suitable for both classification and regression problems. It provides a flexible tree-like structure, where internal nodes, branches, and leaf nodes represent the features of a dataset, the decision rules, and the outcomes, respectively. Two kinds of nodes are in a decision tree, i.e., a decision node, including multiple branches for making decisions, and a leaf node for outputting a decision without further branches.

The principle of CART is described as follows. Without losing generality, X and Y are given as the input vector and the response vector, respectively. The training set is defined as D = { ( x 1 , y 1 ) , ( x 2 , y 2 ) , … , ( x N , y N ) } , where x i = ( x i ( 1 ) , x i ( 2 ) , … , x i ( n ) ) is known as feature vector, n represents the number of features, i = 1 , 2 , … , N , N is sample size.

The feature space is divided by the heuristic method. In each division, each value of the features in the current set will be investigated one by one, and the optimal one will be selected as the segmentation point according to the least-squares error criterion. The j-th feature of the training set is denoted as x ( j ) , the value of which is s. x ( j ) and s are regarded as segmentation variable and segmentation point, respectively. Suppose two regions: R 1 ( j , s ) = { x | x ( j ) ≤ s } and R 2 ( j , s ) = { x | x ( j ) > s } . Solve the following equation to obtain optimal j and s:

where c 1 and c 2 denote the fixed output values of the two regions. The above formula can also be written as follows:

where c ˆ 1 = 1 N 1 ∑ x i ∈ R 1 ( j , s ) y i and c ˆ 1 = 1 N 2 ∑ x i ∈ R 2 ( j , s ) y i .

After finding the optimal segmentation point (j, s), the input space is divided into two regions. Repeat the above division process for each region until the termination condition is reached. Based on the above process, a least-squares decision tree is constructed.

2.3 Lasso

Lasso regression [23,24], sometimes called L1 regularization of LR, is a modified form of LR. A regularization parameter as punishment is used to control model complexity. It is a biased estimation for dealing with complex collinearity data and obtains a more accurate model by constructing a penalty function. Suppose a regression function:

Its loss function can be expressed as follows:

where X ∈ R n × m is the input vector, Y ∈ R m is the response vector, n and m are the feature number and the sample size, respectively. The emergence of Lasso is to tackle two problems: over-fitting of LR and the occurrence of X transpose multiplied by irreversible X in course of solving θ by regular methods. To this end, a regularizer is introduced into the loss function, which is as follows:

where λ is a regularization parameter. If λ is set too big, it will cause an under-fitting result; otherwise, it will cause over-fitting.

In addition, the loss function can also be matrixed into the following form:

2.4 EN

EN, the mixture of Ridge and Lasso regressions, is a linear regression method with L1 and L2 norms as the priori regularizers [25]. A model with only a few non-zero and sparse parameters is learned by this combination, which is just like Lasso. But it still maintains some regular properties like Ridge regression. The matrix of the loss function of EN is described as follows:

2.5 KNN

KNN is one of the most basic and simplest machine learning algorithms [26]. Its idea is very simple: Each of the n-dimensional input vectors corresponds to a point in the feature space, and the output is a category label or a prediction value. When it is used for regression prediction, the average of those target values ( y i ) of k nearest samples is taken as the prediction value y ˆ of the new sample [27]. The formula is as follows:

2.6 SVR

Based on kernel functions and ε insensitive function, SVR first proposed by Cortes and Vapnik [28] has been applied in many fields successfully [29,30], such as character recognition [31], drug design [32], and combinatorial chemistry [33]. The input and output of the i-th sample ( i = 1 , 2 , … , l ) are denoted as x i ∈ R n and y i , respectively. The solution of the nonlinear SVR could be obtained via the following equation:

s.t. 0 ≤ α i ≤ C , i = 1 , 2 , … , l

where α i and α i ∗ represent Lagrange coefficients; C is the penalty factor; and K ( ⋅ ) denotes the kernel function. The regression function f ( x ) is

where N _SV denotes the number of SVs.

2.7 AdaBoost

In view of multiple variations of AdaBoost, its R ² regression is illustrated as an example here [34]. Denoting T = { ( x 1 , y 1 ) , ( x 2 , y 2 ) , … , ( x m , y m ) } as the sample set. The iteration number of the weak learner is K.

The initial weight of the sample set is expressed as follows:

Supporting k = 1 , 2 , … , K :

1. Gain a weak learner G k ( x ) by training the sample set with weights D k .

2. Calculate the maximum error of the training set.

   (12) E k = max | y i − G k ( x i ) | , i = 1 , 2 , … , m .

3. Calculate the relative error for each sample. Exponent error is expressed as follows:

   (13) e k i = 1 − exp − | y i − G k ( x i ) | E k .

4. Calculate the regression error rate as follows:

   (14) e k = ∑ i = 1 m ω k i e k i .

5. Calculate the coefficient of the weak learner as follows:

   (15) α k = e k 1 − e k .

6. Update the weight distribution of the sample set as follows:

where the normalization factor Z k is expressed as follows:

Finally, the strong learner is constructed as follows:

where g ( x ) , k = 1 , … , K , is the median value of all α k G k ( x ) .

2.8 RF

RF [35] is an ensemble of B trees { T 1 ( X ) , … , T B ( X ) } , where X = { x 1 , … , x p } is a p-dimensional input vector [36]. The ensemble produces B outputs { Y ˆ 1 = T 1 ( X ) , … , Y ˆ B = T B ( X ) } , where Y ˆ b , b = 1 , … , B , is a prediction value of the output variable generated by the b-th tree. The outputs of all trees are aggregated to generate one final prediction, Y ˆ . In regression, it is the average prediction of all the individual trees.

The training set, D = { ( X 1 , Y 1 ) , … , ( X n , Y n ) } , where X i , i = 1 , … , n , is an input vector and Y i is an output vector. The training process is as follows:

1. A bootstrap sample with replacement is randomly drawn from the training data including n samples.

2. For each bootstrap sample, a tree will be grown by the following process: at each node, the optimal segmentation is chosen among the subsets, including random m try (rather than all) variables. Here, m try is the only tunable parameter in RF. The tree is grown to its maximum size (i.e., until no further splits are possible) and not pruned back.

3. Repeat the above steps until B (a sufficiently large number) trees are grown.

4. Average the predictions of the individual trees to obtain the final prediction as follows:

2.9 KNN-AdaBoost, LR-AdaBoost, and RF-AdaBoost

Ensemble learning algorithms can significantly promote the generalization ability of learning systems by training a certain number of weak learners and combining them into a strong learner. KNN-AdaBoost, LR-AdaBoost, and RF-AdaBoost are all ensemble algorithms that improve AdaBoosts by training KNN, LR, or RF as weak learners, respectively [37,38,39]. The training process is the same as AdaBoost, and only one difference can be found in the first step. KNN, RF, or LR algorithms are selected as the weak learners G k ( x ) , and the sample set is trained with the weights D k in iterations.

2.10 GBR

GBR fits the negative gradient by iterative process [40,41]. Considering a regression question, n samples, ( x 1 , y 1 ) , ( x 2 , y 2 ) , … , ( x n , y n ) given are fitted into a function F ( x ) with minimal error.

Set the fitted values y ˆ = F ( x ) , square error function L ( y , y ˆ ) = ( y − y ˆ ) 2 / 2 , and the total error for all samples J = ∑ i = 1 n L ( y i , y ˆ i ) .

The target value y i of the sample is known, and the final result of error J is a numerical scalar. Error J is a function of the n-variable ( y ˆ 1 , y ˆ 2 , … , y ˆ n ) for a total of n samples.

Calculate the gradient of the function J as follows:

It can be seen from the above derivation that the residual error is equal to the negative gradient y − y ˆ = − ∇ J .

The weak models are established through iteration, gradually enhanced (boosting) and combined into a strong model.

2.11 ERT

ERT, proposed by Geurts et al. [42] in 2006, is very similar to the RF algorithm, which is composed of many decision trees [43,44]. The two major differences between the two algorithms are described as follows:

1. The bagging model is applied in the RF algorithm. While all samples are trained to generate each decision tree in ERT, i.e., each decision tree is constructed with the same training samples.

2. A random subset of RF needs to find the best bifurcation attribute. While the bifurcation value of ERT is obtained completely at random, so as to realize the bifurcation of a decision tree. The ERT model is trained with random features and threshold values. Therefore, the training process of ERT is much faster than RF.

4.1 Dataset

The dataset consists of two parts: the training set for modeling and the test set for validation. After excluding the outliers caused by the abnormal state of BF (such as blowing down, record fault, or overhaul) and data with missing values, the training set with 326 samples was established by collecting the historical data throughout a year from an industrial BF (an internal volume of 2,000 m³) with a sampling interval of 1 day. The test set includes 87 samples collected from the same BF in the first 3 months of next year. The disciplinarian of BF burden distribution remained the same during the period of data collection [17].

The dependent variable of the dataset is CCR (kg·t⁻¹). The data for every CCR is the average of a day. According to the experience of BF experts, 36 process variables (shown in Table 1) were chosen as the independent variables of the dataset. The quantity and quality of data are two key factors influencing industrial optimization and constructing a reasonable model. A common thought is that the amount of data needed to construct a plausible data-driven model is at least three times of feature number. The quantity of data in the work meets the above requirement. The quality of the data depends on the spatial coverage of the target variable and the uncertainties associated with the data [45]. Normally, data with a normal distribution is in favor of establishing a reasonable model. Data uncertainty, such as experimental error and input error, could influence the quality of the data. To further analyze the quality of the data, the distribution curves of CCR and the independent variables were plotted in Figures S1 and S2 of Supporting Information, respectively. It can be seen from the two figures that all variables are basically consistent with normal distribution, and the stability of the data is higher because the data show relative concentration. Most of the data measurement instruments have random and systematic errors, so the data collected from the BF system has some degree of uncertainty. On the whole, the quality of the data set used is suitable for data mining and industrial optimization, although the BF data have some inevitable errors [17,19].

4.2 Selecting variables

Before establishing the models, variable selection can not only reduce the dimension of feature space to further decrease the risk of over fitting but also better remove some variables unrelated to the target variable and noise interference. Meanwhile, it can also greatly reduce training time and further improve the prediction accuracy and generalization performance of the model [17,46]. Apparently, the model’s performance will be harshly affected, and the model is too complex if the variables irrelevant to the output are not deleted before training a model.

In the work, the primary feature screening was conducted using the Pearson correlation coefficients between the features. First, the max relevance min redundancy (mRMR) method [47] was employed to sort the features. The sorting result is represented in Figure S3. Second, the Pearson correlation coefficients between the features were calculated. Its matrix is shown in Figure S4. If the correlation score between the features is higher than 0.9, the feature with the lower mRMR score will be omitted. After PI_Si, R _CL, and SS_B were omitted, the 33 features were retained by the primary screening.

GA based on the learners and RFE based on the base_estimators were used to further optimize the feature sets, respectively. GA, first put forward by Holland [48], is a global optimization algorithm of random search that simulates the processes of inheritance and evolution in natural conditions. It has highly collateral, random, and adaptive features, is an effective method to remove from local optima present on the response surface, and can solve a wide variety of optimization problems with the requirement of no knowledge about the response surface or gradient present on it [49].

RFE is a method that works by selecting data features recursively based on the smallest feature value [50]. In the RFE concept, RFE based on the base_estimators works by eliminating irrelevant features in each iteration, namely the lowest-weight feature. This method is divided into three stages as described in ref. [51]. The two feature selection processes are related to the learners or the base_estimators and are automatically performed during the training of the learners or the base_estimators.

Figure 2 illustrates the implementation processes of variable selection of GA based on the seven learners, respectively. The learners employed are CART, KNN, SVR, KNN-AdaBoost, GBR, LR-AdaBoost, and RF-AdaBoost. Genetic algebra and the best feature set for each GA-learner are listed in Table 2. Figure 3 shows the processes of feature selection of RFE based on the six base_estimators, namely, MLR, Lasso, EN, AdaBoost, RF, and ETR. Table 3 illustrates the best feature set for each RFE – base_estimator. In the variable screening process, the optimal feature sets will be found when R is the maximum or RMSE is the minimum. It can be seen from Tables 2 and 3 that the best feature set is different with different learners or base_estimators, and the feature number for different optimal feature set varies greatly. The least feature number is 4, and the most is 18.

Figure 2

The procedures of variable selection in GA-learners: (a) GA-CART, (b) GA-KNN, (c) GA-SVR, (d) GA-KNN-AdaBoost, (e) GA-GBR, (f) GA-LR-AdaBoost, and (g) GA-RF-AdaBoost.

Figure 3

The procedures of feature selection in RFE based on the base_estimators: (a) MLR, (b) Lasso, (c) EN, (d) AdaBoost, (e) RF, and (f) ETR.

The top three most frequency features selected by the various GA-learners are P, U _IO, and Y _Coal in the feature subsets. P refers to the amount of ore smelted by unit coke, whose higher value may result in lower CCR. U _IO is the amount of iron ore consumption for generating 1 t of pig iron, and the larger U _IO content may lead to the higher CCR. In terms of the mechanism of BF operation, Y _Coal has large contributions to CCR and their relationship is related to the quality and dosage of pulverized coal and coke.

As shown in Figure 4, the scatter plots between the top important features and CCR are generated for further analysis based on the training data. Figure 4(a) indicates that the larger P tends to result in a smaller CCR value and has a negative correlation with CCR. Meanwhile, Figure 4(b) shows the positive correlation between U _IO and CCR, and Figure 4(c) reveals the positive one. Therefore, the relationships between the three features and CCR are consistent with the mechanism of BF operation.

Figure 4

Relationships between the top three most important features: (a) P, (b) U _IO, and (c) Y _Coal and CCR.

4.3 CCR prediction

The parameters of the models were optimized with the grid search technique and TFCV. The results were compared by three evaluation criteria, namely RMSECV, MAE_CV, and R ²_CV, where “CV” represents “cross validation.” The three methods with the top performances were discovered by TFCV and used to establish the models. The generalization capabilities and robustness of the models were evaluated by using RMSE, MAE, R, and R ².

4.3.1 Model validation

The TFCV results of six conventional and seven ensemble methods constructed with the optimal features obtained from feature selection are shown in Table 4. Compared with the other methods, KNN-AdaBoost has an inferior prediction performance, for it achieved the largest RMSECV and MAE_CV and the lowest R ²_CV in TFCV. MLR, SVR based on RBF, and LR-AdaBoost exhibit the best predictive properties among all the technologies. LR-AdaBoost is an ensemble learning algorithm that improves AdaBoost with linear regression as weak learner. According to the excellent performance of the MLR and LR-AdaBoost models, a conclusion can be drawn that a certain linear relationship exists between the target variable and the independent variables. Due to sufficiently high R ²_CV and very low RMSECV and MAE_CV, the three models have excellent generalization performance. Therefore, the further optimization of other methods would not be considered.

Table 4

The TFCV results of the thirteen methods

Methods	RMSECV	MAE_CV	R 2_CV
MLR	1.668	1.079	0.973
CART	5.630	3.860	0.757
Lasso	3.018	2.365	0.928
EN	4.776	3.635	0.819
KNN	4.233	2.890	0.882
SVR	1.878	1.158	0.975
AdaBoost	5.073	3.808	0.803
RF	3.711	2.741	0.885
KNN-AdaBoost	6.404	4.882	0.674
GBR	3.459	2.578	0.904
LR-AdaBoost	2.276	1.638	0.953
ERT	4.053	2.976	0.867
RF-AdaBoost	4.226	2.931	0.858

Figure 5 reveals the box plots of RMSECV and R ²_CV for six conventional methods. In Figure 5(a), the technologies are ranked from low to high according to the center of RMSECV: MLR, SVR, Lasso, KNN, EN, and CART. Especially, the center values of RMSECV distribution of the MLR and SVR models are far less than other methods. As shown in Figure 5(b), the zones of R ²_CV for MLR, SVR, and Lasso are much narrower than those of the other methods, and their center values are much higher than those of the others. According to the above analysis, MLR, SVR, and Lasso are more reliable and more precise for predicting CCR than the other techniques.

Figure 5

Box plots of (a) RMSECV and (b) R ²_CV for six conventional methods.

Figure 6 represents the box plots of RMSECV and R ²_CV of TFCV for seven ensemble methods. In Figure 6(a), the methods are ordered according to the central values of RMSECV in ascending order: LR-AdaBoost, GBR, RF, ERT, RF-AdaBoost, AdaBoost, and KNN-AdaBoost. The RMSECV range of LR-AdaBoost is the narrowest among the seven methods. In Figure 6(b), the centers of LR-AdaBoost, RF, GBR, ERT, and RF-AdaBoost are much higher than those of the others. The methods with a narrow distribution of R ²_CV include LR-AdaBoost, RF, and ERT. On the whole, LR-AdaBoost is the most accurate and stable technique for predicting CCR among all the ensemble techniques.

Figure 6

Box plots of (a) RMSECV and (b) R ²_CV for seven ensemble methods.

In addition, it can be seen from the above two figures that MLR has the superior prediction performance, while SVR and LR-AdaBoost trail in second and third, respectively.

4.3.2 Model training and testing

Based on the above comparison and analysis of the thirteen methods in TFCV, MLR, SVR, and LR-AdaBoost have the top three performances for CCR prediction. This section will illustrate how they would be employed to establish the models and how well the models performed in external validation.

By setting a regression equation, the MLR model can estimate the dependent variable with multiple independent variables and then predict the dependent variable. The independent variables were standardized because of their different units. The larger the value of the standardized regression coefficient, the greater the influence of the independent variable on the dependent variable. The MLR model constructed with the optimal feature set is shown as follows:

(26) CCR = − 4.16 PI Ti + 9.05 P Top − 1.64 R Ox + 0.82 M 10 − 0.19 C 25 + 0.24 Ic − 0.05 SS + 287.4 SC − 1.88 FB − 0.41 R Gas + 0.25 U IO − 9.05 S B − 8.76 P 1B − 9.5 P 2B − 9.61 M B + 9.25 C BW + 0.21 Y Coal + 348.52 .

Parameter tuning of SVR and LR-AdaBoost is a crucial step to develop models with high generalization performance. The three hyperparameters, C, ε, and σ need to be optimized using the grid search method based on leave-one-out cross-validation before SVR modeling [17]. In the work, the range of ε was from 0.01 to 0.1 with a step of 0.01; C changed from 1 to 100 with an interval of 2; σ varied from 0.5 to 1.4 with a step of 0.1. As can be seen from Figure 7, the higher the C and σ, the larger the RMSE. The RMSE value ranges from 1.867 to 5.470. Obviously, the reasonable selection of C, ε, and σ has a great influence on the prediction performance of the SVR model. After conducting the grid search shown in Figure 7, the optimal C, ε, and σ were determined to be 23, 0.01, and 0.7, respectively, when RMSE (equal to 1.867) was lowest.

Figure 7

The variety of RMSE with σ and C in optimizing SVR’s hyperparameters.

The SVR model constructed with the optimal C, ε, and σ is shown as follows:

(27) y = ∑ i = 1 n β i ⋅ exp ( − 0.7 ⋅ | | x − x i | | 2 ) + 0 .64286,

where x is the unknown vector, x i is the support vector of the SVR model, n is the corresponding sample number, and β _i is the Lagrange multiplier of the support vector.

The two hyperparameters of the LR-AdaBoost model, namely n_estimators and learning_rate need to be tuned by using the grid search method. Here, n_estimators is the number of weak learners. Generally, too small n_estimators could cause under-fitting, otherwise could cause over-fitting. Furthermore, learning_rate is the weight reduction factor of the weak learner. In this work, n_estimators varied from 1 to 10 with a step of 1 and from 20 to 100 with a step of 10. Learning_rate changed from 0.1 to 1 with an interval of 0.1. The result indicated that n_estimators and learning_rate were equal to 2 and 0.1, respectively, when RMSE was at a minimum. The optimization process of the parameters is shown in Figure 8, from which it can be found that n_estimators have a more influential impact on RMSE than learning_rate. When n_estimators are less than 20, a smaller RMSE can be obtained.

Figure 8

The variety of RMSE with learning_rate and n_estimators in optimizing LR-AdaBoost’s hyperparameters.

The MLR, SVR, and LR-AdaBoost models were employed to predict the samples in the test set. The evaluation criteria of the three models, RMSE, MAE, R, and R ² are listed in Table 5. The subscripts “TR” and “TE” represent “training” and “test,” respectively. It can be seen from the table that MLR has the lowest RMSE_TR and MAE_TR, and the highest R _TR and R ² _TR. Moreover, it also performs best in the external test. In model training, SVR is superior to LR-AdaBoost. The reverse is true for the performance indicators generated by the external test of the two technologies. This shows that LR-AdaBoost has stronger extrapolation performance than SVR. Figure 9 represents the scatter plots of the three models. It can be seen from the figure that they could primely fit the experimental data since the predictive data are closely distributed along the diagonal direction. In particular, the training result of SVR is slightly better than LR-AdaBoost, but LR-AdaBoost in external validation is more reliable than SVR. It can be seen from the above analysis that the three methods all possess high prediction accuracy and practical values.

Table 5

The prediction results of the MLR, SVR, and LR-AdaBoost models

Indicator	MLR	SVR	LR-AdaBoost
RMSETR	1.337	1.812	1.855
MAETR	1.037	1.084	1.123
R TR	0.993	0.988	0.986
R 2 TR	0.987	0.976	0.973
RMSETE	1.593	2.464	1.863
MAETE	1.328	1.887	1.508
R TE	0.989	0.981	0.987

Figure 9

The scatter plots of the (a) MLR, (b) SVR, and (c) LR-AdaBoost models.