A Deep Neural Network Approach towards Performance Prediction of Bituminous Mixtures Produced Using Secondary Raw Materials

Abstract

With the progressive reduction in virgin material availability and the growing global concern for sustainability, civil engineering researchers worldwide are shifting their attention toward exploring alternative and mechanically sound technological solutions. The feasibility of preparing both cold and hot asphalt mixtures (AMs) for road pavement binder layers with construction and demolition wastes (C&DWs) and reclaimed asphalt pavement (RAP) partially replacing virgin materials like limestone aggregates and filler has already been proven. The technical suitability and compliance with technical specifications for road paving materials involved the evaluation of mechanical and volumetric aspects by means of indirect tensile strength tests and saturated surface dry voids, respectively. Thus, the main goal of the present study is to train, validate, and test selected machine learning algorithms based on data obtained from the previous experimental campaign with the aim of predicting the volumetric properties and the mechanical performance of the investigated mixtures. A comparison between the predictions made by ridge and lasso regression techniques and both shallow (SNN) and deep neural network (DNN) models showed that the latter achieved better predictive capabilities, highlighted by fully satisfactory performance metrics. DNN performance can be summarized by R² values equal to 0.8990 in terms of saturated surface dry void predictions, as well as 0.9954 in terms of indirect tensile strength predictions. Predicted observations can be thus implemented within the traditional mix design software. This would reduce the need to carry out additional expensive and time-consuming experimental campaigns.

1. Introduction

In recent years, machine learning (ML) has emerged, radically changing several engineering areas and showing significant potential in process optimization and predictive analysis [1,2,3,4,5]. This has led to the growing interest of both the academic and industrial communities in the application of machine learning approaches for modeling both the physical and mechanical characteristics of construction materials [6,7,8]. Within the framework of asphalt mixtures for road pavement construction, several research studies have shown that these soft-computing techniques can effectively improve prediction accuracy and simultaneously decrease the costs associated with experimentations [9,10,11,12,13,14,15,16]. ML approaches can be integrated within conventional mix design, improving standard procedures by making them more efficient and less dependent on natural resources [17,18]. The best parameters for mechanical characterization of mixture behavior are experimentally determined, but changing even one variable requires costly new tests to be carried out. Therefore, a reliable predictive model could accelerate mixtures’ volumetric and mechanical characterizations, reducing the need for additional experiments and minimizing resource consumption and waste generation [19]. This is a particularly relevant aspect in an era when reducing environmental impact has become a global priority [20,21,22,23,24]. Within the context just outlined, the implementation of secondary raw materials represents an additional challenge, given their variability in terms of mechanical behavior [25,26,27,28,29,30,31,32]. However, machine learning supported by the latest advanced techniques can provide an innovative approach to address such variability, reducing the need for expensive and time-consuming experimental campaigns [33,34,35,36]. By way of example, Pattanaik et al. [37] developed a linear regression model known as ridge regressor (RR) to predict the abrasion loss of open-graded friction course mixtures prepared with EAF steel slag aggregates. Another predictive soft-computing technique based on linear regression is the lasso regressor (LR), an enhancement of the ridge regressor that implements a slight adjustment to the cost function to be minimized [38]. In the past few years, ML modeling based on artificial neural networks (ANNs) has also gained growing recognition among the scientific community. Many applications in the road engineering field can be found within the relevant literature due to their layered architecture and their ability to identify even highly nonlinear relationships between variables [39,40,41,42,43]. However, their complex structure requires sensitivity analyses to be performed after modeling in order to better understand the influence each variable had in achieving the final predictions. The primary objectives of the current research can be placed within the outlined framework, as they include the design of performance predictive models implementing several cutting-edge soft-computing techniques.

To this end, different linear regression techniques (LR and RR) and both shallow (SNN) and deep neural network (DNN) models were fine-tuned to process data obtained from the experimental campaign described in Rondinella et al. [44], thus allowing an ML-based procedure to be developed for predictive purposes. Both the volumetric and mechanical performance of the investigated mixtures were predicted by ML models trained on input information. These included gyratory revolutions, mixing technology, and the percentage contents of C&DWs, RAP, water, cement, emulsion, and total bitumen.

