A Study on Agricultural Commodity Price Prediction Model Based on Secondary Decomposition and Long Short-Term Memory Network

Abstract

In order to address the significant prediction errors resulting from the substantial fluctuations in agricultural product prices and the non-linear features, this paper proposes a hybrid forecasting model based on variational mode decomposition (VMD), ensemble empirical mode decomposition (EEMD), and long short-term memory networks (LSTM). This combined model is referred to as the VMD–EEMD–LSTM model. Initially, the original time series of agricultural product prices undergoes decomposition using VMD to obtain a series of variational mode functions (VMFs) and a residual component with higher complexity. Subsequently, the residual component undergoes a secondary decomposition using EEMD. All components are then fed into an LSTM model for training to obtain predictions for each component. Finally, the predictions for each component are linearly combined to generate the ultimate price forecast. To validate the effectiveness of the VMD–EEMD–LSTM model, empirical analyses were conducted for one-step and multi-step forecasts using weekly price data for pork, Chinese chives, shiitake mushrooms, and cauliflower from China’s wholesale agricultural markets. The results indicate that the composite model developed in this study provides enhanced forecasting accuracy.

1. Introduction

Agricultural products, as essential commodities in the daily consumption of Chinese residents, directly influence people’s economic well-being. The high or low prices of agricultural products are closely linked to individuals’ economic circumstances. Because of the characteristics of significant price volatility, nonlinearity, and non-stationarity, frequent price fluctuations add complexity to the challenges faced by producers, consumers, and decision-makers. Therefore, constructing accurate forecasting models for agricultural products is of significant importance for government macroeconomic regulation, balancing agricultural market supply and demand, and improving the living standards of residents.

Extensive research has been conducted both domestically and internationally on predicting agricultural product prices. Currently, the primary forecasting methods encompass traditional econometric models, intelligent models, and hybrid models.

Traditional econometric models are built on a solid theoretical foundation and offer strong interpretability. Examples of such models include autoregressive moving average (ARMA), autoregressive integrated moving average (ARIMA), grey model (GM), and generalized autoregressive conditional heteroskedasticity (GARCH) models [1,2,3,4,5,6]. Their advantage lies in having a clear internal structure, with each parameter having a well-defined interpretation. However, they have limitations when it comes to addressing nonlinear issues and are unable to capture the non-stationary, nonlinear, and multiscale characteristics of time series data [7,8].

In recent years, artificial intelligence model technology has rapidly advanced. Given that these models do not need to adhere to statistical assumptions, they are well-suited to address non-stationary and nonlinear problems. Typical models include the backpropagation neural network (BP), recurrent neural networks (RNN), artificial neural networks (ANN), and long short-term memory networks (LSTM), among other machine learning models [9,10,11,12]. Compared to traditional econometric models, intelligent models often achieve higher accuracy, but it can be challenging to find a single model that performs well in all scenarios [13]. At the same time, mainstream intelligent modeling methods struggle to describe the interdependencies between time series data [14]. They also require a large amount of data for training, and insufficient training data can lead to issues like poor fit and low accuracy in predictions [15]. Consequently, researchers have been focusing on improving and optimizing artificial intelligence optimization algorithms [16,17,18,19]. While these optimized intelligent models offer better predictive accuracy, the process can be intricate and the runtime can be lengthy.

The third category of hybrid modeling methods combines the advantages of both traditional econometric models and intelligent models. These models incorporate prior assumptions and data processing for prediction problems, thereby enhancing the model’s fitting capability. They have found widespread application in the prediction of agricultural product prices [20,21]. For instance, Ling et al. [22] developed a GM-VAR combination model to forecast the prices of various livestock products, and their predictions outperformed single forecasting models. Wu et al. [23] combined variational mode decomposition (VMD) with bald eagle optimization in an LSTM (long short-term memory) framework, creating the VMD-IBES-LSTM model for predicting seafood prices. However, this model overlooked essential information within the residuals obtained after VMD decomposition.

To further enhance predictive accuracy, a complex system theory approach that combines decomposition techniques with forecasting models has been proposed by Wang Shouyang [24]. This approach can effectively overcome the limitations of traditional econometric models and artificial intelligence models [25]. In fact, in most combination forecasting models, signal processing methods are used to decompose time series, and intelligent models are employed to predict the components after decomposition [26]. Typical sequence decomposition methods include empirical mode decomposition (EMD), which has a significant advantage in handling non-stationary and nonlinear complex signals as it does not rely on arbitrary basis functions. For example, Pankaj and others [27] combined the EMD method with support vector regression (SVR) to create an EMD-SVR model for predicting the wholesale price index (WPI) of chili peppers, and they found that the results were better than a standard SVR model. However, the EMD method is prone to mode mixing issues. Fang et al. [28] combined the ensemble empirical mode decomposition (EEMD) method with SVM, NN, and ARIMA models to forecast various agricultural products, and they found that the predictions were better than those of individual models. Nevertheless, the EEMD method struggles to decompose components with similar frequencies and faces issues related to envelopes and undershoots, which limits its decomposition results [29]. As an improved decomposition technique, variational mode decomposition (VMD) can adaptively decompose the effective components corresponding to each center frequency, leading to higher decomposition accuracy. VMD is more effective in feature selection for prediction models [30].

Given the limitations of traditional econometric models in addressing nonlinear problems and the inability of mainstream intelligent models to capture the interdependencies among time series data, we opt for the LSTM model for forecasting. This choice is motivated by the LSTM model’s effectiveness in capturing long-term trends and dependencies within time series data while adapting to different time scales. Considering the issues of mode mixing in EMD decomposition and the inability of EEMD to separate components with similar frequencies, we select the VMD method, which offers higher decomposition precision, to perform a single decomposition on the original price series.

