Experimental and Machine Learning Approach to Investigate the Mechanical Performance of Asphalt Mixtures with Silica Fume Filler

Abstract

This study explores the potential in substituting ordinary Portland cement (OPC) with industrial waste silica fume (SF) as a mineral filler in asphalt mixtures (AM) for flexible road pavements. The Marshall and indirect tensile strength tests were used to evaluate the mechanical resistance and durability of the AMs for different SF and OPC ratios. To develop predictive models of the key mechanical and volumetric parameters, the experimental data were analyzed using artificial neural networks (ANN) with three different activation functions and leave-one-out cross-validation as a resampling method. The addition of SF resulted in a performance comparable to, or slightly better than, OPC-based mixtures, with a maximum indirect tensile strength of 1044.45 kPa at 5% bitumen content. The ANN modeling was highly successful, partly due to an interpolation-based data augmentation strategy, with a correlation coefficient R_CV of 0.9988.

1. Introduction

Asphalt mixtures (AM) are the most used construction material for flexible road pavements worldwide; they are designed to withstand traffic loading, ageing, and water damage [1]. Such mixtures of coarse aggregates, sand, and bitumen are integrated by mineral filler, which leads to porosity reduction and bitumen stiffening, resulting in improvement of mechanical response against rutting and shoving [2,3,4,5,6,7]; however, overdoses of mineral filler reduce the aggregate content, and a reduction in strength is observed [8].

The effective utilization of various mineral fillers such as lime [9], hydrated lime [10], limestone powder [11], fly ash [12] and cement [13] to improve the strength of asphalt mixtures has been investigated in various studies. In particular, Portland cement has demonstrated better effectiveness as a filler material [14]; however, on account of high CO₂ emissions associated with cement manufacturing, researchers are constantly looking for alternate environment-friendly filler materials [15].

Silica fume (microsilica), an industrial waste product, has an adverse impact on the environment if not disposed of safely; it is an extremely fine silica nanocrystalline polymorph formed as a by-product of silicon and ferrosilicon alloy processing in electric arc furnaces (EAFs). Silica fume (SF) has been used in the concrete industries to improve compressive strength, hardening of concrete, and as a partial replacement for cement [16,17,18,19]. In fact, the material is involved in the pozzolanic reaction with the calcium hydroxide crystals, producing additional cementing material (C-S-H) and eliminating areas of stress concentrations prone to failure initiation [17,20,21,22]. Furthermore, SF is used in road pavement subgrade to reduce swelling pressure and increase the bearing capacity of the pavement’s subgrade [23,24,25]; however, a minimal investigation in necessary to study the use of this material as filler in AM. Al-Hdabi et al. investigated the effect of silica fume as mineral filler additive with waste fly ash on the indirect tensile strength of AM [26]. The study revealed that the addition of SF enhances the initial mechanical properties, long-term performance, and durability. Abutalib and Abutalib et al. studied the application of silica fume (0.25–1.0%) to improve the rheological properties and to reduce the oxidation rate of the asphalt binder [27,28]. Al-Taher et al. carried out experimental work to quantify the effect of SF on various mechanical properties of asphalt concrete [29]. The effect of SF on the mechanical strength of reclaimed asphalt concrete was investigated by Larbi et al. [30]. Choudhary et al. reported that silica fume can effectively improve the mechanical and rheological properties of asphalt concrete and can be used as an alternative filler material [31].

This research investigates the potential of using SF, a sustainable filler material, in place of OPC in asphalt mixtures for the base course of road pavements. The study focuses on evaluating the impact of SF as a mineral filler for heavy load requirements in the Indian context. The mechanical response of the asphalt mixtures was evaluated using Marshall stability (MS), the mean Marshall stability ratio (MMSR), and indirect tensile strength (ITS) and water sensitivity (WS) tests for different proportions of silica fume, OPC, and bitumen.

As well as providing a comprehensive mechanical and volumetric characterization of asphalt mixtures for base courses made with SF rather than OPC as filler, this study aims to develop a machine learning-based methodology for predicting the mechanical and volumetric parameters of such mixtures. To this end, the experimental data were analyzed using artificial neural networks (ANNs) to create predictive models. Developing accurate predictive equations for the mechanical and volumetric response of asphalt mixes can enhance traditional mix design, which relies on costly and time-consuming laboratory testing. Numerical prediction can reduce the number of experiments required to optimize the mix without sacrificing reliability of the results. Large datasets are usually required to properly train a neural model using conventional approaches [32,33]. However, in this study, a nonconventional approach called leave-one-out cross-validation was used to adequately train the neural network model based on the limited experimental data available, maximizing the reliability of the ANN model. This approach was taken to address the relatively small amount of data in the laboratory dataset.

2. Materials and Methods

2.1. Bitumen

Tiki Tar Industries (Baroda) Limited, Gujarat, India, provided the conventional VG-30 grade bitumen used in this study. The binder’s properties were assessed through viscosity, flash point, penetration, softening point, solubility, ductility, and specific gravity tests.

Table 1. Properties of bitumen.

Test Parameter	MoRTH Limits	Test Results	Test Method
Absolute viscosity (poises) at 60 °C	2400–3600	2855	IS 1206 (P-2)
Kinematic viscosity (cSt) at 135 °C	(Min) 350	392	IS 1206 (P-3)
Flash point Cleveland open cup (°C)	(Min) 250	304	IS 1448 (P-69)
Penetration (1/10 mm) at 25 °C, 100 gm, 5 s	(Min) 45	49	IS 1203
Softening point (°C)	(Min) 47	48	IS 1205
Matter soluble in trichloroethylene (% by mass)	(Min) 99	99.45	IS 1216
Viscosity ratio at 60 °C	(Max) 4.0	1.3	IS 1206 (P-2)
Ductility (cm) at 25 °C after TFOT	(Min) 40	75	IS 1208
Specific gravity (g/cm3)	0.97–1.02	0.987	IS 1202