The structure of this paper can be outlined as follows: the Section 2 provides a brief summary of the experimental campaign and the mathematical framework for the machine learning algorithms analyzed; the Section 3 describes the results obtained from the predictive modeling procedures, providing accuracy and sensitivity information; the Section 4 outlines the most important conclusions, highlighting the applicability constraints of the developed predictive models and pointing out the possibilities for future developments.

2. Materials and Methods

2.1. Raw Materials

The asphalt mixtures investigated here are blends of coarse aggregates (either limestone virgin aggregates, aggregates from construction and demolition waste, or RAP), fillers (either virgin limestone filler or recycled filler from construction and demolition waste), and binders (either neat bitumen or modified bitumen for hot AM, bitumen emulsion and cement for cold AM). The virgin limestone aggregates and limestone filler were supplied by a crushing and sieving plant located in the Campania region (Italy). The main technical properties of the aggregates are shown in

Table 1. Main properties of the inert fractions used to blend asphalt mixtures.

Property	Standard	Unit	Limestone 10–20 mm	Limestone 6–12 mm	Limestone Sand	Limestone Filler	CD&W1 2/16 mm	CD&W2 Filler
Bulk specific gravity	EN 1097-6 [46]	g/cm3	2.685	2.686	2.689	2.737	2.925	2.934
Los Angeles value	EN 1097-2 [47]	%	20.6	20.1			36
Flakiness index	EN 933-3 [48]	-	8	11			13
Rigden voids	EN 1097-4 [49]	%				41.44		53.82
Sand equivalent	EN 933-8 [50]	%			95.3	92.0		81.1

, along with their test standards. The RAP was derived from the milling of the wearing course and binder layers of an existing distressed asphalt pavement located in the Campania region (Italy); the size designation according to UNI EN 13108-8 [45] and the main classification parameters are shown in

Table 2. Main properties of the reclaimed asphalt.

Property	Standard	Unit	RAP
Designation	EN 13108-8 [45]	-	16 RA 0/10
Binder content	EN 12697-1 [51]	%	4
Bulk specific gravity of the aggregates	EN 1097-6 [46]	g/cm3	2.52
Flakiness index	EN 933-3 [48]	-	10

.

The recycled coarse aggregates and filler from C&DW were supplied from a recycling plant; it separates the inert fractions from metals and other foreign matter, then crushes and sieves the remaining inert fractions to obtain new secondary construction materials. In detail, the coarse aggregates from C&DW (C&DW1) were designated as 2/16 mm size, while the filler from C&DW (C&DW2) was entirely passing through a sieve size equal to 1 mm. All the relevant technical properties of the secondary aggregates are shown in Table 1.

As concerns the binders, the main properties of the neat bitumen with a penetration class 50/70 and the 5% SBS-modified bitumen used to blend hot AM are shown in

Table 3. Main properties of the binders: (a) neat and modified bitumen, (b) bitumen emulsion, and (c) pozzolanic cement.

(a)
Property	Standard	Unit	Neat Bitumen	Modified Bitumen
Penetration at 25 °C	EN 1426 [52]	dmm	68	52
Softening point	EN 1427 [53]	°C	46	87
Dynamic viscosity at 135 °C	EN 13702 [54]	Pa s	0.25	0.77
(b)
Property	Standard	Unit	Bitumen Emulsion
Water content	EN 1428 [55]	%	40
pH value	EN 12850 [56]	-	4.2
Settling tendency at 7 days	EN 12847 [57]	%	5.8
Softening point after water evaporation	EN 1427 [53]	°C	49.5
(c)
Property	Standard	Unit	Pozzolanic Cement
Initial setting time	EN 196-3 [58]	min	112
Compressive strength at 2 days	EN 196-1 [59]	MPa	27.8
Compressive strength at 28 days	EN 196-1 [59]	MPa	61.2
Volume constancy	EN 196-3 [58]	min	0.52

