Prediction of Grain Porosity Based on WOA–BPNN and Grain Compression Experiment

Abstract

The multi-field coupling of grain piles in grain silos is a focal point of research in the field of grain storage. The porosity of grain piles is a critical parameter that affects the heat and moisture transfer in grain piles. To investigate the distribution law of the bulk grain pile porosity in grain silos, machine learning algorithms were incorporated into the prediction model for grain porosity. Firstly, this study acquired the database by conducting compression experiments on grain specimens and collecting data from the literature. The back propagation neural network (BPNN) algorithm was optimized using three metaheuristic algorithms (genetic algorithm (GA), particle swarm optimization (PSO), and whale optimization algorithm (WOA)). Five machine learning models (GA–BPNN, PSO–BPNN, WOA–BPNN, BPNN, and random forest (RF)) were developed to predict the grain porosity using three input parameters (vertical pressure, grain type, and moisture content). The five models were assessed using four evaluation metrics: coefficient of determination (R²), root mean square error (RMSE), mean absolute error (MAE), and mean absolute percentage error (MAPE), to determine the best porosity prediction model. Finally, the generalization ability of the best prediction model was verified using the results of the grain cell box experiment on wheat piles. The results indicated that the WOA–BPNN model was the best prediction model with an R² value of 0.9542, an RMSE value of 0.0079, an MAE value of 0.0044, and an MAPE value of 1.1467%. The WOA–BPNN model demonstrated strong generalization ability, confirming the feasibility of using this model to predict grain porosity. It also established an expression for the relationship between wheat porosity and the vertical pressure of the grain pile. This study presents a machine learning prediction method for determining the porosity of grain piles. The obtained porosity distribution law serves as a crucial basis for conducting comprehensive multi-field coupling analysis of grain piles and offers theoretical support for safe grain storage.

1. Introduction

Grain is the material basis for human survival and development, and grain security is crucial to national economic development and social stability. Grain storage is an important link to guarantee grain security and plays a significant role in the country’s macroeconomic control [1]. During the storage process, grain is prone to local condensation, heat, mold, and other issues in the grain pile. These problems can lead to significant grain storage losses and pose a threat to grain storage safety [2]. The mechanical ventilation of grain silos is an effective control measure to ensure grain security storage. The complex and variable structure of the pores in grain piles affects the airflow rate, airflow direction, and ventilation efficiency during the ventilation process. Therefore, understanding the structural characteristics of the pores is crucial for comprehending the heat and moisture transfer law of the grain pile and guaranteeing grain security storage [3,4].

The porosity of bulk grain is closely related to the pressure it is subjected to. Under the influence of gravity, the pressure at each point within the grain pile varies due to factors such as internal friction (friction between grains) and external friction (friction between grains and the silo wall). Consequently, the distribution of porosity in the grain pile is non-uniform. This study examined the extent of porosity in the bulk grain under vertical pressure, aiming to integrate it with the bulk pressure of the grain pile in the silo and determine the distribution pattern of porosity in the grain pile within the silo.

Existing research has primarily focused on characterizing the pore structure of porous media materials like sandstone, soil, and coal. However, there have been fewer studies conducted on characterizing the pore structure of grain piles. Thus far, numerous methods have been developed to measure or evaluate the pore structure, which are mainly categorized into experimental and numerical simulation methods. The experimental method is the most intuitive and effective approach for directly analyzing the physical properties of porous media structures and determining their porosity on a macroscopic scale. Early experimental methods include mercury intrusion porosimetry [5] and the nitrogen adsorption method [6,7]. Commonly used methods for measuring the pore structure of bulk materials include the gas pycnometer method [8] and the gas displacement method [9]. In 1967, Thompson et al. [8] measured the porosity of oats, wheat, grain sorghum, soybeans, rye, barley, alfalfa, and corn using the gas pycnometer method. In 2020, Mai et al. [9] developed a porosity measuring instrument for the constant volume of bulk materials based on the a and the gas displacement method. This instrument achieved automatic control of the measurement process, real-time display, and storage of the detection parameters as well as a measurement error of 0.3%. The system detection demonstrated high stability and accuracy. With advancements in measurement technology, researchers are now able to physically image porous media using high-precision scanners. Neethirajan et al. [10] utilized CT scanning imaging to reconstruct three-dimensional images of wheat and pea grain piles under simulated silo conditions. They quantitatively described the pore network structure of the grain porous medium. Ge et al. [11] investigated the impact of vertical pressure on the pore structure of grain piles using grain slice imaging and image processing techniques. The results indicated that the pore area and pore diameter of the piled grain gradually decreased with increasing vertical pressure. Currently, numerous researchers are extensively utilizing numerical simulation methods to investigate porous media and model them with the assistance of rapidly advancing computer technology [12,13,14].