shows the physical properties of the bitumen used in the study. With respect to the bitumen characterization, strong similarities can be found between the Indian standards (IS 1203 and IS 1205) and European standards (EN 1426 and EN 1427) for the evaluation of penetration and softening point, respectively. Regarding the viscosity measurements, both the Indian (IS 1206 P-2 and P-3) and European standards (EN 12595 and EN 12596) are based on the same capillarity approach.

2.2. Aggregate

The investigation used crushed quartz aggregates, and several laboratory experiments were conducted to evaluate their properties. These experiments included cleanliness, grain size distribution, bulk specific gravity, Los Angeles abrasion, soundness, flakiness and elongation index, impact strength, and water absorption. The physical properties of the aggregates are presented in

Table 2. Properties of aggregate.

Test Parameter	MoRTH Limits	Test Results	Test Method
Cleanliness (Dust) (%)	(Max) 5	3	IS 2386 Part I
Bulk specific gravity (g/cm3)	2–3	2.68	IS 2386 Part III
Percent wear by Los Angeles abrasion (%)	(Max) 35	10.6	IS 2386 Part IV
Soundness loss by sodium sulphate solution (%)	(Max) 12	3.4	IS 2386 Part V
Soundness loss by magnesium sulphate solution (%)	(Max) 18	3.7	IS 2386 Part V
Flakiness and elongation index (%)	(Max) 35		IS 2386 Part I
20 mm		27.93
10 mm		32.13
Impact strength (%)	(Max) 27		IS 2386 Part IV
20 mm		4.15
10 mm		5.91
Water absorption (%)	(Max) 2	1.67	IS 2386 Part III

and Figure 1. The grading curve of the aggregate was compared with the reference envelope prescribed by the Indian Ministry of Road Transport and Highways (MoRTH, 2013) for flexible road pavements, which resulted in the classification of the material as grade II [34].

2.3. Filler

In India, ordinary Portland cement (OPC) is the conventional mineral filler that is typically used in asphalt mixes for heavy-load road pavements. In this study, the effectiveness of silica fume (SF) as a mineral filler was evaluated by comparing it to OPC grade 43. According to the Indian Standard IS 8112-1989, OPC grade 43 must provide a designed strength of at least 43 MPa at 28 days.

The processed silica fume used in the study was obtained in powder form and air-dried from Safew Tech system in Indore (India), while OPC was sourced from a local supplier. Both materials underwent assessment of their properties, including density, specific surface area, grain size distribution, methylene blue value analysis, mineralogical and chemical composition, German filler, and pH value. The German filler test, as prescribed by the National Asphalt Pavement Association (1995) [35], was used to determine the porosity of the fillers, which is inversely proportional to its German filler value. The silica fume grading curve was compared with the required size of mineral filler prescribed by MoRTH [34] and was found to lie within the recommended range for asphalt mixture (

Table 3. Grain size curve of silica fume as per MoRTH requirement for filler material.

Sieve Size	Passing Percentage (%)
	MoRTH Upper Limit	MoRTH Lower Limit	Considered Silica Fume
0.6 mm	100	100	100
0.3 mm	100	95	96.8
0.075 mm	100	85	90.3

).

Table 4. Properties of fillers.

Property	OPC	SF
Specific gravity (g/cm3)	3.04	2.2
MBV (g/kg)	3.00	3.85
German filler (g)	85	94
Fineness modulus (FM)	4.96	1.96
Average particle size (µm)	10.12	0.243
Surface area (m2/g)	1.75	16.45
pH	12.90	6.98
Loss on ignition (%)	1.8	2.0
Particle shape	Granules and sub-angular particles with rough texture	Spherically shaped and very fine
Mineralogical composition (XRD)	Alite, Belite, Calcite, Quartz, Portlandite	Quartz and Calcite
SiO2 (%)	21.43	93.5
CaO (%)	66.58	0.89
Al2O3 (%)	3.01	0.08
MgO (%)	1.39	0.82
SO3 (%)	2.26	-
Na2O (%)	0.03	0.4
K2O (%)	0.02	3.52
TiO2 (%)	0.21	-
MnO (%)	0.20	0.06
Fe2O3 (%)	4.68	-
TiO2 (%)	-	0.05

presents the various physical and chemical properties of silica fume.

2.4. Active and Passive Adhesion

The aggregate–bitumen interface significantly affects the moisture resistance of AM. This aspect has been evaluated in accordance with ASTM D3625 with varying percentages of filler [36]. Active adhesion is the time required to accurately coat the aggregate with bitumen during the mixing operation. Bitumen, filler, and aggregate were preheated at a temperature of 170 °C and then mixed with optimal proportions. The total time elapsed between the addition of the target binder content and 100% coating of the aggregate is measured in seconds [37]. The time required to completely coat the aggregate was noted and reported as active adhesion. The shorter elapsed time represents good active adhesion.