a; Table 3b and Table 3c summarize, respectively, the main properties of the bitumen emulsion and Portland cement used to blend cold AM.

2.2. Asphalt Mixtures

An experimental campaign was carried out to identify the optimum binder content of 7 alternative AMs for the binder layer of a road pavement. Each mixture underwent gyratory compaction at the desired number of gyrations [60] using a variable binder content until the optimum mix composition was identified. An exhaustive description of both the carried-out experimental campaign and the prepared hot and cold AMs can be found in Rondinella et al. [44]. However, the abbreviations used to identify the individual mixtures are provided in

Table 4. Summary of the abbreviations used to identify the produced mixtures.

Mixture	Description
HMAmod	hot AM made up of 100% limestone aggregates and variable SBS polymer-modified bitumen content in the range [4.5%–5.5%]
HMAC&DW1	hot AM made up of 40% construction and demolition waste aggregates (C&DW1), 60% limestone aggregates, and a neat bitumen 50/70 content in the range [7.0%–7.5%]
HMAmodC&DW1	hot AM made up of 40% construction and demolition waste aggregates (C&DW1), 60% limestone aggregates, and an SBS polymer-modified bitumen content in the range [6.0%–7.0%]
CMA	conventional cold AM made up of 76% RAP and 24% limestone aggregates with 4% water, 0.5% cement, and variable bitumen emulsion content, in the range [3%–5%]
CMAC&DW1	cold AM with 30% RAP, 30% construction and demolition waste aggregates (C&DW1), and 40% limestone aggregates with 7% water, 5% cement, and variable bitumen emulsion content, in the range [3%–5%]
CMAC&DW2_1	cold AM with 76% RAP, 20% limestone aggregates, and 4% filler from construction and demolition waste aggregates (C&DW2) with 6% water, 2.5% cement, and variable bitumen emulsion content, in the range [3%–5%]
CMAC&DW2_2	cold AM with 30% RAP, 64% limestone aggregates, and 6% filler from construction and demolition waste aggregates (C&DW2) with 15% water, cement content in the range [6.5%–7.5%], and variable bitumen emulsion content, in the range [3%–5%]

since they will represent the categorical variables of the subsequent predictive models.

All the mix compositions underwent the experimental characterization of two main technical control parameters, complying with the main Italian technical specifications for road asphalt materials: the air void content (EN 12697-8 [61]), determined according to the saturated surface dry (SSD) bulk density (EN 12697-6 [62]) and the mathematical maximum density (EN 12697-5 [63]), and the indirect tensile strength (ITS) at 10 °C (EN 12697-23 [64]).

2.3. Database Statistics

The analyzed dataset consisted of 70 observations for each evaluated variable, namely percentage content of water (WC), cement (CC), emulsion bitumen (EBC), total bitumen (TBC, expressed as the cumulative amount of fresh, emulsified, and rejuvenated bitumen contained in RAP), construction and demolition wastes (C&DW1, C&DW2), and reclaimed asphalt pavement (RAP), along with gyratory revolutions (GRs), the SSD void content (SSDV), and the indirect tensile strength (ITS) at 10 °C. Finally, as previously explained, a categorical variable was used to identify the mixing technology adopted to prepare the corresponding mixtures.

To provide a brief summary of the processed variables, descriptive statistics referring to the analyzed dataset are shown in

Table 5. Descriptive statistics referring to the analyzed dataset.

