Demand forecasting of green metal materials using non-equidistant grey prediction with robust nonlinear interval regression analysis.

The main advantages of magnesium alloys are that they are lightweight, easy to recycle, and have high vibration absorption. These unique characteristics make magnesium alloys important green metal materials for manufacturing, especially for the automotive and 3C products industries. The developing trends of these related industries can be recognized by forecasting the demand for magnesium alloys. This study develops grey prediction power models to forecast the demand for such a promising green metal material. Grey prediction is an appropriate technique because available data regarding the demand for magnesium alloys are not in line with any statistical assumptions. In particular, because outliers might cause a deterioration of forecasting performance, a robust nonlinear interval regression analysis with neural networks is applied to detect outliers by estimating data intervals. Then, a power model is applied to the newly generated non-equidistant data sequence without outliers. Residual modification is further considered here to improve the forecasting performance of the power model. The forecasting abilities of the proposed grey residual modification models are verified using actual magnesium alloy demand data. The experimental results for ex-post testing show that the mean absolute percentage errors of the proposed models that can work on non-equidistant data were minimal among all considered models.

Introduction

Magnesium is a lightweight engineering material, and magnesium alloys are regarded as a twenty-first century green metal material. Compared with other metal materials, magnesium alloys have a few advantages: they are lightweight, have high impact resistance, high corrosion resistance, and are easy to recycle. These advantages are the reason why magnesium alloys have been widely used for the development and manufacture of various products, including the components of automobiles, aircraft, bicycles, sports equipment, and the metal cases of computer, consumer electronic, and computer (3C) products (Hong, 2016; Lu, 2014). In 2016, as a result of major demand from the automotive industry, 32% of the worldwide magnesium was consumed in the die casting of magnesium alloys (Lee et al., 2017). Next, as the main determinant of fuel consumption is vehicle weight, lightweight components are essential for the automotive and aerospace industries. As national authorities make efforts to significantly mitigate the effects of CO2 emissions and promote waste recycling, using magnesium materials in the manufacture of products, such as vehicles and 3C products, can have a very positive influence on environmental protection and energy conservation. Therefore, magnesium alloy demand forecasting can help these industries understand developing trends to quickly recognize market configurations and set up investment strategies.

In contrast to artificial intelligence techniques and time series models, grey system theory has been applied to many real-world problems (Gao et al., 2021; Rajesh, 2020). Grey prediction has drawn attention because of the need to use limited samples to characterize an unknown system without complying to any statistical assumptions (Liu & Lin, 2010). This characteristic motivated us to apply a power model of a first-order grey model with one variable (GM(1,1)), abbreviated to GM(1,1)-P, to predict the demand trends for magnesium alloys. Indeed, the GM(1,1)-P model has been widely applied to real-world problems. The literature on the applications of GM(1,1)-P is briefly described in Table 1.

Table 1

Summary of the literature on the applications of GM(1,1)-P

Author(s)	Application
Wang et al. (2011)	Qualified discharge rate of industrial wastewater in China
Pao et al. (2012)	CO2 emissions, energy consumption and economic growth in China
Lu et al. (2016)	Traffic flow forecasting
Chen et al. (2019)	Short-term forecast of the passenger volume of Taiwan High Speed Rail
Wu et al. (2019)	Short-term renewable energy consumption
Xiao et al. (2020)	Biomass energy consumption
Ma et al. (2019)	Tourist income prediction
Şahin (2020)	Cumulative number of confirmed cases of COVID-19
Zheng et al. (2020)	Hydroelectricity consumption of China

In light of the applicability of the GM(1,1)-P model to real-world problems, it is interesting to apply the GM(1,1)-P and its variants to forecast the magnesium alloy demand. Two issues are particularly addressed in this paper. First, GM(1,1)-P often works well on annual or equidistant data; however, there is always the possibility that outliers might worsen the performance of the prediction model (Huang et al., 1998). To avoid outliers adversely affecting forecasting performance, robust interval regression analysis can be used to detect outliers by gauging the intervals consisting of upper and lower bounds (Hladík & Černý, 2014; Wang et al., 2015). Thus, it is interesting to develop the GM(1,1)-P that can work on a newly generated non-equidistant data sequence after removing those outliers to better forecast magnesium alloy demand.

Next, it has been demonstrated that residual modification is available to promote forecasting accuracy of an original prediction model (Liu et al., 2017). Therefore, it is reasonable to incorporate the mechanism of residual modification into GM(1,1)-P. Based on these features, this study aims to develop GM(1,1)-P and its residual version that works on non-equidistant data sequences, generated using robust nonlinear interval models to remove outliers.

To deal with the first aforementioned issue, we use two neural networks (NNs) to build a nonlinear interval model that is robust against outliers. This approach is adopted here because several studies have demonstrated the effectiveness of using two NNs on the construction of nonlinear interval models. The related literature is summarized in Table 2. As for the second issue with respect to residual modification, the functional-link net (FLN) is considered because Hu (2017) and Hu et al. (2019) indicate the effectiveness of the FLN on the combination of a grey prediction model with its residual version. The main advantage of FLNs is relieving the restrictions arising from commonly used sign estimation methods, such as Markov-chain estimation (Hsu, 2003) and genetic programming (Lee & Tong, 2011).