Passive adhesion was calculated by visual inspection using the boiling water test according to ASTM D 3625-96. The passive adhesion was determined through visual inspection. A 250 g AM was prepared and stored at room temperature to cool the mixture down to a temperature of 85–100 °C. The cooled bitumen aggregate mixture is boiled in water for 10 min and after boiling it is cooled again to room temperature. The water is decanted from the boiled mix and then dried on absorbent paper after boiling the aggregate–bitumen mix for 10 min. The percentage of bitumen that remained spread on the aggregate was reported as passive adhesion. The percentage of bitumen that remained coated on the aggregate was studied by a team of 5 technicians using a magnifying glass. The average value obtained for the percentage of bitumen remain coated on aggregates is reported as passive adhesion.

The results of passive and active adhesion are presented in Figure 2. Results show a better active and passive adhesion of silica fume with respect to OPC, which should ensure higher moisture damage resistance to SF bituminous mixes.

2.5. Asphalt Mixes

To isolate the effect of various filler materials on the physical-mechanical properties of asphalt mixes, the study kept aggregate type, grain size distribution, and bitumen type constant across all mixtures. Three different volumes of silica fume (4.0%, 6.0% and 8.0%) and the same amount of OPC were used for comparison. Each dosage of SF and OPC was prepared with 4 different bitumen contents (4.5%, 5.0%, 5.5%, and 6.0% by volume of mix), resulting in 24 distinct asphalt mixtures analyzed.

The aggregate, SF, OPC, and bitumen were heated separately at 150 °C before mixing. Dry samples of the necessary amount of aggregate and fillers were mixed at 170 °C in a planetary mixer with a speed of 140 ± 5 rpm. Then, different percentages of bitumen were added to create asphalt mixture specimens according to Indian MoRTH guidelines. The entire laboratory procedure was performed in triplicate, and the experimental findings were reported with standard deviation.

2.6. Marshall Stability Strength Test

The specimens for the subsequent physical-mechanical tests were prepared using the Marshall method, which is based on the impulsive compaction principle. Although other laboratory techniques, such as the gyratory compaction, are widely used in the USA within the SUPERPAVE method (NAPA, 1995) [38], the Marshall compaction is still commonly used in many road laboratories due to its simplicity, low cost, and extensive database available in the literature, despite the fact that it leads to volumetric properties of the compacted specimens quite different from those observed in road pavement as a result of rollers compaction [8,37,39,40,41,42,43]. Marshall testing was carried out in accordance with ASTM D6927 [44] and the results were reported as Marshall stability, flow, and quotient.

2.7. Volumetric Properties

The desirable physical and mechanical characteristics of an asphalt mixture are closely related to the volume of its binder, mineral filler, and coarse aggregates. Hence, the volumetric properties of asphalt mixtures significantly impact the durability and strength of road pavements. In accordance with Indian Ministry of Road Transport and Highways (MoRTH) standards, various volumetric properties, such as maximum theoretical specific gravity, dry density, percentage of air voids (Vv), voids in mineral filler, and voids filled with bitumen, were determined for both silica fume and OPC mixtures. The experimental procedures involved in MoRTH for the evaluation of the volumetric properties are basically equivalent to the experimental procedures described in EN 12697-6 and EN 12697-8.

2.8. Long-Term Aging (LTA)

LTA is usually considered to evaluate the long-term hardening of asphalt pavement during its in-service life. Three specimens per mixture were cured in an oven at 80 ± 5 °C for 5 days to simulate the long-term hardening produced by 10 years of service life [45] and subsequently tested for Marshall stability, according to the standard ASTM D6927. Before placing the specimen for curing, it was clamped to avoid complete collapse. The long-term aging effect was evaluated in terms of the mean Marshall stability ratio (MMSR). The equation to calculate MMSR is given in Equation (1):

$$MMSR=\frac{MSA}{MS}$$

(1)

where MSA is the Marshall stability after aging and MS is the standard Marshall stability of the mixes. Although long-term aging can be determined using various strength indices, in this study, the Marshall strength was considered as described by Al-Hdabi [46].

2.9. Indirect Tensile Strength Test

In accordance with ASTM D6931, an indirect tensile strength (ITS) test was conducted on cylindrical specimens of the Marshall type after conditioning at 25 °C for 4 h. The test was carried out at a constant strain rate of 50.4 mm/min [47].

2.10. Water Sensitivity Test

The evaluation of water sensitivity was carried out according to the modified Lottman test [48]. WS was expressed by the percentage reduction in ITS evaluated before and after moisture conditioning performed in accordance with AASHTO T-283 standards [49]. The 63.5-mm-thick and 150-mm-diameter cylindrical specimens were prepared using a modified Marshall hammer to obtain an air voids content of 7 ± 0.5%. These compacted specimens were subsequently stored at room temperature for 24 ± 3 h before being moisture conditioned. Specimens were saturated for 10 min at 50 kPa absolute vacuum pressure and then left submerged for 10 additional minutes in a vacuum container. A degree of saturation equal to 75 ± 5% was maintained. Cured specimens were covered by a thin plastic film and subsequently sealed in a plastic bag containing 10 ± 0.5 mL of water. Sealed specimens were then frozen at −18 ± 3 °C for 24 hrs. After removal from the plastic bag, the frozen specimens were placed in a water bath for 24 h at 25 ± 5 °C. Following specimens’ conditioning, ITS was evaluated according to ASTM D6931 [47].

3. Methodology

3.1. Artificial Neural Networks