In recent years, due to numerous breakthroughs in computer science, many scholars have utilized artificial intelligence to predict issues in various engineering fields [15,16,17,18]. Machine learning techniques have been applied to solve numerous practical problems in various aspects of grain security and storage. Chen et al. [19] utilized machine vision and back propagation neural network technology to classify five corn varieties. The experimental results demonstrated an average classification accuracy of up to 90%. Jia et al. [20] proposed a hyperspectral classification model for detecting corn seed mildew with an accuracy of 92%. This model offers technical support and new insights for early mildew detection as well as for the selection and breeding of corn seeds. Duan et al. [21] proposed a grain temperature prediction model based on support vector regression. They analyzed the correlation between the grain temperature data and the corresponding weather forecast data, conducted principal component analysis on the weather data, reduced the dimensionality of the model input, and utilized the weather information to predict the average temperature of the grain piles. The predicted root mean square error (RMSE) was 1.26 °C. Duysak et al. [22] utilized the KNN method to estimate the quantity of grains in the silo, yielding practical results with an accuracy of 96.71%. Meanwhile, Wu et al. [23] introduced a swift approach for evaluating the quality of maize kernels through image analysis and support vector machine, achieving a recognition accuracy of 97.44%. Zeymer et al. [24] evaluated five different prediction models—cubist, random forest, artificial neural network, support vector machine with polynomial kernel, and support vector machine with radial kernel—using machine learning techniques to simulate the loss of dry matter in stored soybeans. The results indicated that the random forest model yielded the most accurate predictions. Anami et al. [25] utilized a machine learning approach to detect adulteration in bulk paddy. They used the color and texture features of the image and developed a BPNN classification model trained on reduced features derived from the principal component analysis. This model achieved an adulteration classification accuracy of 93.31%. Several studies have investigated the application of machine learning for pest detection, mold detection, and monitoring storage temperature in grain storage. However, there is limited research on the use of machine learning for predicting grain porosity.

In this paper, a machine learning algorithm was incorporated into the prediction model for grain cell porosity. Three influencing factors of grain porosity—grain type, moisture content, and vertical pressure—were selected as the input variables. A self-developed experimental device for the compression of grain specimens was used to obtain experimental data on the porosity of grain piles under vertical pressure. Three hybrid algorithms, GA–BPNN, PSO–BPNN, and WOA–BPNN as well as two single algorithms, BPNN and random forest, were employed to predict the porosity. The models were assessed using various evaluation indices to identify the optimal grain porosity prediction model. The generalization ability of the models was ultimately verified using the grain cell box experiment.

2. Porosity Experiment Determination Method

2.1. Materials

Four common grains found in granaries were chosen for this experiment: wheat, corn, soybean, and paddy. The initial wet basis moisture content of the grains was measured as 10.70%, 11.58%, 10.14%, and 12.93%, respectively, after drying at 105 °C for 72 h, in accordance with the ASAE standards [26]. Broken grains, immature seeds, and impurities were filtered out using a sieve. Then, 100 seeds were randomly selected for measuring their triaxial dimensions with vernier calipers for all the filtered grains [27], and their mean values were calculated. According to the code for the Inspection of Grain and Oils—Determination of Test Weight [28] and for the Inspection of Grain and Oil—Determination of Relative Density of Grain and Oilseeds [29], the material parameters of each type of grain were measured, as shown in

Table 1. Material parameters for each grain type.

Grain Type	Initial Bulk Density ${\mathit{\rho }}_0$/(kg/m3)	Particle Density ${\mathit{\rho }}_{\mathit{s}}$/(kg/m3)	Grain Length L/mm	Grain Width W/mm	Grain Thickness T/mm
Wheat	806.8	1383.75	6.53	3.41	2.91
Corn	735.0	1256.68	11.84	8.62	4.72
Soybean	717.1	1219.87	7.72	7.58	7.46
Paddy	599.0	1164.49	7.38	3.07	2.13

Note: Initial bulk density ${\rho }_0$ is the density of the grain pile (cell) when it is not under vertical load; particle density ${\rho }_s$ is the ratio of the mass of a single grain particle to its volume.

.

This paper investigated the effects of grain type, moisture content, and vertical pressure on grain porosity. Distilled water was sprayed to formulate three different moisture contents for each grain. The samples were sealed in polyethylene bags and stored at 4 °C for 48 h to allow the moisture to equilibrate. The moisture content of various grain samples is presented in