Table 2

Summary of the literature on the use of NNs to build nonlinear interval models

Author(s)	Subject
Cheng and Lee (2001)	Fuzzy regression with radial basis function network
Chen and Jain (1994)	A robust back propagation learning algorithm
Huang et al. (1998)	Robust interval regression analysis using multi-layer perceptrons (MLPs)
Ishibuchi and Tanaka (1992)	Fuzzy regression analysis using MLPs
Ishibuchi et al. (1993)	Fuzzy regression analysis using NNs with interval weights
Ishibuchi and Nii (2001)	Fuzzy regression using fuzzified neural networks
Nasrabadi and Hashemi (2008)	Robust fuzzy regression analysis using NNs
Hao (2009)	Interval regression analysis using support vector networks
Hu (2009)	Robust nonlinear interval regression analysis using functional-link nets
Hu (2014)	Robust nonlinear interval regression analysis using MLPs

The remainder of this paper is organized as follows. The robust nonlinear interval regression analysis using NNs and grey prediction are introduced in Sects. 2 and 3, respectively. The proposed GM(1,1)-P along with its residual model that work on non-equidistant sequences are presented in Sect. 3. For forecasting magnesium alloy demand, we investigate the forecasting accuracy of the different grey prediction models in Sect. 4. Finally, the conclusions are presented in Sect. 5.

Robust nonlinear interval regression analysis using NNs

Imprecise observations can be represented by estimated data intervals (Xie et al., 2014; Zeng et al., 2014). Nonetheless, an unsatisfactory interval can be obtained when a data sequence contain outliers. In Sect. 2.1, Eqs. (1), (2), (3), (4), (5), and (6) describe nonlinear interval regression analysis using NNs, and Eqs. (7), (8), (9), (10), and (11) in Sect. 2.2 describe the development of the robust nonlinear interval models against outliers.

Nonlinear interval regression analysis

In practice, a nonlinear interval model can be determined using two NNs, denoted by NN and NN, from a given original data sequence x = (x1, x2,…, x). As a result, (t1, x1), (t2, x2),…, and (t, x) comprise a data set for the two NNs (NN and NN) where (t, x) denotes the p-th input–output pattern at time t. Also, let and be the upper and lower limits obtained from NN and NN, respectively, such that the difference between and denotes the width of the gauged interval. A nonlinear optimization problem is formulated as.

subject to

To approximately satisfy the constraint condition, Ishibuchi and Tanaka (1992) proposed an NN-based approach in which two MLPs are used to set up the nonlinear interval model, and each MLP has a single output, a single input, and one hidden layer with five hidden units. Note that what the NN this study addresses is MLP. With a weighting parameter , a cost function E is formulated to discover the upper limit:

in whichwhere ω ∈ (0, 1). To determine the lower limit, cost function E is formulated as

in which

Robust nonlinear interval models

A robust nonlinear interval model can be built on the basis of the aforementioned NN-based approach, where the quality index, , with respect to (t, ) at time t beyond the estimated interval should be considered to help the MLP training resist outliers. The quality index is given as follows:where c is a cutoff point, and is a very small positive real value. During the model construction, r > c means that (t, x) falls in the rejection region. To determine the upper limit against outliers, the absolute residual, r, with respect to (t, x) is defined as follows (Chen & Jain, 1994; Huang et al., 1998):

In this approach, data in the rejection region impose little effect on the construction of NN.

To determine c, we renumber the data by sorting r in descending order, such that (t, ) is renumbered as one as it has min{r1, r2,…, r}. is assigned to c, where q (0 ≤ q ≤ 1) indicates the tolerable limit on the rate of outliers beyond the true interval. With , cost function E is defined as (Hu, 2009, 2014)

To determine the lower limit, r is defined as

In a similar manner, with q being the tolerable limit on the rate of outliers beneath the true interval, c is given by as well. With , cost function E is defined as The data in the rejection region also impose little effect on the construction of NN.

By defining the fitness function as E (E), the connection weights of NN (NN) can be determined using the toolbox in MATLAB to implement a real-valued genetic algorithm (GA). The estimated data intervals can be obtained from the best chromosome whose fitness value is a minimum among the straight generations. As a result, any data outside the data intervals can be treated as outliers. The details of applying MLPs to robust nonlinear interval regression analysis using GAs are omitted here for simplicity, but the reader is referred to Hu (2014) for further information.

Grey prediction

Equidistant grey prediction

GM(1,1)-P and its residual model, both of which use equidistant data sequences are described here prior to the introduction of the grey prediction models that use non-equidistant data sequences. Section 3.1 including Eqs. (12)-(20) introduces the traditional GM(1,1)-P first. Then, the residual modification model, NGM(1,1)-P, consisting of GM(1,1)-P and its residual version, is described in Sect. 3.2 including Eqs. (21)-(25).