ANNs are mathematical models falling in the class of non-linear parametric functions [50,51,52,53,54]. Such models are the result of the weighted and biased connection of logistic regression units (called artificial neurons), organized in several sequential layers: an input layer that receives the features vector, one or more hidden layers, and an output layer that produces the network outcome. The connections, which establish the complexity and computing power of the ANN, link neurons belonging to different layers so that the information can travel from the input to the output layer in only 1 direction. Nevertheless, the activity of each neuron is ruled by a non-linear activation function that determines whether the processed output should be transmitted or interrupted. Thanks to this network architecture, ANNs can identify the relationship or pattern that links input predictors to target variables: connection weights and biases are set by a supervised training process that aims to minimize the prediction error of the experimental targets.

This study focuses on shallow neural networks (SNNs), i.e., 3-layer perception networks, which have been shown to solve arbitrarily well any multidimensional input-target fitting problem by providing a sufficient number of neurons in their only hidden layer [55]. The proposed SNN consists of a 3-neuron input layer (1 neuron for each input feature), an $N$-neuron hidden layer whose processed output is passed to a non-linear activation unit, and a 7-neuron output layer associated with a linear activation function. The number of neurons $N$ and the activation function $FUN$ characterizing the hidden layer were identified through a grid search procedure, choosing the combination that yielded the lowest score on the performance function. In this study, the integer range $X_N=\{1,\cdots ,100\}$ was used for selecting the number of hidden neurons and the following activation functions were considered: the exponential linear (ELU), the hyperbolic tangent (TanH), and the logistic sigmoid (LogS) units (Figure 3). The input features were the bitumen content (% by volume of mix), percentage of filler (%), along with a categorical variable that distinguishes the filler materials (values: 0 for Silica Fume, 1 for OPC). The output variables were air void content (%), the voids in the mineral aggregate (%), the Marshall stability before and after aging (kN), the Marshall quotient (kN/mm), the indirect tensile strength (kPa), and the water sensitivity of asphalt mix (%). Each variable belonging to both input and target vectors was standardized before being processed by the network, i.e., all variables were rescaled to have 0 mean and unity standard deviation. Such an adjustment served to reduce computational time and improve the efficiency of the neural model.

3.2. ANN Training and Regularization

The supervised training process identifies connection weights and biases that minimize the difference between the ANN output $\widehat{\mathit{y}}$ and the experimental target $\mathit{y}$, corresponding to the input feature vector $\mathit{x}$. This process is divided into 2 distinct phases: a forward and a backward pass. In the forward stage, the training feature vector $\mathit{x}$ is inputted to the network, and the neuron activations produce the output $\widehat{\mathit{y}}$. After that, a backward comparison is made between the computed output $\widehat{\mathit{y}}$ and the experimental target vector $\mathit{y}$ by a performance function $F\left(\widehat{\mathit{y}},\mathit{y}\right)$, also called loss function, with the aim of defining the corrections to the weights and biases of the network. The training process involves the use of a learning rule that defines the update of network parameters $\mathit{W}$ (the matrix of weights and biases) according to the value assumed by the performance function for a fixed number of iterations $E$. The mean sum of the error squares (or mean square error, MSE) is commonly accepted as $F\left(\cdot \right)$ function and its gradient with respect to $\mathit{W}$, calculated by means of a backpropagation algorithm [52], is used in the learning rule so that the network parameters are updated to minimize the loss value.

For a better understanding of the training process implemented in this study, the analytical expression of the learning rule is presented in Equation (2):

$${\mathit{W}}^{e+1}={\mathit{W}}^e-{\left({\nabla }^{2}F({\mathit{W}}^e)\right)}^{-1}\nabla F({\mathit{W}}^e),e\in \{1,\cdots ,E\}$$

(2)

Such equation, assuming the Levenberg–Marquardt (LM) backpropagation algorithm [56], becomes:

$${\mathit{W}}^{e+1}={\mathit{W}}^e-{\left[{\mathit{J}}^T\left({\mathit{W}}^e\right)\mathit{J}\left({\mathit{W}}^e\right)+{\mu }_{\mathrm{e}}\mathit{I}\right]}^{-1}{\mathit{J}}^T\left({\mathit{W}}^e\right)\mathit{v}\left({\mathit{W}}^e\right),e\in \{1,\cdots ,E\}$$

(3)

where ${\mathit{W}}^e$ is the matrix of weights and biases at iteration $e$, $\mathit{J}$ is the Jacobian matrix of the training loss $F\left(\cdot \right)$ with respect to ${\mathit{W}}^e$, $\mathit{I}$ is the identity matrix, and $\mathit{v}\left({\mathit{W}}^e\right)=\widehat{\mathit{y}}\left({\mathit{W}}^e\right)-\mathit{y}$ is the vector of network errors. ${\mathit{W}}^{e+1}$ are the updated values of network parameters to be used in the forward pass of iteration $e+1$. While the direction towards the minimum is determined by the gradient ${\mathit{J}}^T\left({\mathit{W}}^e\right)\mathit{v}\left({\mathit{W}}^e\right)$, the scalar $\mu$ determines the step size taken in that direction at each iteration and, as a result, the convergence rate. In particular, a high $\mu$ value corresponds to a small learning step and a very long convergence time; on the contrary, a low $\mu$ allows the convergence speed to increase but the algorithm could jump over the minimum. In order to achieve faster convergence and, in parallel, avoid undesirable local minima, the parameter $\mu$ is varied during training: the value ${\mu }_{\mathrm{e}+1}$ corresponds to ${\mu }_{\mathrm{e}}$ multiplied by ${\mu }_{inc}>1$ (or ${\mu }_{dec}<1$) if the performance index has increased (or decreased) between iterations $e-1$ and $e$. In addition, if the parameter $\mu$ becomes too large, ${\mu }_{\mathrm{e}+1}>{\mu }_{\mathrm{m}\mathrm{a}\mathrm{x}}$, the LM algorithm is stopped.