Variable	Description	Average	Std	Min	25%	50%	75%	Max
GR [–]	Gyratory Compaction Revolutions	165.14	8.80	160.00	160.00	160.00	175.00	180.00
C&DW1 [%]	Aggregates from Construction and Demolition Waste	12.43	17.56	0.00	0.00	0.00	30.00	40.00
C&DW2 [%]	Filler from Construction and Demolition Waste	2.11	2.68	0.00	0.00	0.00	5.50	6.00
RAP [%]	Reclaimed Asphalt Pavement Content	32.20	29.85	0.00	0.00	30.00	76.00	76.00
WC [%]	Water Content	6.13	5.85	0.00	0.00	6.00	13.00	15.00
CC [%]	Cement Content	2.86	2.95	0.00	0.00	2.50	6.12	7.50
EBC [%]	Emulsion Bitumen Content	2.63	2.02	0.00	0.00	3.00	4.00	5.00
TBC [%]	Total Bitumen Content	6.57	0.78	4.50	6.21	6.81	7.31	7.50
SSDV [%]	Saturated Surface Dry Voids	6.48	3.42	2.06	3.94	5.52	9.39	13.35
ITS at 10 °C [kPa]	Indirect Tensile Strength at 10 °C	1237.25	1036.81	202.08	407.80	738.77	2260.18	3529.01
Categorical [–]	Hot and Cold Mixing Technologies	–	–	–	–	–	–	–

. Details about the 25th, 50th, and 75th percentiles, as well as minimum, maximum, average, and standard deviation values, are provided. The following machine learning operations will have as their main objective the development of a model that can reliably predict both SSDV and ITS at 10 °C values on the basis of the remaining investigated variables. This approach could allow for a comprehensive characterization of the investigated mixtures from both volumetric and mechanical perspectives, simplifying the integration of predicted values into traditional practices for road pavement mix design.

2.4. Ridge and Lasso Regression

Denoting by $y$ the target variable involved in model prediction and by $\{x_1,\dots ,x_n\}$ the set of independent variables, a linear regression is aimed at finding $\{{\beta }_0,\dots ,{\beta }_n\}$ parameters so that ${\beta }_0+{\sum }_{i=1}^n{\beta }_i{x}_i$ results as a reliable approximation of $y$. Since a dataset is composed of $m$ observations for each $x_i$, denoted by $x_{i,j}$, ${\beta }_i$ parameters are selected in order to identify the minimum of the residual sum of squares (RSS), expressed as (Equation (1)):

$$\mathrm{R}\mathrm{S}\mathrm{S}={{\sum }_{j=1}^m\left({\beta }_0+{\sum }_{i=1}^n{\beta }_i{x}_{i,j}-y_i\right)}^2$$

(1)

Hastie et al. [65] formally introduced the so-called forward stepwise selection method [66]. It consists of assuming ${\beta }_0$ is equal to the average of $y$, while ${\beta }_i$ are all set equal to 0. Subsequently, $i$ and the corresponding ${\beta }_i$ are iteratively selected in order to produce a model capable of reducing RSS to the lowest value. Once chosen, the value of ${\beta }_i$ remains unchanged for the remainder of the process. Consequently, each following model is differentiated from the previous one only by a single $\beta$ parameter. Finally, n + 1 different models are produced, each resulting in an iteratively lower RSS.

Ridge regression (RR) represents a particular linear regression model since, even before tuning $\beta$ parameters, a pre-processing normalization procedure is applied to predictor variables [67]. This procedure consists of evaluating the mean and the standard deviation of each predictor and then, for each of them, subtracting the estimated mean from its respective values and dividing by the estimated standard deviation. This is performed to ensure that the same scale is representative of all values and that the standard deviation is always equal to 1. Furthermore, a penalty to be applied to ${\beta }_i$ is introduced by RR models, making those with lower values of regression coefficients preferable. This penalty is aimed at lowering the values of the least important predictors while keeping the higher ${\beta }_i$ associated with the most important ones. The cost function (CF_RR) to be minimized during an RR is slightly modified by virtue of the introduced penalty, as follows (Equation (2)):

$${\mathrm{C}\mathrm{F}}_{\mathrm{R}\mathrm{R}}=\mathrm{R}\mathrm{S}\mathrm{S}+\alpha {\sum }_{j=1}^m{{\beta }_j}^2$$

(2)

with $\alpha$ representing the positive penalty parameter that determines the trade-off between both minimizing error and the size of regression coefficients. Different values of $\alpha$ result in the design of different models.