GM(1,1)-P model

Given an equidistant data sequence  = (,,…,) composed of n data samples, a new sequence  = (,,…,) is generated from by a one-order accumulated generating operation (1-AGO) as (Liu & Lin, 2010)

Since monotonically increases, the following first-order Bernoulli equation can be used as the mathematical form of the GM(1,1)-P aswhere a denotes the developing coefficient, b represents the control variable, and r represents the power exponent.

The time response function with respect to Eq. (13) is (Chen, 2008; Chen et al., 2008)where r ≠ 1. Then, a and b can be estimated by a linear regression model consisting of n–1 grey difference equations with the following form:where the background value is given by

in which 0 ≤  ≤ 1. The interpolation coefficient is usually set to 0.5 for convenience, although such a setting may not be optimal.

In turn, with the ordinary least squares (OLS) method, a and b can be derived as

in which

Finally, the predicted value of can be obtained by the inverse 1-AGO (1-IAGO):

To optimally construct GM(1,1)-P, we apply the well-known optimization modeling software, Linear Interactive and General Optimizer (LINGO), to solve the optimization problem by minimizing the mean absolute percentage error (MAPE) as

subject to r ≠ 1,  = , Eq. (3), and Eqs. (5), (6), (7) and (8).

That is, MAPE is the sum of the APE of x (k = 1, 2,…, n) divided by the sample size. The computational steps of setting up the optimized GM(1,1)-P model are as follows.

Step 1 Generate accumulated generating sequence from the equidistant data sequence.

Step 2 Construct GM(1,1)-P with LINGO to obtain the optimal parameters, including the developing coefficient (a), control variable (b), interpolation coefficient (), and power exponent (r).

Step 3 Obtain the predicted sequence with 1-IAGO.

It should be noted that GM(1,1)-P, also called the nonlinear grey Bernoulli model, is identical to GM(1,1) when r = 0. In particular, the power exponent makes GM(1,1)-P much more useful than GM(1,1) (Dang et al., 2016; Lu et al., 2016).

Residual modification of GM(1,1)-P model

After constructing the GM(1,1)-P model, a sequence of absolute residual values,  = (, ,…, ), is generated where

A residual model is established with respect to as GM(1,1)-P using such that the corresponding time response function is

In common with GM(1,1)-P, a and b can be derived by the OLS method. A new predicted value can be obtained by modifying the predicted value produced by GM(1,1)-P as

The resultant prediction model consisting of GM(1,1)-P and its residual version is called RGM(1,1)-P.

Traditionally, s is 1 or − 1 (Liu & Lin, 2010). Even so, to optimally determine s, since it has been indicated that an FLN without hidden layers is effective for building a residual modification model (Hu et al., 2019), we set up an FLN to better integrate GM(1,1)-P and its residual version. Let t denote the time k with respect to . The input pattern for the FLN is an enhanced pattern, v = (t, sin(πt), cos(πt), sin(2πt), cos(2πt)). Let w = (w1, w2, w3, w4, w5) be a vector of connection weights. Then, the actual output of v can be obtained as:where θ is the bias and ˗1 ≤ y ≤ 1. As a result, a new predicted value generated from the proposed model can be computed as follows:

in which  > 0 and s =  y. Thus, the new value of ranges from to .

To set up the RGM(1,1)-P, we use the toolbox in MATLAB to develop a real-valued GA, using MAPE as the fitness function, by which the connection weights of the FLN can be optimally determined. The relevant parameters can be obtained from the best chromosome whose fitness value is a minimum among all the straight generations. As a result, the computational steps to construct the RGM(1,1)-P are briefly summarized as follows.

Step 1 Build the optimized GM(1,1)-P using the equidistant data sequence with LINGO.

Step 2 In common with the construction of the optimized GM(1,1)-P model, use the equidistant residual sequence to establish the optimized residual GM(1,1)-P model with LINGO.

Step 3 Obtain the predicted sequence using an FLN to optimally combine the predicted values obtained from the optimized GM(1,1)-P and its optimized residual model. The optimized FLN is derived with the GA toolbox in MATLAB.

Non-equidistant grey prediction

After removing outliers with the robust nonlinear interval model, the original equidistant data sequence becomes non-equidistant (Li et al., 2012). We use NGM(1,1)-P to denote a GM(1,1)-P model that can work on non-equidistant data. In our design, we first construct NGM(1,1)-P and then the residual NGM(1,1)-P. Therefore, NGM(1,1)-P is first introduced in Sect. 4.1 including Eqs. (26), Eq. (27) and Eq. (28). Generation of the proposed residual combination model, called RNGM(1,1)-P, consisting of NGM(1,1)-P and its residual version, is then described in Sect. 4.2 including Eqs. (29), (30), (31), and (32).

NGM(1,1)-P model