At the end of the $E$ iterations of the training process or when ${\mu }_{max}$ is reached (a good indication that the algorithm has truly converged), the optimal weights and biases are kept fixed while the test feature vector is processed just in the forward manner to define the model’s loss index on novel data.

However, under the conditions set out above, a neural model may overfit the training data: such a situation is generally caused by too large values of connection weights, which increase the slope of the network function. To avoid overfitting and to improve the performance score on test data, a regularization technique was implemented in the current study setup [55]. It involved modifying the sum squared error by adding a term as the sum of squares of the network weights in order to penalize network complexity and force the resulting function to be smooth. In this way, the ANN optimization objective becomes:

$$F_{opt}\left(\widehat{\mathit{y}}\left({\mathit{W}}^e\right),\mathit{y},{\mathit{W}}^e\right)=\beta {‖\widehat{\mathit{y}}\left({\mathit{W}}^e\right)-\mathit{y}‖}_2^2+\alpha {‖{\mathit{W}}^e‖}_2^2$$

(4)

where the operator ${‖\cdot ‖}_2^2$ represents the 2-norm applied to the network’s parameters ${\mathit{W}}^e$ and errors $\mathit{v}\left({\mathit{W}}^e\right)=\widehat{\mathit{y}}\left({\mathit{W}}^e\right)-\mathit{y}$. $\alpha$ and $\beta$ are the regularization parameters that control the complexity of the network solution: the ratio $\alpha /\beta$ assumes values in the interval [0, 1] where the bigger is the ratio and the smoother is the ANN response. However, the regularization raises the issue of setting the $\alpha$ and $\beta$ parameters properly: if the ratio $\alpha /\beta$ is too large compared to the optimum, the network will create an overly smooth interpolation through the training data, and on the contrary, if such ratio is too small, the network will be more likely to overfit the data.

In this study, David MacKay’s approach [57], which involves applying Bayesian statistics to the training of the neural network, was used to optimize the regularization parameters. Therefore, $\alpha$ and $\beta$ were automatically re-estimated every iteration of the training algorithm with the aim of restricting the connection weights to be small by seeking to minimize the loss function of Equation (4).

Conversely, the hyperparameters that define the functioning of the LM algorithm ($\mu$, ${\mu }_{inc}$, ${\mu }_{dec}$, ${\mu }_{max}$ and $E$) were kept at their default values as assigned in the implementation offered by the well-known MATLAB^® ANN Toolbox: the maximum number of training epochs $E$ was set to 1000, $\mu \_dec$ to 0.1, $\mu \_inc$ to 10, the highest $\mu$ value that stops training is 1 × 10¹⁰, and the initial $\mu$ is 0.001. As previously mentioned, the number of neurons $N$ in the SNN hidden layer and the activation unit type $FUN$ were identified by means of a grid search approach that in the current study setup aims to find the lowest averaged test error, penalized by the number of SNN parameters, of 12 fitted networks (with the same structure for each combination of $X_N=\{1,\cdots ,100\}$ and $X_{FUN}=\{ELU,TanH,LogS\}$) as a result of a leave-one-out cross-validation partitioning.

3.3. Leave-One-Out Cross-Validation

The standard practice of splitting the available dataset into 2 random subsets for training and testing (also called hold-out method) may result in biased performance evaluations due to the different distribution of data within such splits (the test error rate can be highly variable depending on precisely which observations are included in each set), along with the risk of missing some relevant trends in training data. These effects are particularly marked when the dataset is relatively small [58]. Conversely, the leave-one-out cross-validation (LOOCV) method [59] repeatedly performs a partition of the experimental dataset into a training set containing all but 1 observation and a test set that contains only that observation (Figure 4).

Therefore, if the available set of observations consists of $n$ data, the ANN model is trained on $n-1$ observations and the fitted model is used to predict the response for the excluded observation. Such procedure is repeated for $n$ splits: the first training set contains all but observation 1, the second training set contains all but observation 2, and so on. Each fitted model produces a squared error $M{SE}_i$ computed as the squared difference between the prediction ${\widehat{\mathit{y}}}_{\mathit{i}}$ and the experimental value ${\mathit{y}}_{\mathit{i}}$ of the excluded observation: ${MSE}_i={‖{\widehat{\mathit{y}}}_{\mathit{i}}-{\mathit{y}}_{\mathit{i}}‖}_2^2$. The test error of the ANN model is then estimated by averaging the $n$ resulting MSE’s:

$${MSE}_{CV}=\frac{1}{n}\sum _{i=1}^n{MSE}_i$$

(5)

Thereby, the LOOCV allows the model to be trained on an almost entire dataset, resulting in the absence of a designated test set, and provides an unbiased estimate of the prediction error rate [60]. Furthermore, by removing the randomness in the hold-out splits, LOOCV always returns the same results.

Nevertheless, the model presented in this study aims to describe the physical properties and mechanical response of 2 groups of mixtures characterized by different filler materials. For this reason, a cross-validation was implemented by excluding 1 observation for each type of filler used, defining 12 stratified couples of observations. In other words, the model was fitted on 22 observations and tested on an independent data couple, repeating the procedure 12 times. In this way it was possible to reduce the computational cost of fitting 24 models ($n=24$ is the size of available dataset) and to decrease the variance of the cross-validation procedure because the overlap between the 12 training sets is lower, and consequently, the model outputs are somewhat less correlated with each other.