From previous studies, it is evident that the predictive performance of composite models surpasses that of individual intelligent models. However, in some composite models utilizing decomposition techniques, the prediction overlooks crucial information within the residual terms post-decomposition, leading to an insufficient extraction of features from the sequence and preventing the model accuracy from reaching its optimum. Therefore, this paper’s contribution lies in the fact that no one has yet applied second-order decomposition to agricultural commodity price forecasting. This study proposes a VMD–EEMD bidecomposition method, using the components obtained after two decompositions as input for LSTM, to establish a VMD–EEMD–LSTM model for predicting agricultural commodity prices.

In addition, the VMD–EEMD–LSTM model proposed in this study represents a hybrid approach that combines frequency domain decomposition with the LSTM model. Unlike traditional hybrid models that involve restructuring the econometric model into a machine learning framework by integrating different types of models, the VMD–EEMD–LSTM model takes a unique approach. Through two decomposition steps, it reduces data complexity and introduces the LSTM model to effectively capture the multi-level structure of time series data. This integration enhances the model’s ability to learn and predict time series, providing greater flexibility and adaptability compared to conventional hybrid models.

The main work of this paper includes:

(i)

In response to the non-stationary and high complexity characteristics of the agricultural product price sequences used in the experiments in this paper, VMD is employed to decompose the original complex sequences into multiple sub-sequences based on their central frequencies. Additionally, a secondary decomposition is performed on the residual components containing complex information.

(ii)

After decomposing the original sequence with VMD, different variational mode functions (VMFs) and residual components (Res) are obtained. The residual components are further decomposed using EEMD, and then all these sequences are input into the LSTM model for prediction.

(iii)

To validate the predictive performance and accuracy of the proposed model, various metrics such as mean absolute error (MAE), mean absolute percentage error (MAPE), and root mean square error (RMSE) are chosen for evaluation. The baseline models and the method proposed in this paper are used to perform one-step and multi-step predictions for weekly agricultural product prices.

2. Materials and Methods

This section mainly describes the construction of the VMD–EEMD–LSTM hybrid model, including the VMD algorithm, EEMD algorithm, LSTM model, and the overall modeling steps.

2.1. Variational Mode Decomposition

Variational mode decomposition (VMD) is an adaptive and completely non-recursive modal variational and signal processing method proposed by Dragomiretskiy and Zosso [31]. Unlike empirical mode decomposition, VMD can determine the number of modal decompositions for a given signal sequence based on the actual situation. It decomposes the original input signal into $k$ variational modal functions (VMFs) that follow central frequencies and have finite bandwidth. By minimizing the sum of bandwidth estimates for all modes, it obtains the corresponding modal component signals $u_k$, where variations in each modal signal are centered around its central frequency ${\omega }_k$. The VMD variational constraint model is as follows:

$$\begin{array}{c}\\ \{\begin{array}{l}\underset{\left\{u_k\right\},\left\{{\omega }_k\right\}}{\min }\left\{{\displaystyle \sum _k\Vert {\mathsf{\partial }}_t\left[\left(\delta \left(t\right)+\frac{j}{\pi t}\right)\ast u_k\left(t\right)\right]e^{-j{\omega }_{k}t}\Vert {}_2^2}\right\},\\ s.t. {\displaystyle \sum _k{u}_k=f(t).}\end{array}\end{array}$$

(1)

In the equation, $\left\{u_k\right\}:=\left\{u_1,u_2,\cdots ,u_k\right\}$, represents the variational mode function (VMF) obtained after decomposition. $\left\{{\omega }_k\right\}:=\left\{{\omega }_1,{\omega }_2,\cdots ,{\omega }_k\right\}$ represents the central frequency corresponding to each variational mode function (VMF). $\delta \left(t\right)$ is the Dirac function, $k$ represents the number of modes, $\ast$ is the convolution operator, $f\left(t\right)$ is the original signal, ${\mathsf{\partial }}_t$ is the gradient operation. $\left(\delta \left(t\right)+\frac{j}{\pi t}\right)\ast u_k\left(t\right)$ is the spectrum of $u_k\left(t\right)$ obtained through the Hilbert transformation. The exponential term $e^{-j{\omega }_{k}t}$ is used to adjust the estimates of each ${\omega }_k$, then integrate the spectra of the modes into the baseband. If you want to solve the constrained optimization problem described above, it needs to be transformed into an unconstrained problem by introducing quadratic penalty factor $\alpha$ and Lagrange multiplier operator $\lambda \left(t\right)$, resulting in the expression of the unconstrained variational problem:

$$L\left(\left\{u_k\right\},\left\{{\omega }_k\right\},\lambda \right)=\alpha {\displaystyle \sum _k\Vert {\mathsf{\partial }}_t}\left[\left(\delta \left(t\right)+\frac{j}{\pi t}\right)\ast u_k\left(t\right)\right]e^{-j{\omega }_{k}t}\Vert {}_2^2+\Vert f\left(t\right)-{\displaystyle \sum _k{u}_k\left(t\right)}\Vert {}_2^2+\langle \lambda \left(t\right),f\left(t\right)-{\displaystyle \sum _k{u}_k\left(t\right)}\rangle .$$

(2)

By finding the extremum of this Lagrangian function, we obtain the expressions for VMF $u_k$ and the center frequency ${\omega }_k$ in the unconstrained variational problem:

$${\widehat{u}}_k^{n+1}\left(\omega \right)=\frac{\widehat{f}\left(\omega \right)-{\displaystyle \sum {}_{i\ne k}{\widehat{u}}_i\left(\omega \right)+\frac{\widehat{\lambda }\left(\omega \right)}{2}}}{1+2\alpha {\left(\omega -{\omega }_k\right)}^2},$$

(3)

$${\widehat{\omega }}_k^{n+1}=\frac{{\displaystyle {\int }_0^{\infty }\omega |{\widehat{u}}_k{}^{n+1}\left(\omega \right)|{}^{2}d\omega }}{{\displaystyle {\int }_0^{\infty }{\widehat{u}}_k{}^{n+1}|\left(\omega \right)|{}^{2}d\omega }}.$$

(4)

The steps of the VMD decomposition method are as follows:

(1)

Initialize $\left\{u_k^1\right\}$, $\left\{{\omega }_k^1\right\}$, $\left\{{\lambda }_k^1\right\}$, and $n$, and set an appropriate number of VMF components $k$;

(2)

Update $u_k$ and ${\omega }_k$ based on Equations (3) and (4);

(3)

Update the value of $\lambda$:

$${\widehat{\lambda }}^{n+1}\left(\omega \right)={\widehat{\lambda }}^n\left(\omega \right)+\tau \left[\widehat{f}\left(\omega \right)-{\displaystyle \sum _k{\widehat{u}}_k^{n+1}\left(\omega \right)}\right]$$

(5)

(4)

Given the precision criterion $\epsilon >0$, if satisfied:

$$\frac{{\displaystyle \sum {}_k}\Vert {\widehat{u}}_k^{n+1}-{\widehat{u}}_k^n\Vert {}_2^2}{\Vert {\widehat{u}}_k^n\Vert {}_2^2}<\epsilon$$

(6)

Then, stop the iteration; otherwise, return to step (2).

In the equation, ${\widehat{u}}_k^n\left(\omega \right)$, $\widehat{f}\left(\omega \right)$, and ${\widehat{\lambda }}^n\left(\omega \right)$ correspond to the Fourier transforms of ${\widehat{u}}_k^n$, $f\left(t\right)$, and ${\lambda }^n$ respectively.

2.2. Ensemble Empirical Mode Decomposition

The empirical mode decomposition (EMD) method can handle non-stationary and nonlinear signals, but it is prone to modal aliasing, where the resulting intrinsic mode functions (IMFs) exhibit mode mixing. To mitigate this issue, Wu and Huang [32] introduced the ensemble empirical mode decomposition (EEMD) technique, which involves injecting a low-amplitude white noise sequence into the signal to assist in the analysis. The decomposition steps are as follows:

(1)

Add a signal $n_i\left(t\right)$ of the same length to the initial signal $x\left(t\right)$, creating a new signal:

$$x_i\left(t\right)=x\left(t\right)+n_i\left(t\right).$$

(7)

(2)

Perform EMD decomposition on the $x\left(t\right)$ obtained in step (1) to obtain various IMF mode components:

$$x_i\left(t\right)={\displaystyle \sum _{j=1}^J{c}_{i,j}\left(t\right)}+r_i\left(t\right).$$

(8)

(3)

Repeat steps (1) and (2) $N$ times, each time adding white noise sequences with different amplitudes during the decomposition to obtain various IMF components:

$$c_{1,j}\left(t\right),c_{2,j}\left(t\right)\dots c_{N,j}\left(t\right),J=1,2,\cdots J.$$

(9)

(4)

Since the statistical average of uncorrelated sequences is zero, the corresponding IMF components are averaged to obtain the final decomposition result of EEMD:

$$c_j\left(t\right)=\frac{1}{N}{\displaystyle \sum _{i=1}^N{c}_{i,j}}\left(t\right).$$

(10)

In the equation, $x_i\left(t\right)$ represents the sequence after adding white noise for the $i$-th time; $c_{i,j}\left(t\right)$ is the $j$-th IMF component obtained after the $i$-th addition of white noise, $r_i\left(t\right)$ is the residual function, $J$ is the number of IMFs.

2.3. Long Short-Term Memory

Long short-term memory (LSTM) is an improved type of recurrent neural network that introduces the concept of a “cell state” and “gates” to effectively address the vanishing or exploding gradient problem in traditional RNNs. LSTM is designed as a sequential structure primarily consisting of a Forget Gate, an Input Gate, and an Output Gate. The LSTM algorithm uses these gate structures to maintain and update the memory cell state.

The information forgetting of the memory cell is determined by the Forget Gate, and the specific formula is as follows:

$$f_t=\sigma \left(W_f\left[{\stackrel{⌢}{y}}_{t-1},x_t\right]+b_f\right).$$

(11)

The Sigmoid layer’s output of the Input Gate serves as the updated information and is determined by the following formula:

$$\begin{gathered} i_t=\sigma \left(W_i\left[{\stackrel{⌢}{y}}_{t-1},x_t\right]+b_i\right), \\ \tilde{C}=\tanh \left(W_c\left[{\stackrel{⌢}{y}}_{t-1},x_t\right]+b_c\right), \\ {C}_t=f_t{C}_{t-1}+i_t{\tilde{C}}_t. \end{gathered}$$

(12)

The Output Gate determines the value of the next hidden state, and the specific formula is as follows:

$$\begin{gathered} o_t=\sigma \left(W_o\left[{\stackrel{⌢}{y}}_{t-1},x_t\right]+b_o\right), \\ {\stackrel{⌢}{y}}_t=o_t\tanh \left(C_t\right). \end{gathered}$$

(13)

In the equations provided, $\sigma$ represents the Sigmoid function; $\tanh$ represents the hyperbolic tangent (tanh) function; $W_i$, $W_c$, $W_f$, and $W_o$ represent weight matrices; $b_i$, $b_c$, $b_f$, and $b_o$ represent bias matrices; ${\stackrel{⌢}{y}}_{t-1}$ and ${\stackrel{⌢}{y}}_t$ are, respectively the inputs from the previous unit and the outputs of this unit; $C_{t-1}$ and $C_t$ are information states of two adjacent units.

2.4. VMD–EEMD–LSTM Prediction Model

The VMD decomposition can break down complex signals into several lower complexity modal components, with each mode representing a frequency component within the data. Through VMD decomposition, the prediction model can better capture the local and frequency-domain features of the data. EEMD is an enhanced empirical mode decomposition method that better handles nonlinear and non-stationary signals through multiple iterations and the addition of random noise. EEMD contributes to capturing the nonlinear dynamic characteristics in the data, enhancing the model’s adaptability to complex patterns. Therefore, using EEMD to decompose the residual term RES after VMD decomposition helps reduce the complexity of the residual components. The LSTM model is employed to process and learn from time series data, effectively capturing the long-term trends and dependencies in the time series while adapting to different time scales. The dual decomposition with VMD and EEMD effectively reduces the complexity of the data, enabling the LSTM model to better learn the features of the time series and thereby improving the predictive accuracy of the VMD–EEMD–LSTM model. The specific modeling steps are as follows:

(i)

By applying the VMD decomposition algorithm to the original agricultural product price time series, lower-complexity VMF components and a residual term Res (obtained by subtracting each VMF component from the original sequence) are obtained.

(ii)

Each VMF component obtained in step (1) is individually fed into an LSTM model for training, resulting in predictions for each VMF component.

(iii)

Using EEMD for the secondary decomposition of the Res obtained in step (1), each IMF component derived from this decomposition is individually input into an LSTM for prediction. The predicted values for each IMF component are then summed to obtain the prediction result for the residual term Res.

(iv)

Combine the prediction results for each VMF component and the residual term Res to obtain the final prediction result for the original sequence.

The modeling process, as shown in Figure 1.

2.5. Evaluation Metrics

In order to evaluate the predictive performance of the VMD–EEMD–LSTM model for agricultural product prices, this paper employs three evaluation metrics: root mean square error (RMSE), mean absolute error (MAE), and mean absolute percentage error (MAPE). These three evaluation metrics are used to assess the errors between the model’s predictions and the actual values, and a smaller value for these metrics indicates higher predictive accuracy. The specific formulas for the three evaluation metrics are as follows:

$$\begin{gathered} RMSE=\sqrt{\frac{1}{n}{\displaystyle \sum _{i=1}^n{\left(y_i-{\stackrel{⌢}{y}}_i\right)}^2}} \\ MAE=\frac{1}{n}{\displaystyle \sum _{i=1}^n|y_i-{\stackrel{⌢}{y}}_i|} \\ MAPE=\frac{1}{n}{\displaystyle \sum _{i=1}^n|\frac{y_i-{\stackrel{⌢}{y}}_i}{y_i}|} \end{gathered}$$

(14)

In the formulas provided, $y_i$ and ${\stackrel{⌢}{y}}_i$ represent the actual (true) value and predicted value of the $i$-th sample, respectively. $n$ represents the total number of test samples.

3. Results

In this section, we conducted the following experiments. Firstly, we selected the optimal $k$ values for VMD decomposition and performed secondary EEMD decomposition on the residuals for pork, leeks, shiitake mushrooms, and cauliflower. Secondly, we determined the appropriate time step for prediction in the VMD–EEMD–LSTM model. Finally, we obtained the ultimate prediction results and evaluated the model performance using metrics such as MAE, MAPE, and RMSE, conducting comparative analysis with other benchmark models.

3.1. Data Source

The empirical analysis in this paper is based on the weekly average pork price. The data is sourced from the National Agricultural Products Wholesale Market of the Chinese Ministry of Agriculture and Rural Affairs. The sample data covers a period from 1 March 2018, to 28 February 2023, with a total of 261 data points. The time step is set at five (Why the time step is chosen as five is demonstrated in Section 3.3.), which means historical data from the previous five weeks is used to predict the price for the sixth week. The first 183 data points are designated as the training sample, and the remaining 73 data points are used as the test sample, resulting in a training-to-test sample ratio of 7:3. The original sequence of pork prices is depicted in Figure 2, and the standard statistical results for the sample data are presented in

Table 1. Conventional statistical analysis of pork, Chinese chives, shiitake mushrooms, and cauliflower price sequences.

Agricultural Product Categories	Average Value	Maximum Value	Minimum Value	Standard Deviation	Kurtosis	Skewness
Pork	29.34	52.38	15.83	11.14	−1.21	0.58
Chinese chives	4.77	7.78	3.0	1.22	−0.66	0.43
Shiitake mushrooms	7.33	11.14	4.14	1.06	1.85	0.30
Cauliflower	4.36	8.01	2.78	0.99	2.86	1.45

. This empirical study was conducted using Matlab R2021a.

3.2. The Process of Secondary Decomposition

3.2.1. VMD Primary Decomposition

Before performing VMD decomposition on the raw pork price series, it is necessary to determine the appropriate number of mode components, denoted as $k$, as failure to do so can affect the effectiveness of the VMD decomposition and subsequently impact the overall predictive performance of the model. If the value of $k$ is too large, it may lead to the overlapping of central frequencies of components, resulting in excessive decomposition. If $k$ is too small, the complexity of the original sequence cannot be sufficiently reduced. Therefore, in this paper, according to the central frequency method [33], we determine the optimal value of k by initially selecting a range of modal numbers $k$ from 3 to 8 for pre-decomposition. We obtain the central frequency of each modal component. When the central frequencies of adjacent modes are close, it is considered as an indication of over-decomposition of the signal, as shown in

Table 2. The central frequencies for different values of k.

K	VMF1	VMF2	VMF3	VMF4	VMF5	VMF6	VMF7	VMF8
3	0.47	98.01	643.83
4	0.36	81.87	348.35	1309.36
5	0.29	70.56	284.28	731.68	3415.09
6	0.29	68.66	268.29	605.97	1352.75	4109.05
7	0.28	66.47	259.82	564.02	1090.77	1936.27	4187.84
8	0.27	64.74	251.67	537.42	994.39	1610.41	2683.42	4288.26

. A comparison of the central frequencies of the last layer of mode components is made, and when these frequencies become relatively stable, the optimal $k$ value can be identified. Setting VMD parameters as follows: $\alpha =2000$, $\epsilon ={10}^{-7}$.

From Table 2, it can be observed that when $k>6$, the central frequencies of the last layer of mode components tend to stabilize. Therefore, the optimal $k$ value is chosen to be 7. Subsequently, the raw pork price series is decomposed using both VMD and EEMD, as illustrated in Figure 3. It is evident that the VMD decomposition results in a smoother sequence, effectively reducing the complexity of the original series.

3.2.2. EEMD Secondary Decomposition

From Figure 3, it can be observed that the residual component “Res” after VMD decomposition still exhibits significant fluctuations, indicating a high level of data complexity, with a calculated sample entropy of 0.8914. To improve the overall predictive performance of the model, EEMD is applied to decompose the residual component “Res” The white noise standard deviation is set to 0.4, and the number of ensemble trials is configured as 100. The residual component “Res” after VMD decomposition is decomposed into 7 IMF (intrinsic mode function) components and 1 residual trend component, as shown in Figure 4.

3.3. Experiments with Different Time Steps

Considering the distinct volatility patterns of various agricultural products, setting an appropriate time step for prediction can help reduce forecasting errors, i.e., using price data from previous weeks to predict prices for the following week. In order to validate the optimal time step for price prediction using the LSTM model for pork, Chinese leeks, shiitake mushrooms, and cauliflower, experiments with different time steps were conducted, keeping the parameters consistent with Section 3.4.1. The experimental results are shown in

Table 3. Experiments with Different Time Steps for Various Agricultural Products.

Agricultural Product Categories	Evaluation Metrics	Time Step
		3	5	7	9	11
Pork	MAE	0.9151	0.8547	0.9417	0.9147	0.9485
	MAPE	0.0371	0.0354	0.0368	0.0356	0.0367
	RMSE	1.2450	1.2268	1.2648	1.2409	1.2590
Chinese chives	MAE	0.3339	0.2899	0.3150	0.3056	0.3116
	MAPE	0.0679	0.0590	0.0633	0.0617	0.0629
	RMSE	0.3912	0.3549	0.3878	0.3726	0.3822
Shiitake mushrooms	MAE	0.4850	0.4310	0.4517	0.4601	0.4716
	MAPE	0.0609	0.0541	0.0552	0.0566	0.0559
	RMSE	0.6038	0.5518	0.5746	0.5749	0.6073
Cauliflower	MAE	0.3514	0.3239	0.3424	0.3365	0.3309
	MAPE	0.0831	0.0755	0.0816	0.0799	0.0783
	RMSE	0.4303	0.4135	0.4306	0.4390	0.4437

. The results revealed that for pork, Chinese leeks, shiitake mushrooms, and cauliflower, the minimum prediction errors were observed at a time step of 5. Therefore, in practical forecasting, the optimal time step for each product was chosen for prediction.

3.4. Price Prediction Based on VMD–EEMD–LSTM

3.4.1. LSTM Hyperparameter Configuration

Training epochs, learning rate, and other factors are important parameters affecting LSTM model accuracy. We experimented with different settings for the LSTM network, including training epochs of 100, 200, and 400, learning rates of 0.01, 0.001, and 0.005, and batch sizes of 32, 64, and 128. After multiple comparisons of validation set errors, we found that with 200 training epochs, a learning rate of 0.005, a batch size of 32, 70 LSTM hidden units, and a learning rate decay factor of 0.2, the model achieved the lowest loss.

3.4.2. Prediction Results Comparison and Analysis

In order to better demonstrate the predictive performance of the VMD–EEMD–LSTM model, this study, while maintaining consistency in the experimental process, concurrently employed RF (random forest), ELM (extreme learning machine), LSTM, EEMD–LSTM, and VMD–LSTM models for predicting pork prices. Different prediction horizons (1, 2, 3) were set for comparison. The error results for different prediction horizons are presented in

Table 4. Different prediction horizon error results.

Prediction Horizon	Evaluation Metrics	RF	ELM	LSTM	EEMD-LSTM	VMD-LSTM	VMD–EEMD–LSTM
One	MAE	1.2864	1.1915	0.8924	0.7977	0.5675	0.4790
	MAPE	0.0511	0.0473	0.0354	0.0344	0.0231	0.0195
	RMSE	1.6566	1.6145	1.2268	0.9611	0.7071	0.5926
Two	MAE	1.9177	1.6274	1.5269	1.1125	0.6512	0.5617
	MAPE	0.0757	0.0639	0.0603	0.0485	0.0264	0.0229
	RMSE	2.3496	2.2520	2.0292	1.3649	0.7950	0.6775
Three	MAE	2.8440	2.5612	2.1314	1.6517	0.7820	0.7075
	MAPE	0.1105	0.1023	0.0850	0.0730	0.0313	0.0286
	RMSE	3.4253	3.222	2.7491	2.0175	0.9408	0.8529

, and the predictions of each model at various prediction horizons are illustrated in Figure 5, Figure 6 and Figure 7.

Based on Table 4 and Figure 5, Figure 6 and Figure 7, it can be observed that within a single model, LSTM outperforms RF and ELM in terms of prediction performance across different prediction horizons. It achieves lower MAE, MAPE, and RMSE values. Taking a prediction horizon of 1 as an example, when compared to RF, LSTM shows a reduction of 30.63% in MAE, 30.72% in MAPE, and 25.94% in RMSE. When compared to ELM, LSTM results in a reduction of 25.10% in MAE, 25.16% in MAPE, and 24.01% in RMSE. This result indicates that within a single model, the LSTM model exhibits superior predictive performance.

Although the LSTM model exhibits better predictive performance than other individual models, when compared to models that incorporate decomposition techniques, there is still a significant gap in the evaluation metric values. This is because decomposition methods can reduce the complexity of the time series, thereby enhancing the predictive accuracy of the models. When the prediction horizon is 1, the proposed VMD–EEMD–LSTM model, compared to RF, shows a reduction of 62.76%, 61.84%, and 64.23% in MAE, MAPE, and RMSE, respectively. In comparison to ELM, the reductions are 59.80%, 58.77%, and 63.30%, while compared to LSTM, they are 46.32%, 44.92%, and 51.70%. When compared to EEMD–LSTM, the reductions are 39.95%, 43.31%, and 38.34%, and in comparison to VMD–LSTM, they are 15.59%, 15.58%, and 16.19%.

For a prediction horizon of 2, the VMD–EEMD–LSTM model, in comparison to RF, exhibits reductions of 70.71%, 69.75%, and 71.17% in MAE, MAPE, and RMSE, respectively. When compared to ELM, the reductions are 65.48%, 64.16%, and 69.92%, and in comparison to LSTM, they are 63.21%, 62.02%, and 66.61%. Compared to EEMD–LSTM, the reductions are 49.51%, 52.78%, and 50.36%, and compared to VMD–LSTM, they are 13.74%, 13.26%, and 14.78%.

At a prediction horizon of 3, the VMD–EEMD–LSTM model, compared to RF, shows reductions of 75.12%, 74.12%, and 75.10% in MAE, MAPE, and RMSE, respectively. When compared to ELM, the reductions are 72.38%, 72.04%, and 73.53%, and in comparison to LSTM, they are 66.81%, 66.35%, and 68.98%. When compared to EEMD–LSTM, the reductions are 57.17%, 60.82%, and 57.72%, and compared to VMD–LSTM, they are 9.53%, 8.63%, and 9.34%.

This indicates that the proposed VMD–EEMD–LSTM model exhibits excellent predictive performance, both in single-step and multi-step forecasting. The trend predictions align well with the original data, demonstrating a good fit.

3.5. Generalizability Experiment of the VMD–EEMD–LSTM Model

To validate the predictive applicability of the model proposed in this study on price sequences of other agricultural products, price data for Chinese chives, shiitake mushrooms, and cauliflower were selected for prediction validation. The sample time range for all datasets is from 26 January 2020, to 25 December 2022, and the data sources are the same as described in Section 3.1, with model parameters consistent with Section 3.4.1. The conventional statistical results for each dataset are displayed in Table 1 and the original price sequences for chinese leeks, shiitake mushrooms, and cauliflower are shown in Figure 8. The optimal $k$ values for Chinese chives, shiitake mushrooms, and cauliflower are selected as discussed in Section 3.2.1, specifically shown in

Table 5. The central frequencies of the original price sequences of Chinese chives, shiitake mushrooms, and cauliflower under different values of $k$ in the decomposition. (a) The central frequencies of Chinese chives price sequences under different values of $k$ in the decomposition. (b) The central frequencies of shiitake mushroom price sequences under different values of $k$ in the decomposition. (c) The central frequencies of cauliflower price sequences under different values of $k$ in the decomposition.

(a)
K	VMF1	VMF2	VMF3	VMF4	VMF5	VMF6	VMF7	VMF8
3	0.79	349.9	3849.5
4	0.70	314.3	1261.4	4096.3
5	0.60	283.1	1160.4	3426.5	4194.0
6	0.26	220.9	566.1	1354.7	3531.5	4274.2
(b)
K	VMF1	VMF2	VMF3	VMF4	VMF5	VMF6	VMF7	VMF8
5	0.10	333.0	1780.5	2529.8	3654.7
6	0.12	339.4	1748.0	2243.8	2904.6	4405.4
7	0.11	313.0	1152.7	1875.4	2673.5	3612.6	4451.9
8	0.10	302.5	1112.1	1790.0	2272.5	2858.8	3655.7	4464.5
(c)
K	VMF1	VMF2	VMF3	VMF4	VMF5	VMF6	VMF7	VMF8
5	0.46	374.7	859.2	1738.4	2976.0
6	0.43	360.2	802.3	1615.5	2625.0	4146.6
7	0.40	342.5	736.0	1328.2	1991.8	2873.3	4248.0
8	0.39	337.2	715.5	1261.0	1895.3	2599.8	3350.8	4383.9

. According to Table 5, the optimal values of $k$ for VMD decomposition are selected as four for leek price sequence, six for shiitake mushroom price sequence, and six for cauliflower price sequence. The VMD decomposition sequences for Chinese chives, shiitake mushrooms, and cauliflower are shown in Figure 9a,c,e, respectively. The secondary decomposition of their residuals (Res) using EEMD is depicted in Figure 9b,d,f.The error results for different prediction horizons are listed in

Table 6. Error comparison for different prediction horizons in Chinese chives.

Prediction Horizon	Evaluation Metrics	RF	ELM	LSTM	EEMD-LSTM	VMD-LSTM	VMD–EEMD–LSTM
One	MAE	0.3500	0.3155	0.2899	0.1978	0.1570	0.1288
	MAPE	0.0709	0.0645	0.0590	0.0389	0.0335	0.0269
	RMSE	0.4228	0.3845	0.3549	0.2562	0.1923	0.1588
Two	MAE	0.5384	0.4909	0.4752	0.2551	0.2384	0.1990
	MAPE	0.1105	0.1030	0.0973	0.0526	0.0507	0.0420
	RMSE	0.6652	0.6018	0.5662	0.3052	0.2953	0.2464
Three	MAE	0.6877	0.6387	0.6022	0.3434	0.2853	0.2550
	MAPE	0.1428	0.1359	0.1246	0.0729	0.0614	0.0551
	RMSE	0.8435	0.7630	0.7350	0.4047	0.3649	0.3149

,

Table 7. Error comparison for different prediction horizons in shiitake mushrooms.

Prediction Horizon	Evaluation Metrics	RF	ELM	LSTM	EEMD-LSTM	VMD-LSTM	VMD–EEMD–LSTM
One	MAE	0.6912	0.6304	0.4310	0.3781	0.1785	0.1443
	MAPE	0.0828	0.0792	0.0541	0.0534	0.0221	0.0185
	RMSE	0.9300	0.7800	0.5518	0.5246	0.2157	0.1718
Two	MAE	0.8739	0.7604	0.7269	0.4993	0.2167	0.1805
	MAPE	0.1062	0.0952	0.0901	0.0700	0.0267	0.0226
	RMSE	1.1139	0.9661	0.9378	0.7330	0.2733	0.2153
Three	MAE	1.0250	0.9748	0.8739	0.6192	0.2457	0.2255
	MAPE	0.1253	0.1192	0.1068	0.0884	0.0307	0.0285
	RMSE	1.2406	1.2242	1.1230	0.9254	0.3170	0.2715

and

Table 8. Error comparison for different prediction horizons in cauliflower.

Prediction Horizon	Evaluation Metrics	RF	ELM	LSTM	EEMD-LSTM	VMD-LSTM	VMD–EEMD–LSTM
One	MAE	0.3904	0.3568	0.3239	0.2597	0.1596	0.1528
	MAPE	0.0978	0.0826	0.0755	0.0669	0.0409	0.0389
	RMSE	0.4893	0.4714	0.4136	0.2881	0.1960	0.1889
Two	MAE	0.6209	0.5925	0.5834	0.3931	0.2317	0.2248
	MAPE	0.1590	0.1453	0.1444	0.1014	0.0596	0.0575
	RMSE	0.7326	0.7213	0.7033	0.4365	0.2814	0.2744
Three	MAE	0.7387	0.7025	0.7012	0.5324	0.3067	0.3013
	MAPE	0.1925	0.1790	0.1772	0.1362	0.0804	0.0788
	RMSE	0.8589	0.8366	0.8244	0.5917	0.3785	0.3717

. The prediction results for various models (using a prediction horizon of 1 as an example) can be found in Figure 10, Figure 11 and Figure 12.

Combining Table 6, Table 7 and Table 8, in the price prediction of Chinese chives, oyster mushrooms, and cauliflower, when the forecasting step is 1, the model proposed in this study, VMD–EEMD–LSTM, compared to RF, demonstrates improvements in MAE by 63.2%, 79.12%, and 60.86%, and in RMSE by 61.57%, 81.53%, and 61.39%, respectively. Compared to ELM, the improvements in MAE are 59.18%, 77.11%, and 57.17%, and in RMSE are 58.70%, 77.97%, and 59.93% for Chinese chives, oyster mushrooms, and cauliflower, respectively. In comparison to LSTM, the improvements in MAE are 55.01%, 66.52%, and 52.91%, and in RMSE are 54.73%, 68.87%, and 54.63% for Chinese chives, oyster mushrooms, and cauliflower, respectively. Compared to EEMD–LSTM, the improvements in MAE are 34.88%, 61.84%, and 41.16%, and in RMSE are 38.02%, 67.25%, and 34.43% for Chinese chives, oyster mushrooms, and cauliflower, respectively. Lastly, compared to VMD–LSTM, the improvements in MAE are 17.96%, 19.16%, and 4.26%, and in RMSE are 17.42%, 20.35%, and 3.62% for Chinese chives, oyster mushrooms, and cauliflower, respectively.

When the forecasting step is 2, the model proposed in this study, VMD–EEMD–LSTM, compared to RF, exhibits significant improvements in MAE and RMSE for Chinese chives, oyster mushrooms, and cauliflower. The improvements in MAE are 63.04%, 79.35%, and 63.79%, and in RMSE are 62.96%, 80.67%, and 62.54%, respectively. Compared to ELM, the improvements in MAE are 59.46%, 76.26%, and 63.06%, and in RMSE are 59.06%, 77.71%, and 61.96% for Chinese chives, oyster mushrooms, and cauliflower, respectively. Relative to LSTM, the improvements in MAE are 58.12%, 75.17%, and 61.47%, and in RMSE are 56.48%, 77.04%, and 60.98% for Chinese chives, oyster mushrooms, and cauliflower, respectively. Compared to EEMD–LSTM, the improvements in MAE are 21.99%, 63.85%, and 42.81%, and in RMSE are 19.27%, 70.63%, and 37.14% for Chinese chives, oyster mushrooms, and cauliflower, respectively. Lastly, compared to VMD–LSTM, the improvements in MAE are 16.53%, 16.71%, and 2.98%, and in RMSE are 16.56%, 21.22%, and 2.49% for Chinese chives, oyster mushrooms, and cauliflower, respectively.

When the forecasting step is 3, the model proposed in this study, VMD–EEMD–LSTM, compared to RF, demonstrates substantial improvements in MAE and RMSE for Chinese chives, oyster mushrooms, and cauliflower. The improvements in MAE are 62.92%, 78.00%, and 59.21%, and in RMSE are 62.67%, 78.12%, and 56.72%, respectively. Compared to ELM, the improvements in MAE are 60.08%, 76.87%, and 57.11%, and in RMSE are 58.73%, 77.82%, and 55.57% for Chinese chives, oyster mushrooms, and cauliflower, respectively. Relative to LSTM, the improvements in MAE are 57.66%, 74.20%, and 57.03%, and in RMSE are 57.16%, 75.82%, and 54.91% for Chinese chives, oyster mushrooms, and cauliflower, respectively. Compared to EEMD–LSTM, the improvements in MAE are 25.74%, 63.58%, and 43.41%, and in RMSE are 22.19%, 70.66%, and 37.18% for Chinese chives, oyster mushrooms, and cauliflower, respectively. Lastly, compared to VMD–LSTM, the improvements in MAE are 10.62%, 8.22%, and 1.76%, and in RMSE are 13.70%, 14.35%, and 1.80% for Chinese chives, oyster mushrooms, and cauliflower, respectively.

By calculating the reduction rates of MAP, MAPE, and RMSE between the proposed model and other comparative models, this study further verifies the superiority of the proposed model over other baseline models. The results indicate that the original agricultural product prices become more predictable after undergoing a single decomposition using VMD and EEMD. The secondary decomposition of VMD residual components effectively reduces data complexity and improves predictive results, bringing them closer to actual values. Simultaneously, by calculating the reduction rates of MAP, MAPE, and RMSE between the proposed model and other baseline models in multi-step forecasting, it is demonstrated that the VMD–EEMD–LSTM model maintains good performance in multi-step forecasting. This reaffirms the effectiveness and superiority of the VMD–EEMD–LSTM model.

3.6. The Comparative Experimental Model Analysis Involves Comparing and Analyzing Various Models

Mohanty [34] et al. employed the DTR (decision tree regressor) model for crop yield and price forecasting. In comparison to the DTR model, the VMD–EEMD–LSTM model holds advantages in effectively handling signal components with different frequencies and amplitudes, showcasing enhanced adaptability to non-stationary and multi-scale time series data. The decompositions provided by VMD and EEMD contribute to a more profound understanding and capture of structural information within time series data. When combined with the LSTM model, it efficiently captures long-term dependencies in time series, resulting in outstanding performance in time series prediction tasks.

On the other hand, the DTR model is well-suited for data with relatively simple nonlinear relationships, offering good interpretability. Therefore, the proposed VMD–EEMD–LSTM model in this study exhibits superior efficacy in addressing predictive challenges associated with complex data sequences.

Zhang [35] et al. employed the VMD–ELM method for forecasting the prices of rice, wheat, and soybean meal, achieving favorable results. However, they overlooked the residual components’ continued complexity post VMD decomposition. It was observed that additional processing, such as secondary decomposition of the residual components, could enhance predictive accuracy when combined with the model. In contrast, the model proposed in this study takes into account the crucial information embedded in the residual components after VMD decomposition. By conducting a secondary decomposition of these residuals, the complexity of the data is significantly reduced, leading to improved predictive results. Additionally, through a comparative analysis of results at different prediction steps, this study effectively demonstrates the forecasting advantages of the VMD–EEMD–LSTM model.

Through comparative analysis, it has been observed that the proposed VMD–EEMD–LSTM model in this study exhibits excellent adaptability in handling non-stationary and multi-scale time series data. Furthermore, after undergoing a secondary decomposition with VMD–EEMD, the model effectively reduces the complexity of the data. The combination of this enhanced feature extraction with the LSTM model’s ability to capture long-term dependencies within time series leads to a significant improvement in predictive accuracy for the hybrid model.

4. Conclusions

To enhance the accuracy of agricultural product price forecasting, this study proposes a combined model prediction method based on secondary decomposition and long short-term memory (LSTM) networks. Through empirical analysis of pork, Chinese chives, shiitake mushrooms, and cauliflower prices, the following conclusions are drawn:

(i)

Based on the decomposition approach, the price sequences of pork, Chinese chives, shiitake mushrooms, and cauliflower are decomposed, and their residual components undergo EEMD secondary decomposition. Then, each subsequence is individually predicted using LSTM models, and the final prediction is obtained by combining the predictions of these sub-sequences. Furthermore, the model proposed in this study maintains high predictive accuracy across different prediction horizons.

(ii)

VMD decomposition yields higher accuracy than EEMD decomposition, but the residual components after VMD decomposition exhibit higher complexity. Performing a secondary decomposition on these residuals results in improved model accuracy.

(iii)

Comparing the predictions of the VMD–EEMD–LSTM model with those of the one-time decomposition model, it is evident that the VMD–EEMD–LSTM model demonstrates stronger predictive capabilities, thereby validating its effectiveness and applicability in agricultural product price forecasting.

Although the combined model constructed in this study demonstrates good predictive capabilities, it relies solely on historical data as model inputs. Agricultural product prices are influenced by various factors, such as temperature fluctuations, holidays, and consumer consumption levels, all of which can affect prices. Therefore, to further enhance model accuracy, integrating these influencing factors with the model proposed in this study could be a valuable avenue for improvement.