A non-equidistant sequence, (x(0)(t1), x(0)(t2),…, x(0)(t)) (n ≤ m), is obtained by removing outliers from the original sequence, where x(0)(t) denotes the actual value at time t (i ≤ n). Let , where i = 2, 3,…, n. Of course, is not a constant. Let  = (,,…,) = (x(0)(t1), x(0)(t2),…, x(0)(t)). There are four main differences between GM(1,1)-P and NGM(1,1)-P.

First, GM(1,1)-P and NGM(1,1)-P are set up by equidistant and non-equidistant sequences, respectively. Next, instead of 1-AGO,  = (x(1)(t1), x(1)(t2),…, x(1)(t)) = (,,…,) for NGM(1,1)-P can be generated by the following AGO as (Dang et al., 2016; Liu et al., 2017)

Then, the time response function obtained by solving the following whitening equation iswhere r ≠ 1 and k = 2, 3,…, n. Finally, with the following IAGO, the predicted value for NGM(1,1)-P with respect to is given by

Being analogous with the construction of GM(1,1)-P, the computational steps of setting up the optimized NGM(1,1)-P model are as follows.

Step 1 Remove outliers from the equidistant data sequence using the robust nonlinear interval model to detect outliers.

Step 2 Generate the accumulated generating sequence from the non-equidistant data sequence.

Step 3 Construct NGM(1,1)-P with LINGO to obtain the optimal parameters, including the developing coefficient, control variable, and interpolation coefficient.

Step 4 Obtain the predicted sequence with IAGO.

Residual modification of NGM(1,1)-P model

We set up the residual model of NGM(1,1)-P using the absolute residual value, , which is defined as

The prediction model consisting of NGM(1,1)-P and its residual version is called RNGM(1,1)-P. There are four main differences between RGM(1,1)-P and RNGM(1,1)-P as well. First, RGM(1,1)-P and RNGM(1,1)-P are set up using equidistant and non-equidistant residual sequences, respectively. Next, RNGM(1,1)-P follows the specific AGO where

Then, the corresponding time response function iswhere r ≠ 1. In common with NGM(1,1)-P, a and b can be derived by the OLS method. Lastly, with the specific IAGO, the predicted value for RNGM(1,1)-P can be generated by

Figure 1 depicts a flowchart of the construction the proposed RNGM(1,1)-P model. The computational steps to set up the RNGM(1,1)-P are briefly summarized as follows.

Fig. 1

Flowchart for constructing the RNGM (1,1)-P model

Step 1 Build the optimized NGM(1,1)-P using the non-equidistant data sequence with LINGO.

Step 2 In common with the construction of an optimized NGM(1,1)-P model, use the non-equidistant residual sequence to establish the optimized residual NGM(1,1)-P model with LINGO.

Step 3 Obtain the predicted sequence with an FLN to optimally combine the predicted values obtained from the optimized NGM(1,1)-P and its optimized residual model. The GA toolbox in MATLAB is applied to optimize the connection weights of the FLN.

Experiments

Experiments were conducted for magnesium alloy demand forecasting to make comparisons between the different grey prediction models using the equidistant sequences and those using the non-equidistant sequences. The former include the optimized GM(1,1) and optimized GM(1,1)-P, and the latter includes the proposed NGM(1,1), NGM(1,1)-P, and RNGM(1,1)-P. It is reasonable to include GM(1,1) with equidistant sequences and NGM(1,1) with non-equidistant sequences because GM(1,1) is the most commonly used grey prediction model (Liu et al., 2017). Several distinctive prediction models, the adaptive NGM(1,1) (ANGM(1,1)) (Chang et al., 2014) and FLNGM(1,1) (Hu (2017)), were included as well. It should be noted that, when a grey prediction model is developed on the basis of GM(1,1), “GM(1,1)” can be incorporated into the corresponding acronym by convention. The considered prediction models include:

NGM(1,1): A non-equidistant version of GM(1,1);

GM(1,1)-P: A power model of the GM(1,1);

NGM(1,1)-P: A non-equidistant version of GM(1,1)-P;

RGM(1,1)-P: A power model that integrates GM(1,1)-P with residual modification;

RNGM(1,1)-P: A power model that integrates NGM(1,1)-P with residual modification;

ANGM(1,1): Adaptive NGM(1,1) proposed by Chang et al. (2014) has shown its superiority over some prediction models for non-equidistant data sequences on some real cases including synthetic control charts and color filter manufacturing.

FLNGM(1,1): FLNGM(1,1) was originally developed for equidistant data, which uses an FLN to optimally integrate the optimized GM(1,1) and its residual model. For simplicity, the details of FLNGM(1,1) are omitted here, but the study of Hu (2017) can be consulted for more information.

NFLNGM(1,1): NFLNGM(1,1) denotes a non-equidistant FLNGM(1,1);

NN: The NN we use is a multi-layer perceptron (MLP) that has a single output, a single input, and one hidden layer with five hidden units.