Among the possible combinations of $N$ and $FUN$, a grid search approach was performed to select the combination that maximized the difference between the effective number of parameters (i.e., $\gamma =\{0,\cdots ,\delta \}$ is a measure of how many weights and biases in the neural network are effectively used in reducing the loss function [55] divided by the total number of network weights and biases), $ENP$, and the averaged test error, ${MSE}_{CV}$, of the twelve fitted models:

$$\begin{array}{cc}\hfill L_{CV}\left({\delta }_i,{\alpha }_i,\widehat{\mathit{y}}\left({\mathit{W}}_{\mathit{i}}\right),\mathit{y},{\mathit{W}}_{\mathit{i}}\right)& =ENP-{MSE}_{CV}=\hfill \\ \hfill ={\displaystyle \sum _{i=1}^{12}}& \frac{{\delta }_i-2{\alpha }_i{\mathrm{t}\mathrm{r}\left({\nabla }^{2}F\left({\mathit{W}}_{\mathit{i}}\right)\right)}^{-1}}{{\delta }_i}-{\displaystyle \sum _{i=1}^{12}}{‖{\widehat{\mathit{y}}}_{\mathit{i}}\left({\mathit{W}}_{\mathit{i}}\right)-{\mathit{y}}_{\mathit{i}}‖}_2^2\hfill \end{array}$$

(6)

where ${\mathit{W}}_{\mathit{i}}$ is the matrix of weights and biases of the i-th SNN model trained in the leave-one-out cross-validation procedure and ${\delta }_i$ is the total number of parameters in the same i-th SNN. The trace of Hessian matrix can be computed by making use of the Jacobian matrix of the training set errors and the regularization parameters ${\alpha }_i$ and ${\beta }_i$, as described in Hagan et al. [55]. This approach in selecting the proper SNN model for the given fitting problem becomes important when a regularization technique is applied in order to reduce the computational cost: in fact, as the complexity of the network increases (i.e., in this study, as the number of hidden neurons increases), the effective number of parameters used to reduce the error function can be much less than the total number of weights and errors in the neural network, given that many of these parameters are forced to be very small or even 0 during the training process. In other terms, if a large network is selected, it will not necessarily overfit the data (thanks to the regularization process) but it will certainly require more computational time to calculate the network response which may not be too different from that of a smaller network. Therefore, the use of Equation (6) in the grid search process allows the SNN to be less complex but more efficiently identified. In fact, such a network is characterized by the widest gap between $ENP$ and ${MSE}_{CV}$:

$${\left\{N,FUN\right\}}_{best}=\underset{\left(N,\mathit{FUN}\right) \in \left(X_N,X_{\mathit{FUN}}\right)}{\mathrm{arg\; max}}{L}_{CV}\left(N,FUN\right)$$

(7)

3.4. Data Augmentation

The number of experimental observations available strongly affects the performance of a machine learning model: if this number is very small, as often happens when dealing with asphalt mixtures, it is not always possible to develop a neural model able to represent the physical complexity of the problem. However, there is the possibility to increase the diversity of a training set and improve the performance of a machine learning model by augmenting the data already available. This practice offers the opportunity to generate training sample variations without changing the semantics of the original data [60] when the collection of additional data is difficult due to time or cost limitations. The main challenge of data augmentation is to avoid generating any data that distorts the meaning of the original ones. Unlike computer vision applications in which it is possible to rotate, reverse, or flip the image data to augment the training set, in the case of input-output fitting problems there are no generalized data augmentation strategies, and it is easier to alter the original data without being able to identify the distortion.

Nevertheless, the recent study by Oh et al. has shown that cubic interpolation can be a numerical data augmentation method characterized by robustness against variations of data trend information and ease of use [61]. For these reasons and given the difficulties in preparing a very large amount of Marshall specimens with different bitumen and filler contents (both for cost and testing time), the training sets from the leave-one-out cross-validation procedure were augmented in the current study setup by means of the modified Akima cubic Hermite interpolation, also called the “makima” interpolation, which represents a MATLAB^® modification of Akima’s formula (1974) to avoid excessive local undulations (e.g., the makima interpolating surface is characterized by fewer undulations compared to cubic splines). The makima algorithm for a 2-D grid interpolation produces a piecewise function of cubic polynomials with continuous first-order derivatives using the values of the experimental grid points that are nearer (in both directions) to each query point [58,62]. In particular, the 2-D query points that identify the interpolated data added to each LOOCV training set are shown in Figure 5 (red cross marker) along with the experimental grid points (black plus sign marker). It is worth pointing out that the choice of the interpolation method and the query points is subjective, and these features could be considered a kind of hyperparameter. However, bearing in mind that the response of asphalt mixtures is typically non-linear and, at the same time, that variations in the response can be described by smooth curves, the makima interpolation represents a good choice because it allows overshoot of the most general cubic spline to be avoided and the physical non-linearity of the problem with respect to the linear interpolation to be considered. Instead, the selected configuration of the query points (Figure 5) was identified in trying to match the amount of information added to the size of the experimental dataset.

4. Results and Discussion

Figure 6 displays the volumetric properties of the asphalt mixtures containing SF and OPC. As expected, increasing the filler percentage resulted in lower air voids for both SF and OPC mixtures at each bitumen content. The use of SF as a replacement for OPC slightly decreased air voids in the asphalt mixtures. This could be attributed to the finer particle size of SF compared to OPC (Table 4), making it more effective in reducing air voids for each combination of filler and bitumen content. Air voids in both OPC and SF mixtures remained within the acceptable range of 3–5%, as per the Indian standard.