The only difference between RR and lasso regression (LR) [38] is represented by the slightly revised cost function (${\mathrm{C}\mathrm{F}}_{\mathrm{L}\mathrm{R}})$, expressed as follows (Equation (3)):

$${\mathrm{C}\mathrm{F}}_{\mathrm{L}\mathrm{R}}=\mathrm{R}\mathrm{S}\mathrm{S}+\alpha {\sum }_{j=1}^m\left|{\beta }_j\right|$$

(3)

The advantage of this adjustment is the LR capability to nullify some regression coefficients when the value of $\alpha$ is high, thus yielding smaller models. Therefore, it can be concluded that $\alpha$ represents a hyperparameter of both RR and LR models, so it will require to be optimized by means of suitable optimization procedures.

2.5. Artificial Neural Network Modeling

Artificial Neural Networks (ANNs) represent mathematical models whose architecture is similar to a biological nervous system. In this sense, it is composed of several interconnected artificial neurons, typically stratified and arranged in a sequence of distinct layers, known as the input layer, one or more hidden layers, and the output layer [68]. Each neuron is linked to every neuron of the other layers, but it is independent of the neurons of the same layer. These inter-neuronal connections are associated with certain weights, either positive or negative, that can be adjusted during the training procedure according to a specified learning algorithm and a predetermined loss function [69]. Each neuron of the input layer is associated with a single input variable, known as a predictor. Neurons belonging to the hidden layer process the information coming from the input one, and their number is a typical hyperparameter that needs to be optimized. The number of hidden layers determines whether the architecture should be called shallow (SNN, composed of a single hidden layer) or deep (DNN, composed of more than one hidden layer) and can be optimized as a different hyperparameter. Finally, the output layer is associated with the predicted target variables. Both SNNs and DNNs are capable of dealing with non-linear relationships since hidden layers can be equipped with a non-linear activation function. Conversely, output layers usually employ linear activation functions.

Typically, during the iterations of the training procedure, the neural model’s output is predicted according to Equation (4):

$$f\left(X\right)=W_2\left[f_A\left(W_{1}X\right)\right]$$

(4)

where $X$, $f_A$, $W_1$, and $W_2$ represent the predictors vector, the hidden-layer activation function, the weights matrix of input-hidden connections, and the weights matrix of hidden-output connections, respectively.

Weights matrices progressively evolve according to the implemented learning algorithm, and the most scientifically sound are described within the relevant literature [70,71,72]. During the past few years, ML models that are neural network-based have demonstrated their capabilities of reliably approximating strongly nonlinear relationships, leading to remarkable performance in several engineering-related fields, especially in road pavement design and construction [73,74].

2.6. Hyperparameter Optimization

To determine the portion of the dataset to be used to train the models and that to be used to test them, observations were randomly shuffled and then divided into training and testing subsets containing 55 and 15 observations out of 70, respectively. Therefore, roughly 80% of the available observations were used to properly train the models, whereas roughly 20% of them were used to evaluate their performance. The pre-processing normalization described in Section 2.3 was implemented in each model, both to enhance the predictive accuracy and to reduce the time needed to reach convergence. A k-fold cross-validation method was also implemented in each model to fairly evaluate the training and validation performance of each developed model, to fine-tune their hyperparameters, and also to avoid the occurrence of inconvenient overfitting phenomena. k-fold cross-validation consists of partitioning the training set into k distinct partitions. Among these, k − 1 serves as training set, while the remainder makes up the so-called validation set. These steps are iteratively performed k-times, thus generating k-validation scores. Averaging the k-scores obtained during this sort of pre-test practice will provide a more accurate estimation of the models’ actual predictive capabilities. k-value was chosen equal to 5, in accordance with the most scientific sound literature references [75]. Finally, the best-performing hyperparameters were identified by means of an exhaustive grid search performed for each soft-computing technique. Both with respect to RR and LR, the only parameter that could be optimized was the penalty $\alpha$, and it was searched in a range from 10^–4 to 10⁵ on a log scale. Regarding the ANN model, 5 different hyperparameters were optimized. They were represented by the number of hidden layers, the number of neurons composing each hidden layer, the activation function implemented within hidden layers, the solver algorithm, and the maximum number of training iterations searched within the ranges described in