As for parameter specifications involving the construction of the proposed model, when using the GA to set up the MLPs introduced in Sect. 2.2, we assigned values of 1000, 50, 0.8, and 0.01 to the total generations, population size, crossover, and mutation probabilities, respectively. In addition, ω in Eq. (4) and in Eq. (7) are equal to 0.2 and 10–5, respectively. To construct the proposed model, in Eq. (25) was specified as 20. Although several factors could affect the performance of the GA, there are no optimal values that can be assigned to its parameter specifications (Osyczka, 2003).

Case I: Magnesium die castings demand of China

China has consumed most of the magnesium in the world since 1999, and it has provided at least 80% of all magnesium materials used worldwide since 2008. With a demand of 70% from the automotive industry and 20% from the 3C products industry, about 35% of the magnesium has been consumed in die casting production with magnesium alloys (Lee et al., 2017). Magnesium die casting demand forecasting can thus reflect magnesium alloy demand and can indicate the developing trends of the automotive and 3C products industries in China.

We used the historical data of the demand for magnesium die casting that was reported and analyzed by Hsu (2012), Hong (2015a, 2015b), and Lee et al. (2017). Data from 2002 to 2012 were used for model fitting, and data from 2013 and 2014 were used for ex-post testing. The upper and lower limits are shown in Fig. 2 after performing the robust nonlinear interval regression analysis with q = q = 0.1. There are two outliers outside of the identified data intervals. The optimized NGM(1,1), NFLNGM(1,1), AGM(1,1), optimized NGM(1,1)-P, and RNGM(1,1)-P were constructed using the non-equidistant sequence by removing data in 2007 and 2009. The forecasting results are described in Table 3. Table 4 shows the forecasting accuracy obtained by the NN, GM(1,1), FLNGM(1,1), GM(1,1)-P, and RGM(1,1)-P that worked on the equidistant data sequence.

Fig. 2

Outlier detection of training data for Case I

Table 3

MAPEs of different grey prediction models working on non-equidistant data for Case I (unit: 1000 tons)

Year	Actual	NGM(1,1)		NFLNGM(1,1)		ANGM(1,1)		NGM(1,1)-P		RNGM(1,1)-P
		Predicted	APE	Predicted	APE	Predicted	APE	Predicted	APE	Predicted	APE
2002	0.61	0.61	0.00	0.61	0.00	0.61	0.00	0.61	0.00	0.61	0.00
2003	1.02	2.42	136.89	1.07	4.87	2.40	135.54	0.97	4.95	1.02	0.02
2004	1.8	2.82	56.67	1.43	20.75	2.82	56.76	1.80	0.14	1.80	0.24
2005	2.29	3.29	43.73	2.67	16.80	3.31	44.71	2.90	26.43	2.90	26.81
2006	5.1	3.84	24.68	4.51	11.66	3.89	23.69	4.16	18.49	4.18	18.02
2007	9.2
2008	6.85	4.86	29.08	7.17	4.65	4.97	27.46	6.08	11.21	6.16	10.08
2009	5.77
2010	7.66	6.62	13.61	8.61	12.45	6.85	10.52	8.17	6.66	8.33	8.71
2011	9.2	8.32	9.58	9.45	2.68	8.70	5.48	9.14	0.71	9.35	1.59
2012	9.73	9.71	0.22	9.22	5.19	10.21	4.96	9.42	3.16	9.68	0.52
2013	10.46	11.33	8.33	9.68	7.44	11.99	14.66	9.47	9.51	9.77	6.62
2014	11.27	13.23	17.35	11.43	1.39	14.09	24.98	9.30	17.48	9.65	14.41
2015	9	15.44	71.50	13.73	52.59	16.54	83.80	8.97	0.34	9.36	3.96
2016	10.4	18.01	73.22	16.43	57.94	19.43	86.81	8.51	18.14	8.94	14.07

“Predicted” refers to predicted values, and “APE” denotes the absolute percentage error

Table 4

MAPEs of different prediction models working on equidistant data for Case I (unit: 1000 tons)

Year	Actual	NN		GM(1,1)		FLNGM(1,1)		GM(1,1)-P		RGM(1,1)-P
		Predicted	APE	Predicted	APE	Predicted	APE	Predicted	APE	Predicted	APE