Table 2. Moisture content of various grain samples.

Grain Type	Initial Moisture Content MC0/%	Actual Moisture Content after Formulation MC/%
Wheat	10.70	10.70
		11.26
		12.89
Corn	11.58	11.58
		12.34
		13.47
Soybean	10.14	8.58
		10.14
		13.43
Paddy	12.93	11.69
		12.93
		14.52

. The purpose of formulating multiple sets of moisture contents for each grain is to simulate the occurrence of different moisture contents in the grain in the silo. Different moisture contents will affect the porosity of the grain.

2.2. Experimental Apparatus

The experimental device used in this study for compressing grain specimens is depicted in Figure 1. It was a modified version of the lever-type consolidator commonly utilized in the geotechnical field, capable of applying a maximum vertical pressure of 1200 kPa. The specimen loading chamber as constructed out of high-strength aluminum alloy and Plexiglas, with internal dimensions of 120 mm × 120 mm × 50 mm. The metal loading plate had high stiffness, and a flexible rubber mat was attached to its underside [30]. A water inlet was opened at the bottom of the side wall of the specimen loading chamber, which was connected to the volumetric cylinder via the water injection pipe, and a valve was installed on the water injection pipe. During the experiment, the vertical displacement was measured using a dial indicator. The pressure transducer was positioned at the bottom of the specimen to measure the vertical pressure at that location.

2.3. Experimental Principles

Porosity was determined using the water injection method [31], where the absorption of water by the grain particles can be ignored due to the short injection time (water injected within 1 min).

The equation for the initial grain porosity is:

$$n_0=\left(1-\frac{{\rho }_0}{{\rho }_s}\right)\times 100\%$$

(1)

where $n_0$ is the initial porosity of the grain, ${\rho }_0$ is the initial density of the grain specimen in kg/m³, and ${\rho }_s$ is the density of the grain particles in kg/m³.

The equation for calculating the porosity of grain specimens during loading is:

$$n_i=\frac{\Delta V_i}{A\left(H_0-\Delta H_i\right)}\times 100\%$$

(2)

where $n_i$ is the porosity of the grain specimen after loading at level $i$; $\Delta V_i$ is the total volume change of water in the volumetric cylinder after loading at level $i$, in cm³; $A$ is the bottom area of the specimen loading chamber, in cm²; $H_0$ is the initial height of the grain specimen, in cm; and $\Delta H_i$ is the total compression of the grain specimen after loading at level $i$, in cm.

2.4. Experimental Program

Under vertical loading, the friction on the wall surface of the specimen loading chamber and the grain contact surface is significant. In order to achieve triaxial compression of the specimens, a thin layer of petroleum jelly was applied to the inner wall of the specimen loading chamber before the experiment to reduce the impact of friction. Subsequently, 450 g of the samples was randomly selected and loaded into the specimen loading chamber [32]. The specimens in the loading chamber were loaded at five different pressure levels: 0 kPa, 21 kPa, 42 kPa, 83 kPa, and 125 kPa. Once the dial indicator and pressure transducer readings stabilized, water was injected, and the changes in water volume in the volumetric cylinder and the compression of the grain specimens were recorded. The porosities were then calculated for each pressure level using Equation (2). Each grain was subjected to three groups of experiments with varying moisture contents, resulting in a total of 12 experiments.

3. Database

3.1. Data Sources

In this study, a total of 12 experimental groups were established based on different grain types and moisture content. According to the range of vertical pressure values of the grain pile in the grain silo shown in Figure 2, vertical pressures of 0 kPa, 21 kPa, 42 kPa, 83 kPa, and 125 kPa were selected in five levels to load each group of specimens step by step, and a total of 60 sets of data were obtained. Due to space limitations,

Table 3. Some results of the compression experiments conducted on grain specimens under different vertical pressures.

Grain Types	Moisture Content/%	Vertical Pressure/kPa	Porosity
Wheat	10.70	0	0.414
Wheat	10.70	21	0.405
Wheat	10.70	42	0.398
Wheat	10.70	83	0.385
Wheat	10.70	125	0.377
Corn	11.58	0	0.399
Corn	11.58	21	0.391
Corn	11.58	42	0.384
Corn	11.58	83	0.374
Corn	11.58	125	0.367
Soybean	8.58	0	0.373
Soybean	8.58	21	0.363
Soybean	8.58	42	0.353
Soybean	8.58	83	0.330
Soybean	8.58	125	0.306
Paddy	12.93	0	0.487
Paddy	12.93	21	0.480
Paddy	12.93	42	0.474
Paddy	12.93	83	0.464
Paddy	12.93	125	0.456