Table 6. Summary of the grid search performed to find the best ML hyperparameters.

ML Model	Hyperparameter	Search Range	Selected Value
RR	Penalty parameter	[10−4–105] *	10−3 for SSDV 10−2 for ITS at 10 °C
LR	Penalty parameter	[10−4–105] *	10−4 for SSDV 10−4 for ITS at 10 °C
ANN	Hidden Layers	[1, 2, 3]	3
	Neurons for each hidden layer	[1–30]	24
	Hidden layer activation function	[Identity, Logistic, TanH, ReLU]	ReLU
	Solver	[SGD [70], Adam [76]]	Adam
	Maximum number of iterations	[1000, 5000]	5000

* Hyperparameters searched on a log scale.

, in accordance with literature references [38,69].

An overfitting detection procedure was introduced to further minimize the possibility of overfitting occurring during the training of the neural model. It consisted of observing each subsequent validation score and stopping the training process as a negligible improvement (set as 10⁻⁴) after a predetermined number of consecutive iterations was shown. This stopping criterion was implemented by setting the number of consecutive iterations equal to 20, according to well-established practices [77].

The whole methodology described in this section was fully developed by means of Python 3.9.12 software.

3. Results and Discussion

Machine Learning Modeling Results

Histogram plots displayed in Figure 1 are crucially important to establish the most effective soft-computing technique for achieving predictive modeling goals. One-to-one comparisons between test observations and their corresponding predictions can be observed with respect to both SSDV and ITS values, shown at the top and bottom of the figure, respectively. The IDs of the individual observations that make up the test vector are represented on the x-axis, whereas the corresponding predictions are represented on the y-axis by histograms filled with different colors. Specifically, dark blue, light blue, orange, and dark yellow histograms represent the target values and the predictions produced by DNN, RR, and LR models, respectively. It is worth noting that, based on both SSDV and ITS comparisons, there is an observable consistency between ML models’ predictions and target values. Although the DNN model exhibited a comparatively higher accuracy both in terms of SSDV and ITS predictions, the designed ML models allowed both volumetric and mechanical parameters to be accurately predicted. In this way, a comprehensive characterization of mixture performance could be reliably obtained, thus achieving an important outcome from the engineering perspective. Further considerations about the accuracy of the developed predictive models can be made by considering the results obtained in terms of the six performance metrics implemented (

Table 7. Description of the evaluated goodness-of-fitness metrics.

ML Model	Performance Metric for SSDV
	MAE [%]	MAPE [%]	MSE [%2]	RMSE [%]	R	R2
RR	0.97	24.94	1.53	1.24	0.9540	0.8464
LR	0.98	25.11	1.53	1.24	0.9527	0.8470
DNN	0.75	19.17	1.01	1.00	0.9645	0.8991
	Performance Metric for ITS at 10 °C
	MAE [kPa]	MAPE [%]	MSE [kPa2]	RMSE [kPa]	R	R2
RR	237.48	15.13	118,900.53	344.82	0.9520	0.9063
LR	237.85	15.19	119,220.83	345.28	0.9519	0.9060
DNN	52.22	4.98	5779.01	76.02	0.9988	0.9954

). Specifically, four error metrics and two correlation metrics were analyzed, namely mean absolute error (MAE), mean absolute percentage error (MAPE), mean squared error (MSE), root mean squared error (RMSE), Pearson correlation coefficient (R), and the coefficient of determination (R²). Their mathematical formulations, along with their meaning, were fully described in Rondinella et al. [78].

Remarkable performance was achieved by all developed ML models, both in terms of SSDV and ITS at 10 °C. With respect to the former, MAE ranged from a minimum value of 0.75% (achieved by the best-performing DNN) to a maximum value of 0.98% (achieved by LR). Similarly, the best results in terms of the coefficient of determination R² were achieved by DNN, with a value roughly equal to 0.90, significantly higher than the RR-R², which was roughly equal to 0.85.