2002	0.61	0.35	42.24	0.61	0.00	0.61	0.00	0.61	0.00	0.61	0.00
2003	1.02	0.92	9.82	2.68	163.23	1.05	2.64	1.01	0.99	1.02	0.05
2004	1.8	2.05	13.73	3.08	71.22	1.73	4.03	1.80	0.04	1.79	0.39
2005	2.29	3.62	58.03	3.54	54.49	3.24	41.30	2.77	20.94	2.75	19.93
2006	5.1	5.20	1.97	4.06	20.37	4.75	6.90	3.84	24.72	3.81	25.36
2007	9.2	6.46	29.78	4.66	49.33	6.03	34.40	4.93	46.41	4.89	46.83
2008	6.85	7.35	7.28	5.35	21.89	7.00	2.16	5.97	12.80	5.93	13.45
2009	5.77	7.98	38.29	6.14	6.45	7.66	32.70	6.91	19.81	6.87	19.01
2010	7.66	8.50	11.00	7.05	7.96	8.10	5.70	7.71	0.69	7.68	0.30
2011	9.2	9.05	1.61	8.09	12.03	8.44	8.29	8.35	9.22	8.38	8.95
2012	9.73	9.73	0.01	9.29	4.52	8.84	9.15	8.82	9.32	8.89	8.59
2013	10.46	10.58	1.11	10.66	1.95	9.50	9.19	9.13	12.71	9.22	11.82
2014	11.27	11.55	2.45	12.24	8.61	10.58	6.14	9.29	17.61	9.39	16.65
2015	9	12.51	39.04	14.05	56.12	12.15	34.97	9.30	3.37	9.43	4.77
2016	10.4	13.34	28.28	16.13	55.08	14.20	36.51	9.20	11.50	9.35	10.12

“Predicted” refers to predicted values, and “APE” denotes the absolute percentage error

Table 5 summarizes forecasting accuracy for model fitting and ex-post testing. In terms of MAPE, we can see that the proposed RNGM(1,1)-P outperformed all the other models considered for model fitting and ex-post testing. In addition to MAPE, Table 5 also includes two commonly used measures mean absolute deviation (MAD) and root mean square error (RMSE) defined as (Xiao et al., 2021):

Table 5

Forecasting accuracy for model fitting and ex-post testing for Case I

	NGM(1,1)	NFLNGM(1,1)	ANGM(1,1)	NGM(1,1)-P	RNGM(1,1)-P	NN	GM(1,1)	FLNGM(1,1)	GM(1,1)-P	RGM(1,1)-P
Model fitting
MAPE	34.94	8.78	34.35	7.97	7.33	19.43	37.41	13.39		13.18	12.99
MAD	0.957	0.380	0.922	0.361	0.343	0.771	1.255	0.791		0.895	0.888
RMSE	1.125	0.402	1.059	0.496	0.491	1.186	1.692	1.212		1.467	1.473
Ex-post testing
MAPE	42.60	29.84	52.56	11.37	9.76	17.72	30.44	21.70		11.30	10.84
MAD	4.788	3.493	5.798	0.803	0.790	2.280	3.555	2.373		0.780	0.778
RMSE	5.432	4.039	6.454	1.077	0.889	2.621	4.149	2.635		0.920	0.863

“MAPE”, “MAD”, and “RMSE” refer to the mean absolute percentage error, mean absolute deviation, and root mean square error, respectively

In terms of MAD and RMSE, we can see that the RNGM(1,1)-P outperformed most of the other considered prediction models for ex-post testing, despite that the RNGM(1,1)-P was slightly inferior to RGM(1,1)-P for ex-post testing.

Case II: Magnesium alloy ingot and powder demand of Taiwan

In Taiwan, all magnesium materials are imported from abroad. Magnesium alloy ingots and scrap are the main materials imported for the magnesium alloy industry. Most of the magnesium alloys are used to manufacture the cases of 3C products, such as smart phones, tablet computers, and notebooks, by die casting and thixomolding. Taiwan is the global leader in the die casting of magnesium alloys for 3C products (Hong, 2016). To maintain this competitive advantage, magnesium alloy forecasting can be important in enabling 3C product manufacturers to recognize developing trends in related markets. Accordingly, the imported quantity of magnesium alloy ingots and scrap can be used to indirectly reflect the trends of the magnesium alloy demand of Taiwan.

We used the historical data of the import quantity of magnesium alloy ingots and scrap that was reported and analyzed by Lu (2014) and Hong (2016). Data from 2002 to 2012 were used for model fitting, and data from 2013 and 2015 were used for ex-post testing. In Fig. 3, Two circles indicate outliers outside of the data intervals that were identified after performing the robust nonlinear interval regression analysis with q = q = 0.1. The optimized NGM(1,1), NFLNGM(1,1), AGM(1,1), optimized NGM(1,1)-P, and RNGM(1,1)-P were constructed using the non-equidistant sequence by removing data in 2004 and 2010. Table 6 shows the prediction results. The prediction results obtained by NN, GM(1,1), FLNGM(1,1), GM(1,1)-P, and RGM(1,1)-P that worked on the original equidistant sequence are described in Table 7.

Fig. 3

Outlier detection of training data for Case II

Table 6

MAPEs of different grey prediction models working on non-equidistant data for Case II (unit: tons)

Year	Actual	NGM(1,1)		NFLNGM(1,1)		ANGM(1,1)		NGM(1,1)-P		RNGM(1,1)-P
		Predicted	APE	Predicted	APE	Predicted	APE	Predicted	APE	Predicted	APE