presents some of the results of the compression experiments on grain specimens under various vertical pressures.

According to the Code for Design Reinforced Concrete Silos [33], the vertical pressure of a deep silo (defined as H/D_n ≥ 1.5, where H is the silo height and D_n is the silo diameter) is:

$$p_{\mathrm{v}}=\frac{\gamma R\left(1-e^{-\mu Kz/R}\right)}{\mu K}$$

(3)

The vertical pressure of a shallow silo (defined as H/D_n < 1.5) is:

$$p_v=\gamma z$$

(4)

where $p_v$ is the vertical pressure, $\gamma$ is the gravity density of the storage material, $R$ is the hydraulic radius of the horizontal net section of the silo (for circular silos $R=D_{\mathrm{n}}/4$), $\mu$ is the friction coefficient between the storage material and the silo wall, $K$ is the coefficient of lateral pressure, and $z$ is the depth of the storage material.

According to Equations (3) and (4), the relationship between the vertical pressure of grain in the silo and the depth of the pile is illustrated in Figure 2.

The database for this study collected data from the papers published by Zhang [34], Wang [35], Ge et al. [31], Zhou et al. [36], and Li et al. [37], in addition to data from the experiments. Due to the limitation of the length of the article, the literature-collected data are not listed in this paper. Figure 3 shows a three-dimensional schematic of the relationship between vertical pressure, water content, and porosity. A total of 400 sets of data were obtained through the collection of the experimental and literature data. The data in the database were divided into two parts in a 7:3 ratio, with 280 sets of data used for training and the remaining 120 sets used for testing.

3.2. Data Pre-Processing

During the literature data collection process, variations in the research content and experimental settings among authors as well as the constraints of literature length may result in missing data, differences in magnitude, and varying quantitative formats. These discrepancies can directly or indirectly affect the outcome of model training. Before establishing the grain porosity prediction model, the sample data needed to be pre-processed to facilitate data extraction and improve the accuracy of the model. The sample data were normalized to [0, 1] before model training, and the normalization formula [38] is:

$$x_n=\frac{x_{ij}-\underset{k}{\min }\left(x_{kj}\right)}{\underset{k}{\max }\left(x_{kj}\right)-\underset{k}{\min }\left(x_{kj}\right)}$$

(5)

where $x_n$ is the normalized sample data, $x_{ij}$ is the sample data, $\underset{k}{\max }\left(x_{kj}\right)$ is the maximum value in the sample dataset, and $\underset{k}{\min }\left(x_{kj}\right)$ is the minimum value in the sample dataset.

4. Modeling Prediction of Grain Porosity

4.1. Back Propagation Neural Network

The back propagation neural network (BPNN) was proposed in 1986 by a team of scientists led by Rumelhart and McClelland as a feedforward neural network with a backpropagation of errors. The BPNN consists of three parts: the input layer, output layer, and hidden layer. The number of layers and neurons in the hidden layer can be determined by the dimensions of the input vectors and output vectors. In general, BP neural networks are typically set up with a single hidden layer [39]. BP neural networks can learn and memorize a large number of input-output mapping models without the need for prior mathematical formulas describing such mappings. During the process of network learning, the signal involves both forward propagation and backward propagation. If the actual output of the output layer differs from the expected value, the error is propagated backward. The learning rule employs the steepest descent method, continuously updating the weights and biases of the neural network through iterative operations to minimize the network error. A typical three-layer BPNN topology is illustrated in Figure 4.

The output of the hidden layer of the BPNN is given as:

$$r_j=f_1({\displaystyle \sum _{i=1}^m{\omega }_{ij}}{x}_i+b_1)$$

(6)

where $f_1$ is the activation function of the hidden layer, $x_i$ is the value of the neuron in the input layer, ${\omega }_{ij}$ is the connection weights, and $b_1$ is the bias.

The output of the output layer of the BPNN is provided as:

$$y_k=f_2({\displaystyle \sum _{j=1}^l{\omega }_{jk}}{r}_j+b_2)$$

(7)

where $f_2$ is the activation function of the output layer, $r_j$ is the value of the hidden layer neurons, ${\omega }_{jk}$ is the connection weights, and $b_2$ is the bias.

4.2. Whale Optimization Algorithm

The whale optimization algorithm (WOA) is a novel meta-heuristic algorithm introduced by Mirjalili and Lewis in 2016. It draws inspiration from the hunting behavior of humpback whales [40]. Humpback whales feed by spitting bubbles underneath the fish to create a bubble net, trapping the target fish inside before feeding on them. The predation process of humpback whales is divided into three main stages: encircling the prey, bubble attack, and feeding on the prey.