The VMA of a compacted asphalt mixture is determined by calculating the volume of the intergranular voids between aggregate particles as well as the air voids and the volume of bitumen that is not absorbed by the aggregates. Normally, an inverse bell-shaped curve should be observed for VMA, but the results obtained in this study (Figure 6) only show an increasing trend. This could be due to the fact that the decreasing trend is typically associated with lower bitumen content values than those examined in this study. A similar trend was observed in a study by Tenza-Abril et al. at a binder percentage greater than 4.5% [63].

According to Figure 6, the VMA values for OPC mixtures were slightly higher than those for SF mixtures, with a difference ranging from 0.18% to 1.20%. The VFB, which represents the percentage of VMA filled with bitumen, was slightly higher in the SF mixtures than in the OPC mixtures due to the previously determined air voids and VMA values (Figure 6). However, the differences were small, ranging from 0.18% to 1.19%. All the asphalt mixtures investigated in this study had VFB values within the 65–85% range prescribed by Indian regulations [34].

Experimental results for MS, MQ, and MMSR are shown in Figure 7, and all three Marshall parameters exhibited a bell-shaped curve for all mixtures studied, consistent with relevant literature for both SF and OPC mixes [64]. Increasing filler content resulted in higher values for all three Marshall parameters regardless of the bitumen content, in line with previous studies [59,65,66].

The highest MS values for both SF and OPC mixes were achieved at 5% bitumen and 8% filler, with the SF mix showing a 2.90% higher Marshall stability than the corresponding OPC mix. This suggests that SF can be a viable replacement for OPC in asphalt mixes, which can help reduce the disposal problem of SF waste. AMs with 4% filler and 5% bitumen content meet the MoRTH requirement for conventional bitumen mixes (MS > 10 kN) [34], while filler content raised to 6% or 8% satisfies the acceptance requisite for modified bitumen mixes (MS > 13 kN), highlighting the significant improvement in mechanical response provided by the filler.

Table S1 (Supplementary Materials) indicates that the Marshall flow of all the asphalt mixes examined falls within the MoRTH-specified range of 2 mm to 4 mm, with the lowest flow values being observed at 6% and 8% of SF and OPC. The combination of higher stability and lower flow serves as an empirical indicator of the resistance of asphalt mixtures against permanent deformations at high temperatures, such as the rutting phenomenon. Figure 7 reveals that the Marshall quotient (MQ) values of SF mixes are higher than those of OPC mixtures, ranging from 2.73% to 8.55%, depending on the filler content. Nevertheless, for all the mixtures investigated, MQ results fall within the acceptance range (2–6 kN/mm) specified by the Indian regulation [34]. Therefore, the addition of silica fume industrial wastes is unlikely to produce a highly plastic mechanical response.

Figure 7 shows the mean Marshall stability ratio results, which are similar for both groups of mixtures (SF and OPC) but SF mixtures exhibit slightly higher long-term aging resistance. This improvement ranges from 0.82% to 1.27% and depends on the bitumen and filler contents.

ITS values of SF and OPC mixes are shown in Figure 8. A reasonable bell-shaped curve was obtained for both types of mixes [63], with peak values at 5% of bitumen and a decreasing trend for higher bituminous binder percentages. A significant increment in the indirect tensile strength (up to 34.85%) can be observed using 6% and especially 8% of filler instead of 4% depending on the bitumen content. Comparing OPC and SF mixes for the same percentage of bitumen and filler, slightly higher values can be observed using SF filler with differences within the range 0.06% to 1.64%. Both SF and OPC mixes at 8% filler content have overcome the value of 1 MPa of ITS; therefore, a much higher value than the minimum acceptance requisite, equal to 0.6 MPa, is suggested in the literature [34].

Resistance to water damage is a fundamental property for asphalt mixes in order to ensure an adequate pavement service life. All the mixes investigated have presented a water sensitivity value as expressed by the loss of indirect tensile strength after water immersion equal to 20%, lower than the acceptance requisite prescribed by Indian regulations (Figure 8). Overall, the experimental results outline a better water resistance of the asphalt mixes at 6% and 8% filler contents compared to 4%. Asphalt mixes with SF have presented a slightly better water resistance according to the active and passive adhesion results; however, the differences in OPC asphalts are within the range 0.18% to 1.17% depending on the bitumen percentage.

The best-performing mixture prepared with 5% bitumen and 8% SF was characterized by Vv, VMA, MS, MMSR, MQ, ITS, and WS equal to 3.64%, 15.49%, 13.86 kN, 63.17%, 5.25 kN/mm, 1091.47 kPa, and 10.05%, respectively.

4.1. Machine Learning Results Discussion

As far as data modeling is concerned, the implementation of the LOOCV resampling method into the SNNs training and testing processes resulted in a smooth and accurate regression of the few available observations. In particular, the grid search allowed the optimal model structure to be identified (Figure 9): the SNN, characterized by $N$ = 21 neurons in its hidden layer and $FUN$ = TanH activation function, showed a mean squared error ${MSE}_{CV}$ equal to 0.0032 and a correlation coefficient $R_{CV}$ of 0.9988. This network has a total number of parameters $\delta$ equal to 238, of which 96.4% (=$ENP$) contribute to effectively reduce the loss function. Such results led to the highest $L_{CV}$ value of 0.9609. In the following, this model is referred to as TanH-SNN. The TanH activation unit, which gives the non-linearity required by the problem to the proposed model, has shown greater effectiveness compared to the more recently introduced ELU function (widely used in the computer vision field) for the greater number of parameters (i.e., $ENP$ equal to 0.9641 for TanH, 0.9319 for ELU) that actually contribute to defining the predictive performance of the SNN model, although the higher number of neurons that comprise the TanH-SNN hidden layer increases the calculation time compared to the ELU-SNN (i.e., $N$ is 21 for TanH, 16 for LogS, as shown in Figure 9). No relevant differences were found between the best-performing models with TanH and LogS activation function (Figure 9). All grid search results (for each number of neurons and activation function considered) are provided in the Supplementary Materials (Table S2).