2002	5508	5508	0	5508	0	5508	0	5508	0	5508	0
2003	7257	6397.76	11.84	7263.59	0.09	6581.83	9.30	7267.98	0.15	7287.70	0.42
2004	8410
2005	6687	5883.36	12.02	5852.03	12.49	6068.40	9.25	4698.02	29.74	4648.04	30.49
2006	5347	5406.06	1.10	5329.10	0.33	5590.88	4.56	4249.25	20.53	4154.48	22.30
2007	4246	5110.91	20.37	4989.33	17.51	5294.90	24.70	4242.19	0.09	4103.29	3.36
2008	4155	4831.88	16.29	4665.00	12.27	5014.59	20.69	4332.72	4.28	4146.73	0.20
2009	3144	4568.08	45.30	4363.04	38.77	4749.12	51.05	4491.82	42.87	4254.80	35.33
2010	2503
2011	4006	4200.79	4.86	4058.29	1.31	4378.66	9.30	4836.45	20.73	4507.09	12.51
2012	6476	3859.99	40.40	3872.81	40.20	4034.10	37.71	5273.54	18.57	4837.82	25.30
2013	8387	3649.25	56.49	4004.96	52.25	3820.54	54.45	5623.37	32.95	5097.44	39.22
2014	5605	3450.02	38.45	4320.91	22.91	3618.28	35.45	6017.07	7.35	5384.00	3.94
2015	4006	3261.66	18.58	4797.33	19.75	3426.73	14.46	6455.95	61.16	5695.37	42.17

“Predicted” refers to predicted values, and “APE” denotes the absolute percentage error

Table 7

MAPEs of different prediction models working on equidistant data for Case II (unit: tons)

Year	Actual	NN		GM(1,1)		FLNGM(1,1)		GM(1,1)-P		RGM(1,1)-P
		Predicted	APE	Predicted	APE	Predicted	APE	Predicted	APE	Predicted	APE
2002	5508	5314.62	3.51	5508.00	0	5508.00	0	5508.00	0	5508	0
2003	7257	7885.03	8.65	7289.95	0.45	7193.23	0.88	7373.52	1.61	7018.43	3.29
2004	8410	8157.93	3.00	6665.46	20.74	6224.34	25.99	6414.20	23.73	6273.69	25.40
2005	6687	6573.14	1.70	6094.47	8.86	5549.44	17.01	5809.46	13.12	5601.46	16.23
2006	5347	4970.28	7.05	5572.39	4.22	4903.00	8.30	5346.57	0.01	5346.04	0.02
2007	4246	4084.35	3.81	5095.03	20.00	4283.24	0.88	4963.67	16.90	4654.20	9.61
2008	4155	3786.24	8.88	4658.57	12.12	3707.75	10.76	4633.66	11.52	4239.01	2.02
2009	3144	3836.31	22.02	4259.50	35.48	3263.66	3.81	4342.06	38.11	3943.72	25.44
2010	2503	4081.35	63.06	3894.61	55.60	3223.27	28.78	4080.13	63.01	3771.45	50.68
2011	4006	4444.59	10.95	3560.98	11.11	3924.76	2.03	3842.14	4.09	3781.52	5.60
2012	6476	4894.59	24.42	3255.93	49.72	4971.93	23.23	3624.06	44.04	4028.30	37.80
2013	8387	5420.56	35.37	2977.02	64.50	5772.27	31.18	3422.96	59.19	4519.05	46.12
2014	5605	6015.87	7.33	2721.99	51.44	6422.20	14.58	3236.61	42.25	5200.52	7.22
2015	4006	6670.18	66.50	2488.81	37.87	7143.43	78.32	3063.24	23.53	6003.90	49.87

“Predicted” refers to predicted values, and “APE” denotes the absolute percentage error

In terms of MAPE, Table 7 indicates that, among the prediction models considered here, the proposed RNGM(1,1)-P model demonstrated its superiority in the ex-post testing, despite that the RNGM(1,1)-P was slightly inferior to FLNGM(1,1) and NFLNGM(1,1) for model fitting. Indeed, one should place more emphasis on the results with respect to generalization (Luo et al, 2013). Table 8 summarizes forecasting accuracy for model fitting and ex-post testing. In addition to MAPE, we can see that the RNGM(1,1)-P outperformed most of the other considered prediction models in terms of MAD and RMSE for ex-post testing. The RNGM(1,1)-P outperformed all the other models by RMSE, but it was slightly inferior to NFLNGM(1,1) and ANGM(1,1) by MAD.

Table 8

Forecasting accuracy for model fitting and ex-post testing for Case II

	NGM(1,1)	NFLNGM(1,1)	ANGM(1,1)	NGM(1,1)-P	RNGM(1,1)-P	NN	GM(1,1)	FLNGM(1,1)	GM(1,1)-P	RGM(1,1)-P