(1)

Encircling Prey

Humpback whales can recognize the location of their prey and surround them. The humpback whale group does not know the exact position of the fish in advance, and their communication enables them to share information about the fish’s location. The WOA initially assumes that the humpback whale closest to the prey is the best candidate (i.e., the current optimal solution). The other humpback whales then calculate the distance between their current position and the current optimal solution, and move toward the optimal solution position. This behavior is represented by the following mathematical model [40]:

$$D=\left|C·X^{\ast }(t)-X(t)\right|$$

(8)

$$X(t+1)=X^{\ast }(t)-A·D$$

(9)

where t is the number of iterations; A and C are coefficient vectors; D is the distance between the humpback whale and its prey; X^∗ is the position vector of the optimal solution in the current humpback whale population; and X is the position vectors of humpback whales other than the optimal solution.

The formulas for A and C in Equations (8) and (9) are as follows:

$$A=2a·r-a,\hspace{1em}a=2-\frac{2t}{T_{\max }}$$

(10)

$$C=2·r$$

(11)

where a is the convergence factor, which decreases linearly from 2 to 0 during the iteration; r is a random vector in the range of [0, 1]; and T_max is the maximum number of iterations.

(2)

Attacking by Bubble Net

After encircling the prey, humpback whales spit bubbles toward the school, creating a bubble net to trap the fish. Humpback whales use a narrowing encirclement or spiral updating of their position to trap food inside a bubble net, optimizing their foraging strategy in the solution space. There are two approaches to the mathematical modeling of humpback whale bubble net attacks:

Method 1: Shrinking encircling mechanism. This behavior is achieved by decreasing the value of a in Equation (10). When |A| < 1, individual humpback whales are approaching the humpback whale in the current best position. Under the condition |A| < 1, humpback whales use larger spanning intervals to swim and search for fish when |A| is larger, enabling the WOA to achieve global optimization. Conversely, when |A| is smaller, humpback whales use smaller spacing to swim, enhancing the WOA’s capacity for local optimization.

Method 2: Spiral updating position. First, calculate the distance between the humpback whale and its prey, and then create a spiral equation that represents the relationship between the whale and its prey. This approach can significantly enhance the local optimization capability of the WOA, albeit at the expense of reduced optimization efficiency. The mathematical model that simulates the humpback whale’s spiral update position method is:

$$X(t+1)=D^{\prime }·e^{bl}·\cos (2\pi l)+X^{\ast }(t)$$

(12)

where $D^{\prime }=|X^{\ast }(t)-X(t)|$ is the distance from the whale to the prey, b is a constant that defines the shape of the logarithmic spiral, and l is a random number within the range of −1 to 1.

Humpback whales swim in a shrinking, encircling circle around their prey while following a spiral path. To simulate this simultaneous behavior, it is assumed that the whale’s position is updated during the optimization process by selecting the shrinking encircling mechanism with a 50% probability and swimming along the spiral path with a 50% probability. The mathematical model is as follows.