To verify that the TanH-SNN had returned a smooth regression of the data (without excessive oscillations, not physically acceptable), the seven target variables were predicted, separately for the two filler materials, in increments of 0.05% bitumen content from 4.5 to 6.0% and of 0.10% filler content from 4.0 to 8.0%. In this way, the regression surface of the TanH-SNN model for each target variable was plotted in a space whose axes represent the bitumen content (x-axis), the percentage of filler by volume of the mix (y-axis), and the desired target variable (z-axis). In addition, a contour map was included in the plane defined by the bitumen and filler contents to outline the course of the regression surface. In Figure 10, these surfaces are shown for the Marshall stability parameter before (MS) and after aging (MSA). Regression surfaces of the remaining variables are included in the Supplementary Materials (Figures S1 and S2). The azimuth and elevation angles of the line of sight were changed between the surface graphs of the different target variables for a better understanding of the trends. The blue spheres in Figure 10 represent the experimentally observed data. The prediction surfaces in Figure 10 show an almost perfect adaptation of the model to the experimental data and demonstrate that TanH-SNN has correctly identified the mechanical (and volumetric) response of the mixes: in fact, all the findings reported in Section 4.1 can be verified anew by analyzing these surfaces. As a result, it can be concluded that the modeling procedure followed had a positive outcome on the fitting of the experimental data presented.

The neural model developed can enhance the asphalt mixture optimization by means of the numerical prediction of the mix response for different compositions without the necessity to perform any other laboratory tests.

Finally, to demonstrate that the chosen approach of data augmentation has allowed the numerical modeling to preserve the meaning of the experimental data, Figure 11 shows the contour maps of the percentage variation ($PV$) of the seven target variables between the asphalt mixes with SF and OPC, for the TanH-SNN model, and the makima interpolation of experimental results. By way of example:

$${PV}_{MS}(\%)=\frac{{MS}_{SF}-{MS}_{OPC}}{{MS}_{OPC}}·100$$

(8)

In particular, the two plots are quite similar from a qualitative and quantitative point of view. The TanH-SNN not only succeeds in reproducing the same percentage variation in the experimental grid points (black cross markers in Figure 11) but allows the empirical trend of the $PV$ variable in the domain of interest to be properly represented. Moreover, the results of these maps further confirm that silica fume leads to slightly better mechanical performance when used as filler compared to ordinary Portland cement.

5. Conclusions

This paper presents an integrated experimental and machine learning approach regarding the effective reuse of industrial waste silica fume as a sustainable mineral filler in asphalt mixtures for base courses of road pavements. Comparative laboratory analysis by varying the SF, OPC, and bitumen contents was carried out. The following conclusions can be drawn based on volumetric, durability, and strength tests along with the numerical analysis:

• Asphalt mixtures made with silica fume fulfil all the acceptance requisites prescribed by Indian regulations in terms of volumetric properties, Marshall stability, Marshall quotient, indirect tensile strength, and water sensitivity. In particular, in order to satisfy all acceptance requirements, it is recommended to prepare mixtures with an SF content equal to or higher than 6% and a bitumen content equal to or higher than 5%.
• The analysis carried out to investigate the aging effects demonstrated a good mechanical resistance of the SF asphalt mixes, with MMSR values up to approximately 70% achieved at 5.0% bitumen and 4.0% SF. The use of SF instead of OPC improved the mechanical behavior in terms of MMSR from 0.82 to 1.27% depending on bitumen and filler contents.
• The comparison between mixes with SF and OPC outlines an overall physical–mechanical equivalence, demonstrating the technical feasibility to substitute cement with silica fume without any substantial worsening from a mechanical point of view. Differences in terms of VMA and VFB ranged from 0.18% to approximately 1.20%. In terms of the highest MS, the SF mixture achieved a 2.90% higher value than OPC mixtures, with a slightly higher long-term aging resistance. With respect to ITS, SF mixtures achieved higher values than OPC mixtures, with differences ranging from 0.06 to 1.64%. Finally, differences in terms of resistance to water damage ranged between 0.18 and 1.17% depending on bitumen percentage.
• With respect to the recommended bitumen content (5%), the mix with 6% SF achieves an MS of 13.61 kN which is close to the MS of 13.63 kN obtained for the mix prepared with 8% OPC. Even if SF is characterized by a higher price (average price: 145 $/ton [67,68]) than OPC (average price 105 $/ton [69]), the cost of the mixtures prepared with SF is lower (roughly $19 of SF for cubic meter of mix) than that of OPC mixtures (roughly $24 of OPC for cubic meter of mix). In fact, the cementitious filler-based mix requires a higher amount of cement which, moreover, is characterized by a higher specific gravity.
• Artificial neural networks have allowed accurate modeling of the experimental data related to the main volumetric and mechanical parameters considered in the study as a function of bitumen, SF, and OPC contents.
• For the future development of the modeling procedure, the adoption of an optimization technique for the selection of the ANN structure to reduce computational time is recommended, such as Bayesian statistics.