The differences in the developed models’ accuracies became even more pronounced in terms of the latter. With respect to ITS at 10 °C, DNN outperformed both LR and RR, returning an order of magnitude lower MAE (roughly equal to 52 kPa). DNN performance was significantly better also in terms of R², with a value of 0.99, which is notably higher than the 0.91 achieved by both lasso and ridge regression. It is worth noting that the DNN model allowed simultaneous predictions of both SSDV and ITS to be made. This option was not provided by both linear regressors, thus requiring training, validating, and testing a single model for each predictive variable involved.

Regression plots were also diagrammed in order to provide additional evidence of the developed models’ capabilities (Figure 2). The one on the left refers to SSDV values, while the one on the right refers to ITS ones. Test vector observations are represented along the x-axis, while along the y-axis, DNN, RR, and LR predictions are shown. Finally, a dark blue solid line serves as the benchmark of perfect accuracy, representing an ideal scenario where target observations and predictions perfectly match. Predictions were represented by diamond markers, and their colors were chosen in alignment with the histogram plot of Figure 1. Therefore, light blue crosses identify values predicted by the DNN model, while orange and dark yellow crosses identify RR and LR predictions, respectively. The proximity of the cross markers to the dark blue solid line serves as an indicator of the higher accuracy of the corresponding predictive algorithm. In accordance with the results described in Table 7, the findings showed that the DNN model was characterized by the best consistency in terms of predictive performance. Light blue crosses were mostly close to the perfect accuracy line, thus supporting the remarkable results achieved in terms of low error metrics and high correlation metrics. Considering the results obtained by the previous predictive modeling carried out by Rondinella et al. [44], it can be found that the present DNN model performed better than the formerly developed CatBoost model. DNN predictions were characterized by slightly higher Pearson and determination coefficients, both in terms of SSDV and ITS at 10 °C, highlighting the neural model’s improved reliability. With respect to MSE and RMSE metrics, the neural model developed in this paper showed improvements of about 18% and 10% in terms of SSDV predictions and about 8% and 4% in terms of ITS predictions, compared with the same metrics obtained from the previous CatBoost model.

Since DNN outperformed all its competitors, a SHapley Additive explanation (SHAP) analysis [79] was carried out to make the neural model predictions more interpretable. Such a technique was introduced in 2016, and its goal was the estimation of the contribution of each player in a collaborative game outcome [80]. In this sense, the basic concept behind SHAP analysis is the allocation of a score to each input feature in order to quantify the contribution of that specific feature to the single model prediction. These scores are known as Shapley values (${\varphi }_i$), and they can be evaluated according to Equation (5):

$${\varphi }_i(f,x^{\prime })=\sum _{w^{\prime }\subseteq \{x_1^{\prime },\dots ,x_n^{\prime }\}\backslash \left\{x_i^{\prime }\right\}}{\displaystyle \frac{\left|w^{\prime }\right|!\left(N-\left|w^{\prime }\right|-1\right)!}{N!}}\left[f(w^{\prime }\bigcup x_i^{\prime })-f\left(w^{\prime }\right)\right]$$

(5)

with $w^{\prime }$, $x^{\prime }$, $N$, and $f\left(w^{\prime }\right)$ representing the model’s input features, the feature that needs to be explained, the feature number, and the model prediction values for features $w^{\prime }$, respectively. In this study, two different SHAP analyses were carried out to extract explanations about the developed DNN model, and two summary plots were diagrammed to summarize the obtained information (Figure 3). After assigning the output, features were ranked in a descending importance order from top to bottom, determining the so-called feature importance. It can be observed that in terms of SSDV predictions, RAP, TBC, and GR were the three most impactful input variables. On the other hand, in terms of ITS predictions, C&DW2, EBC, and GR showed the greatest importance. Also, the feature effect was determined in addition to the feature importance. Proceeding from left to right, a gradual shift in colors from blue to red suggests a relationship of direct proportionality between the target variable and the selected feature. Conversely, an inverse transition from red to blue stands for a relationship of inverse proportionality between the couple of variables under consideration. The former could be observed for pairs (SSDV–RAP) and (ITS–GR), while the latter could be observed for (SSDV–TBC), (SSDV–GR), (ITS–C&DW2), and (ITS–EBC).