$$X(t+1)=\left\{\begin{array}{ll}{X}^{\ast }(t)-A·D& \mathrm{if}\hspace{1em}p<0.5\\ D^{\prime }·e^{bl}·\cos (2\pi l)+X^{\ast }(t)& \mathrm{if}\hspace{1em}p\ge 0.5\end{array}\right.$$

(13)

where p represents a random number in [0, 1].

(3)

Searching for Prey

Humpback whales can randomly search for prey in all directions based on the locations of other whales. The purpose of the search is to identify new locations by using the locations of other humpback whales, as mathematically modeled below:

$$D=\left|C·X_{rand}-X\right|$$

(14)

$$X(t+1)=X_{rand}-A·D$$

(15)

where X_rand is the position vector of the randomly selected humpback whale and D_rand is the distance from the randomly selected humpback whale to the prey.

4.3. Construction of Grain Porosity Prediction Model Based on WOA–BPNN

The BPNN is essentially a gradient descent method, which offers significant advantages such as self-adaptation, self-learning, nonlinearity, and strong self-organization ability. However, the BPNN is highly reliant on the initialization parameters during operation. The objective function to be optimized is very complex, which leads to slow convergence and makes it easy to fall into a local optimum.

The WOA is introduced into the BPNN algorithm to address the numerous defects encountered when the BPNN runs independently. The WOA leverages the fast convergence speed and strong global search ability to optimize the weights and biases of the BPNN. In this approach, the weights and biases are treated as the position information of humpback whales and are updated using the whales’ ability to search for prey. This results in the updating of the weights and biases of the BPNN. After numerous iterations, the optimal global weights and biases are obtained, optimizing the BPNN.

The flowchart depicting the grain porosity prediction model based on WOA–BPNN is shown in Figure 5. The specific establishment process is outlined below.

(1) Initialization of a BPNN. Normalize the sample data, determine the topology of the BPNN based on the empirical formula, and initialize the initial weights and biases of the BPNN. The empirical formula is:

$$H_{hidden}=\sqrt{I_{input}+O_{output}}+R$$

(16)

where H_hidden is the number of nodes in the hidden layer, I_input is the number of nodes in the input layer, O_output is the number of nodes in the output layer, and R is an integer whose value range is [1, 10].

(2) Initialize the WOA. In step (1), convert the initial weights and biases into position vectors of WOA. Then, initialize other fundamental parameters of the whale optimization algorithm such as setting the population size (N), the maximum number of iterations (T_max), the initial minimum weight (w₁), the initial maximum weight (w₂), and the convergence factor (a).

(3) Calculate the individual fitness scores. Identify the whale with the highest current fitness value, record its position vector, and use it as the current optimal individual position X^∗(t).

(4) Update the position of each humpback whale in the population. Different position updating strategies are adopted according to the value of A. If |A| ≥ 1, update the position according to Equation (9); if |A| < 1, update the position according to Equation (13).

(5) The optimization algorithm terminates after reaching the maximum number of iterations or meeting the error accuracy requirement. At this point, the current optimal parameters are assigned to the BPNN.

(6) In the BPNN training, the outputs of each layer are utilized to compute the training error and to adjust the weights and biases of each layer. This process is repeated until the maximum number of iterations is reached or the error falls below the bias.

4.4. Model Evaluation Indicators

After establishing the model, four evaluation indices were selected to test the accuracy and quality of the prediction model. These included the coefficient of determination (R²), root mean square error (RMSE), mean absolute error (MAE), and mean absolute percentage error (MAPE). These indices were used to evaluate the performance of the grain porosity prediction model in this paper. Among them, R² is usually used to evaluate the linear fit of the model. The closer the value of R² is to 1, the better the fit of the model. The lower the RMSE, MAE, and MAPE, the greater the prediction accuracy of the model. The corresponding formulas are shown below:

$${\mathrm{R}}^2=1-\frac{{\displaystyle \sum _{i=1}^m{\left(Y_i-\widehat{Y_i}\right)}^2}}{{\displaystyle \sum _{i=1}^m{\left(Y_i-\overline{Y_i}\right)}^2}}$$

(17)

$$\mathrm{RMSE}=\sqrt{\frac{{\displaystyle \sum _{i=1}^m{\left(Y_i-\widehat{Y_i}\right)}^2}}{m}}$$

(18)

$$\mathrm{MAE}=\frac{{\displaystyle \sum _{i=1}^m\left|Y_i-\widehat{Y_i}\right|}}{m}$$

(19)

$$\mathrm{MAPE}=\frac{1}{m}{\displaystyle \sum _{i=1}^m\left|\frac{Y_i-\widehat{Y_i}}{Y_i}\right|}$$

(20)

where m is the total number of samples, i is the number of sample sequences, $\overline{Y_i}$ is the average of the actual values of the samples, $Y_i$ is the actual value of the samples, and $\widehat{Y_i}$ is the predicted value of the samples.

5. Analysis of Prediction Results

5.1. Comparison of Porosity Prediction Results from Five Machine Learning Models

The prediction of grain porosity involves a regression problem, and there are numerous machine learning algorithms that can be utilized. In order to construct a machine learning model with high prediction accuracy, this paper compared three hybrid models (GA–BPNN, PSO–BPNN, and WOA–BPNN) with two single models (BPNN and RF). Figure 6 displays the predicted and actual values of the BPNN, RF, GA–BPNN, PSO–BPNN, and WOA–BPNN models based on the test set. Figure 7 displays the absolute error between the predicted and actual values of the five machine learning models as determined by the test set. From Figure 6, it is evident that there was a significant difference between the predicted and actual values of the BPNN and RF models, while the predicted values of the three hybrid models GA–BPNN, PSO–BPNN, and WOA–BPNN were closer to the actual values. From Figure 7, it can be observed that the absolute errors of the three hybrid models were smaller than those of the two single models.

Table 4. Computational efficiency of the five machine learning models.

Model	Running Time/Second
BP	1.4812
RF	0.2570
GA–BP	7.5487
PSO–BP	4.0655
WOA–BP	37.6912

shows the computational efficiency of the five machine learning models. From Table 4, it can be seen that the random forest model had the shortest computational time (0.2570 s), while the WOA-BP model had the longest computational time (37.6912 s).

Table 5. Values of the evaluation metrics for the five machine learning models.

Model	R2	RMSE	MAE	MAPE/%
BPNN	0.8874	0.0129	0.0088	2.2208
RF	0.8572	0.0155	0.0094	2.3517
GA–BPNN	0.9450	0.0081	0.0049	1.2735
PSO–BPNN	0.9299	0.0102	0.0063	1.6213
WOA–BPNN	0.9542	0.0079	0.0044	1.1467

presents the values of the four evaluation metrics (R², RMSE, MAE, and MAPE) for the five machine learning models. From Table 5, it is evident that the WOA–BPNN model exhibited the highest R² value and the lowest RMSE, MAE, and MAPE values, while the RF model showed the lowest R² value and the highest RMSE, MAE, and MAPE values. Based on the values of the four evaluation metrics, it can be observed that the prediction accuracies of the five machine learning models were ranked from highest to lowest as WOA–BPNN, GA–BPNN, PSO–BPNN, BPNN, and RF. Combining Table 4 and Table 5, it can be observed that while the computational efficiency of the RF model was high, the predictive performance of this model was poor. The computational efficiency of the WOA-BP model was lower than that of the other four models, but the prediction accuracy of the WOA-BP model was high. This high accuracy could meet the requirements for predicting grain porosity in this study. Figure 8 displays the radar charts of the four evaluation metrics, providing a more intuitive way to observe the prediction effect of each model. It can also be observed from Figure 8 that the WOA–BPNN model yielded the most accurate predictions.

5.2. Performance Evaluation of WOA–BPNN Model for Porosity Prediction

Figure 9 illustrates the optimal validation performance of the WOA–BPNN. From Figure 9, it can be seen that the WOA–BPNN model terminated after 80 iterations, with the best validation performance occurring at the 74th iteration. Figure 10 depicts the training status of the WOA–BPNN model. From Figure 10, it can be seen that the training phase stopped after 80 iterations of the WOA–BPNN model, and validation was performed six times. The regression results are depicted in Figure 11, and the correlation coefficient R measures the relationship between the target value and the output value. A larger R indicates a stronger prediction ability of the WOA–BPNN model. As shown in Figure 11d, the correlation coefficient (R) between the actual porosity values in the database and the predicted values of the WOA–BPNN model was 0.9828. This indicates that the WOA–BPNN model demonstrated strong predictive capability.

5.3. Experimental Validation of the WOA–BPNN for Porosity Prediction Model

5.3.1. Grain Cell Box Experiment

To validate the effectiveness of the established WOA–BPNN model in practical applications, the experimental data from the grain cell box experiment were input into the WOA–BPNN model. Subsequently, the predicted porosity values from the WOA–BPNN model were compared with the experimental values. The grain cell box experiment is depicted in Figure 12, with the internal dimensions measuring 600 mm × 600 mm × 600 mm. During the experiment, vertical pressure is applied to the storage material inside the box by a rubber airbag, and an external air pump provides stable gas pressure for the airbag. A 2-mm thick plexiglass plate is placed between the material and the rubber airbag, and the sidewalls around the airbag are coated with a thin layer of petroleum jelly to minimize friction. This setup ensures that the bottom of the airbag can evenly distribute the load on the top surface of the entire grain pile.

The grain cell box experiment was conducted to determine the porosity of wheat as a storage material. The wet base moisture content of the wheat pile was 10.89%, the initial porosity (the porosity of the wheat in the natural stacking state at a vertical pressure of 0 kPa) was 0.429, and the particle density was 1250 kg/m³. The porosity of wheat was measured under vertical pressures of 0 kPa, 25 kPa, 50 kPa, 75 kPa, 100 kPa, 125 kPa, and 150 kPa, resulting in seven sets of experimental data for validation.

5.3.2. Comparison of Predicted and Experimental Results

The experimental data on wheat piles from the grain cell box experiment were input into the developed WOA–BPNN model for prediction. A comparison between the model’s prediction results and the experimental results is presented in Figure 13. As depicted in the figure, the porosity predicted by the WOA–BPNN model was slightly lower than the porosity in the grain cell box experiment. Both the predicted and experimental porosity values gradually decreased with increasing vertical pressure. The model prediction and experimental values were fitted to obtain two curves. The trend of these two curves was essentially the same, and the difference between the two in the longitudinal coordinates was approximately constant. The initial porosity predicted by the WOA–BPNN model and the initial porosity measured in the grain cell box experiment were 0.418 and 0.429, respectively, differing by only 2.56%. When the vertical pressure was 150 kPa, the WOA–BPNN model predicted a porosity of 0.373, while the porosity measured by the grain cell box experiment was 0.389, resulting in a difference of only 4.11%.

According to Equations (17) to (20), the R² between the predicted values of the WOA–BPNN model and the experimental values of the grain cell box was found to be 0.9928. The RMSE was 0.0139, the MAE was 0.0137, and the MAPE was 3.3917%. These results indicate that the WOA–BPNN model in this study had high prediction accuracy and a good training effect. It also indicates that the error between the model’s predicted values and the experimental values was within the permissible range. This suggests that the WOA–BPNN model can be used to predict the porosity in the grain cell box, demonstrating its strong generalization ability. The feasibility of using the WOA–BPNN model to predict the porosity of a grain pile was confirmed by comparing the results of the WOA–BPNN model with those of the grain cell box experiment.

According to the variation rule depicted in Figure 13, the data were fitted using OriginPro version 2021 software to determine the relationship between the porosity of the wheat pile and the vertical pressure. The fitting formula is:

$$n=n_0+a{p}_v^b$$

(21)

where $n$ is the porosity of the grain cell, $p_v$ is the vertical pressure on the top surface of the grain cell in kPa, $n_0$ is the porosity of the grain cell when the vertical pressure is 0 kPa, and a and b are the parameters to be fitted. The parameters obtained by fitting are shown in

Table 6. Values of each fitted parameter in the equation for wheat pile porosity.

Methods for Obtaining Porosity	a	b	n0	R2 (Curve Fitting)
WOA–BPNN model	−0.0011	0.7509	0.418	0.9968
Grain cell box experiment	−0.0014	0.6736	0.429	0.9973

.

From Table 6, it can be seen that the parameter a values of the wheat heap porosity expressions obtained from the WOA–BPNN model and the grain cell box experiment were approximately equal, and the parameter b values did not differ much. This indicates that the fitting curves obtained from these two methods had a similar trend of change. The results demonstrate that the WOA–BPNN model can accurately predict the porosity of the larger-sized (600 mm × 600 mm × 600 mm) grain cell under varying vertical loads.

6. Conclusions

In this paper, three hybrid models (GA–BPNN, PSO–BPNN, and WOA–BPNN) and two classical models (BPNN and RF) were established to predict the porosity of grains using machine learning algorithms with data from the experiments and collected data. The prediction results of the five models were compared using four evaluation indices. The generalization ability of the best prediction model was verified by conducting the grain cell box experiment. The following conclusions can be drawn:

(1) All three hybrid models in this study significantly outperformed the BPNN model and the random forest model in predicting porosity.

(2) The prediction accuracies of the five machine learning models were, in descending order, the WOA–BPNN model, GA–BPNN model, PSO–BPNN model, BPNN model, and RF model. Among them, the WOA–BPNN model exhibited the best overall performance among the three hybrid models and two single models.

(3) The generalization ability of the WOA–BPNN model was verified through the grain cell box experiment, and the model’s prediction results were largely consistent with the experimental results. This indicates that the model can be used to predict the porosity of a grain cell with a larger size (600 mm × 600 mm × 600 mm), demonstrating the strong generalization ability of the WOA–BPNN model.

(4) Based on the prediction results of the WOA–BPNN model and the findings of the grain cell box experiment, an equation describing the relationship between the porosity of wheat and the vertical pressure of the grain pile was derived through additional analysis and data fitting.

In this study, the grain porosity prediction model based on the WOA-BP algorithm, obtained through big data training, could provide the porosity of a specific grain type promptly when the water content and vertical pressure are known. Under the condition that the distribution of the bulk pressure of the grain pile in the silo and the water content are known, the distribution pattern of the porosity of the grain pile in the silo can be determined over time.

The porosity distribution law of the grain pile in a grain silo is fundamental to several key research problems. More accurate results or conclusions can be derived by considering the actual distribution law of porosity. In the study of mechanical ventilation in grain silos as well as the temperature field, humidity field, and multi-field coupling of grain piles, previous studies typically treated the grain piles as a homogeneous medium and considered the porosity as a constant value. However, in fact, the porosity in the grain pile is spatially variable. Based on the grain porosity prediction model obtained in this study, we can first determine the distribution pattern of grain pile porosity in the silo. Based on the distribution law of grain pile porosity in the grain silo, more accurate results or conclusions can be obtained for the aforementioned research issues. For instance, analyzing the distribution of wind resistance in the grain pile during mechanical ventilation can help guide grain custodians in conducting mechanical ventilation to enhance ventilation efficiency and effectiveness.