Model fitting
MAPE	16.91	13.66	18.51	15.22	14.43	14.28	19.85	11.06	19.65	16.01
MAD	833.18	665.26	873.98	740.00	740.36	580.47	920.02	612.79	907.06	790.36
RMSE	1130.55	1042.20	1126.22	1008.02	1041.87	763.50	1172.55	639.35	1112.28	938.91
Ex-post testing
MAPE	37.84	31.64	34.78	33.82	28.45	36.40	51.27	41.36	41.66	34.40
MAD	2012.69	1829.43	1844.48	2408.22	2118.98	2546.83	2737.06	2722.79	2225.40	2353.46
RMSE	2787.44	2577.31	2666.99	2427.92	2278.45	2577.93	3327.58	2739.73	2950.85	2606.37

“MAPE”, “MAD”, and “RMSE” refer to the mean absolute percentage error, mean absolute deviation, and root mean square error, respectively

Statistical analysis

To demonstrate the superiority of the proposed RNGM(1,1)-P, statistical analysis was conducted by using the Wilcoxon signed-rank test. The forecasting accuracy of the RNGM(1,1)-P was compared with that each of the other prediction models considered for ex-post testing. Let the statistic denote the sum of ranks for the last four and three years of Cases I and II, respectively, in which the accuracy of the proposed RNGM(1,1)-P was inferior to that of one prediction model. For instance, compared to NGM(1,1), was equal to the sum of ranks for the seven years in which the accuracy of the RNGM(1,1)-P was inferior to that of the NGM(1,1). As  = 4 which was less than the critical value of 6 at a significance level of 0.1 with respect to a right-tailed test, the RNGM(1,1)-P was significantly superior to the NGM(1,1). It was found that  = 7 which was slightly greater than the critical value at a significance level of 0.1 when comparing with the NFLNGM(1,1), the RNGM(1,1)-P seems to be superior to the NFLNGM(1,1) as well. As a result, the RNGM(1,1)-P significantly outperformed most of the considered prediction models except for the GM(1,1)-P and NN. Even so, the RNGM(1,1)-P performed better than the GM(1,1)-P and NN from the viewpoint of the MAPE for ex-post testing.

Discussion

The prediction accuracy of the proposed RNGM(1,1)-P was examined using real magnesium alloy demand data from China and Taiwan. The experimental results for ex-post testing show that the MAPEs of the proposed model were minimal among all considered models. Statistical analysis further indicated that the RNGM(1,1)-P significantly outperformed the most of the considered prediction models as well. In terms of MAPEs, since the RNGM(1,1)-P outperformed the RGM(1,1)-P that can only deal with the original equidistant sequence, this indicated that the forecasting performance of the RGM(1,1)-P can be improved by removing outliers from data sequences. Altogether, the above results demonstrated the high applicability of the proposed RNGM(1,1)-P in the forecast of magnesium alloy demand.

As mentioned above, MAPE is the primary measure to examine forecasting accuracy of the considered prediction models here. From Tables 5 and 8, we can see that the RNGM(1,1)-P outperformed the other considered models in terms of MAPE for ex-post testing. However, it cannot beat the other considered models in terms of MAD and RMSE. The reason is that all considered grey prediction models were constructed by minimizing MAPE rather than MAD and RMSE. Even so, with MAD and RMSE, the proposed RNGM(1,1)-P is comparable to the other considered models. It should be noted that, for magnesium material demand forecasting, Hu (2020) applied grey relational analysis to access the relative importance of each sample for GM(1,1), and a residual modification without constructing residual models was then provided. Indeed, the contributions of this study are quite different from those of Hu (2020).

Conclusions

Many countries have established high standards for lowering CO2 emissions and energy consumption for automobiles, and constructing lightweight products has become an effective way to achieve these standards. Undoubtedly, magnesium alloys are important metal materials for the automotive industry, as well as for 3C products, aircraft, and medical appliances. As manufacturers need to invest in research and development in the use of magnesium alloys, the development of well-performed demand forecasting models is important and can help manufacturers understand the developing trends of related markets.

The main contribution of this study is to develop the grey residual modification model, namely RNGM(1,1)-P, to forecast the demand of magnesium alloy demand. To remove outliers that might give rise to a deterioration of forecasting accuracy, the robust nonlinear interval models were particularly employed to generate non-equidistant data sequences. It turns out that the proposed RNGM(1,1)-P can work on the resultant non-equidistant data sequences, and it is not necessary for the available data to comply with any statistical assumption.

The findings of this study also highlight several directions for our future work. First, this study did not consider exogenous variables such as the price of magnesium and business indicators in the context of the demand for magnesium alloy. As multivariate models may yield better predictions than single-variable models (Hu, 2020; Ma et al., 2019; Wu et al., 2019), we will investigate the applicability of several models, such as the multivariable models of Moonchai and Chutsagulprom (2020), the cumulative sum operator in grey prediction using integral transformations of Wei et al. (2020), and the background value optimization nonlinear model of Wu et al. (2020), to forecast the demand for magnesium alloy. Next, since there are shortcomings of MAPE outlined by Vandeput (2021), we will investigate how to incorporate the other measures such as MAD and RMSE into the construction of the proposed models to improve forecasting accuracy.