4. Conclusions

In evaluating the feasibility of preparing sustainable AMs using recycled or waste materials, it is crucial to ensure that these innovative mixtures align with the thresholds prescribed by standards requirements, especially in terms of volumetric and mechanical performance. To this end, the present study first evaluated the saturated surface dry voids (SSDVs) and the indirect tensile strength (ITS) at 10 °C of different technological solutions prepared for road pavement binder layers, both with hot and cold mixing methodologies. These mixtures partially replaced virgin aggregate and filler materials like limestone with different sizes of construction and demolition waste (C&DW) as well as reclaimed asphalt pavement (RAP) to match the ninth sustainable development goal by reducing natural resource consumption. Then, several cutting-edge soft-computing techniques, namely ridge and lasso regressions (RRs and LRs, respectively), as well as shallow (SNNs) and deep neural networks (DNNs), were investigated in order to design a reliable machine learning (ML) methodology for predicting both the SSDV and ITS of the investigated sustainable AMs. The accuracy and reliability of the developed ML models were fully characterized through the evaluation of six different performance metrics. On the basis of modeling results, the below considerations were concluded:

• Considerable improvements in ML algorithm predictive capabilities were achieved by means of a comprehensive grid search that allowed optimal hyperparameters to be effectively identified.
• All the soft-computing techniques were trained, validated, and tested using the same data so that the achieved performance could be fairly compared. All the models allowed predictions to be made in terms of SSDV and ITS on the basis of compositional variables, gyratory revolutions, and a categorical variable that distinguished the technology used to mix AMs. In terms of SSDV, the DNN outperformed linear regression models, showing MAE and R² values equal to 0.75% and 0.8991, respectively. These results were comparatively higher with respect to the performance achieved by the simpler regressors, whose best results in terms of MAE and R² were equal to 0.97% (achieved by RR) and 0.8470 (achieved by LR), respectively.
• The DNN also showed outstanding performance with respect to indirect tensile strength predictions. The MAE value (52.22 kPa) was an order of magnitude lower with respect to the results achieved by linear regression models (roughly 237 kPa). Also, the R² value (0.9954) achieved by the DNN was significantly better than the determination coefficients achieved by both LR and RR, whose best value reached a maximum of 0.9063.
• The outlined DNN model also performed slightly better than the former CatBoost model, previously developed by Rondinella et al. [44] based on the same experimental campaign. Comparing the results obtained by the best current predictive model and the best previous one, the MSE and RMSE metrics of the DNN model developed in this manuscript showed improvements roughly equal to 18% and 10% in terms of SSDV predictions and roughly equal to 8% and 4% in terms of ITS predictions with respect to the same metrics achieved by the former CatBoost model.
• Focusing on the predictions made by the DNN model and splitting the sensitivity analysis for each predicted variable, the results obtained from the SHAP analysis showed that in terms of SSDV, the most impactful variables were the percentage contents of RAP and total bitumen as well as the gyratory revolutions, the former demonstrating a direct proportionality with the target variable, while the latter and the third an inverse one. On the other hand, in terms of ITS, the most impactful variables were the percentage contents of C&DW2 and emulsion bitumen, as well as the gyratory revolutions; the first two variables also demonstrated an inverse proportionality with the target variable, while the third demonstrated a direct one.

Machine learning modeling results proved to be encouraging. Some aspects were not considered since they were out of the scope of the present research. However, these will certainly be the target of future developments and will cover: (i) a deeper analysis of AMs containing increased proportions of secondary raw materials like C&DW and RAP in order to further align with the sustainable development goals and circular economy standards; (ii) the mechanical characterization of AMs by investigating different parameters such as fatigue or permanent deformation resistance; and (iii) the design of advanced ML methodologies implementing innovative and potentially more robust soft-computing techniques.