Compact and Interpretable Strength--Ductility Prediction of Heat-Treatable Aluminum Alloys Using Physics-Retained Explainable Multi-Target Machine Learning

1. Introduction

Aluminum alloys are widely used in aerospace, automotive, marine, and advanced structural applications owing to their high specific strength, good corrosion resistance, formability, and suitability for property control through heat treatment Dillibabu et al. (2023); Wahid et al. (2020); Rambabu et al. (2017). Among different aluminum alloy systems, heat-treatable alloys are particularly important because their mechanical response can be tailored through the combined effects of alloy chemistry, solution treatment, quenching, and artificial aging MacKenzie (2016); Siddesh Kumar et al. (2022). The final balance between yield strength (YS), ultimate tensile strength (UTS), and elongation (El) is governed by complex interactions among alloying elements, precipitation behavior, solute distribution, and heat-treatment parameters Girgis et al. (2023). Therefore, identifying suitable composition–process combinations that provide high strength without severe loss of ductility remains a challenging problem in aluminum alloy design.

Traditional alloy development relies heavily on experimental trial and error, thermodynamic calculations, and expert metallurgical knowledge. Although these approaches are valuable, they are often time-consuming and costly, particularly when the design space contains several alloying elements and processing variables. Even small variations in elements such as Cu, Mg, Zn, Si, Zr, Ti, and Li, or changes in solution treatment temperature, aging temperature, and aging time, can significantly modify the resulting microstructure and mechanical properties Wang et al. (2007). In this context, machine learning (ML) has emerged as a useful tool for accelerating materials design by learning composition–process–property relationships from experimental datasets Vasudevan et al. (2021).

Recent studies have demonstrated the potential of ML models, including random forests, support vector regression, gradient boosting, and neural networks, to predict the mechanical properties of aluminum alloys Choi et al. (2026); Yang et al. (2024); Pichlmann et al. (2025). However, high prediction accuracy alone is insufficient for practical alloy design. A useful data-driven framework should also identify compact and physically meaningful descriptors, provide interpretable feature-property relationships, and quantify predictive uncertainty. This is particularly important for literature-curated materials datasets, where variability may arise from differences in experimental procedures, heat-treatment conditions, sample preparation, and testing protocols Pichlmann et al. (2025); Lookman et al. (2017). Therefore, an explainable and reliability-aware ML framework is required for robust strength–ductility prediction and candidate screening in heat-treatable aluminum alloys.

1.1. Research Gap and Motivation

Although several ML-based studies have been reported for predicting the properties of aluminum and other metal alloys, important methodological gaps remain. Some existing studies use all available descriptors without an explicit compact feature-selection stage, whereas others apply correlation-, importance-, or optimization-based selection after constructing statistical, elemental, or physics-derived descriptors Choi et al. (2026); Yang et al. (2024); Pichlmann et al. (2025); Shah et al. (2025). Such approaches can identify predictively useful variables, but they do not necessarily guarantee that processing variables known a priori to govern the material response will remain available to the final model. For heat-treatable aluminum alloys, solution treatment temperature, aging temperature, and aging time influence solute dissolution, precipitation kinetics, and the transition among under-aged, peak-aged, and over-aged conditions Zhang et al. (2024). The specific gap addressed here is therefore not the general absence of physics or feature selection in previous studies, but the absence of an explicit domain-constrained selection rule that preserves these prespecified heat-treatment descriptors throughout subset construction while allowing the remaining composition descriptors to be ranked from the training data.

In addition, strength and ductility are frequently predicted using independent target-wise models, even though YS, UTS, and El are interrelated mechanical responses. In alloy design, increasing strength often reduces ductility; therefore, the simultaneous prediction of these properties is important for evaluating the strength–ductility balance Cecchel and Cornacchia (2025). In addition, many studies focus mainly on point prediction and provide limited discussion on uncertainty-aware candidate recommendation. Without uncertainty estimation, model predictions may appear overly confident, particularly for candidate compositions or heat-treatment conditions that differ from those well represented in the training data Zhang et al. (2022).

1.2. Contributions of This Work

To address these limitations, this study proposes a physics-retained, explainable, and uncertainty-aware multi-target machine-learning framework for predicting the strength–ductility response of heat-treatable aluminum alloys. Unlike unconstrained filter, embedded, or optimization-based feature-selection procedures, the proposed strategy structurally retains solution treatment temperature, aging temperature, and aging time in every evaluated subset, while the remaining composition descriptors are ranked using a normalized combination of target-averaged statistical relevance and mean embedded nonlinear importance. The term “physics-retained” is used in this specific sense: the framework does not embed governing physical equations into the learning algorithm, but uses established metallurgical knowledge to constrain the admissible feature subsets. Its scientific contribution lies in integrating and evaluating this domain-constrained selection strategy with simultaneous prediction of YS, UTS, and El, a matched full-descriptor comparison, target-wise model interpretation, and ensemble-conformal uncertainty assessment. Selected literature-based modeling strategies were implemented as adapted in-house baselines using the present dataset and computational protocol; they are not presented as direct reproductions because the original studies used different datasets, alloy systems, and experimental settings.

The main contributions of this work are summarized as follows:

• A domain-constrained physics-retained hybrid feature-selection strategy was developed in which the three prespecified heat-treatment descriptors were mandatory in every evaluated subset, while the remaining composition descriptors were ranked using target-averaged statistical relevance and mean embedded nonlinear importance. The resulting selected pipeline used 12 rather than 17 descriptors.
• A multi-target machine-learning framework was developed for simultaneous prediction of YS, UTS, and El, enabling strength–ductility assessment within a unified modeling pipeline.
• The 12-feature physics-retained pipeline was compared with a matched 17-feature pipeline across the same 25 outer evaluations, with variability, paired differences, corrected confidence intervals, and statistical significance reported.
• SHAP-based explainability was employed to characterize target-specific model attributions for alloying elements and heat-treatment variables, with correlated-feature and causal-interpretation limitations explicitly acknowledged.
• An ensemble-conformal uncertainty estimation strategy was incorporated to generate sample-specific target-wise marginal prediction intervals and support reliability-aware candidate screening.

The remainder of this paper is organized as follows. Section 2 reviews related works on ML-based aluminum alloy property prediction, feature selection, explainability, and uncertainty estimation. Section 3 describes the dataset, benchmark models, proposed physics-retained feature-selection strategy, model development, SHAP interpretation, and uncertainty-aware candidate screening procedure. Section 4 discusses the experimental results, benchmark comparison, selected features, target-wise prediction performance, SHAP analysis, uncertainty intervals, and candidate screening. Finally, Section 5 concludes the study and outlines future research directions.

2. Related Work

Machine learning has become an important tool in materials informatics for accelerating alloy-property prediction, composition–process optimization, and candidate screening. In aluminum alloys, the mechanical response is governed by coupled interactions among alloy chemistry, processing route, heat-treatment condition, and microstructure. Therefore, recent studies have increasingly used data-driven models to predict yield strength (YS), ultimate tensile strength (UTS), elongation (El), hardness, and other performance indicators from composition, processing, microstructural, or physics-based descriptors. The related literature can be broadly grouped into four categories: ML-based prediction of aluminum alloy properties, feature engineering and feature selection, explainable artificial intelligence, and uncertainty- or reliability-aware prediction.

2.1. Machine Learning for Aluminum Alloy Mechanical-Property Prediction

ML has been widely explored to predict the mechanical properties of aluminum alloys because it can capture nonlinear composition–process–property relationships. Soofi et al. examined ML-assisted alloy design using a small dataset of wrought aluminum alloys and showed that random forest models could predict multiple mechanical and technological properties from compositional and temper information Soofi et al. (2022). Bhat et al. improved property prediction by partitioning the aluminum alloy data into feature-similar classes and training class-based regressors, thereby enhancing interpretability and predictive performance Bhat et al. (2023). Similarly, Fu et al. used literature data on alloy composition, grain size, hardness, and tensile strength and reported strong tensile-strength predictions using XGBoost Fu et al. (2023). Li et al. combined BP, GA-optimized BP, and RBF neural networks to predict and optimize the mechanical properties of 7xxx aluminum alloys Li et al. (2022).

Recent studies on aluminum alloys have also incorporated process variables, microstructure, and physically meaningful descriptors. Choi et al. developed a data-driven inverse-design framework for T5-aged HPDC aluminum alloys using Ridge regression, random forest, XGBoost, and multilayer perceptron models to predict YS, UTS, and El Choi et al. (2026). Yang et al. proposed a data-extension-based process–composition–property framework for die-casting aluminum alloys and used support vector regression to predict YS, UTS, and El Yang et al. (2024). Murakami et al. used deep learning with microstructural images to predict tensile properties of aluminum alloys with varying silicon content, showing that microstructural-image diversity strongly affects prediction performance Murakami et al. (2025). Shah et al. developed an explainable ML framework for predicting tensile properties of aluminum alloys processed by laser powder bed fusion using composition, process parameters, grain size, elemental descriptors, feature screening, and SHAP interpretation Shah et al. (2025). Pichlmann et al. used an LLM-based data-extraction pipeline together with physics-based descriptors, CALPHAD features, and genetic-algorithm-based feature selection to predict the yield strength of aluminum alloys Pichlmann et al. (2025). In a process-specific study, Mohammed Yusuf and Rathinam predicted microhardness and tensile strength of friction stir-processed AA6061-O using process parameters such as tool rotational speed, travel speed, and number of passes Mohammed Yusuf and Rathinam (2026).

In addition to conventional aluminum alloys, related lightweight-alloy studies further demonstrate the usefulness of ML in linking microstructure and mechanical response. Ai et al. used ML and SHAP analysis to study microstructure–mechanical property relationships in Mg–Al alloys Ai et al. (2025), while Liu et al. combined molecular dynamics simulations with ML models to predict Young’s modulus and ultimate tensile strength of graphene-reinforced aluminum nanocomposites Liu et al. (2022). Although these works are not direct benchmarks for heat-treatable aluminum alloy datasets, they support the broader role of ML in capturing complex relationships among composition, processing, microstructure, and mechanical properties.

2.2. Feature Engineering and Feature Selection

Feature engineering and feature selection are important in materials-property modeling because the prediction quality depends strongly on how composition, processing history, microstructure, and physical descriptors are represented. Shah et al. used alloy composition, processing variables, grain size, and elemental properties with PCC and feature-importance-based screening Shah et al. (2025), whereas Pichlmann et al. combined LLM-extracted composition and process information with physics-based elemental and CALPHAD descriptors, followed by genetic algorithm-based feature selection Pichlmann et al. (2025). Owens et al. emphasized that descriptor construction, descriptor transformation, and ML algorithm selection are key stages in structure–property modeling for grain-boundary energy prediction Owens et al. (2025). Broader materials-informatics studies have also shown that physically meaningful and compact descriptors improve prediction and interpretation. Banik et al. evaluated generalized feature importance for elastic-property prediction and highlighted the need for consistent feature ranking across ML models Banik et al. (2024). Lin et al. proposed a physically meaningful feature engineering and sparsification strategy for catalyst design Lin et al. (2024), whereas Kim and Lee used high-throughput alloy experiments and XAI to identify key descriptors governing electrical resistivity Kim and Lee (2025).

Feature selection is also a well-established machine-learning problem, especially for high-dimensional datasets where redundant variables may increase complexity and reduce generalization. Nogales and Benalcázar reviewed several feature selection and extraction methods, including MRMR, sequential selection, Relief-F, decision-tree-based selection, and neighborhood component analysis, showing their effects on accuracy and processing time Nogales and Benalcázar (2023). Metaheuristic and hybrid feature-selection strategies have also been explored in other domains; Penmatsa et al. used bat optimization for wrapper-based feature selection, whereas Mallidi and Ramisetty proposed nature-inspired and embedded-filter feature-selection strategies for selecting compact and informative feature subsets Penmatsa et al. (2021); Mallidi and Ramisetty (2025a,b). Although these studies are not specific to aluminum alloy design, they support the broader relevance of hybrid feature selection that balances relevance, redundancy reduction, and model performance. In the present study, this idea is adapted to predict aluminum alloy properties through a physics-retained hybrid feature-selection strategy that reduces the descriptor space while preserving metallurgically important heat-treatment variables.

2.3. Explainable Artificial Intelligence in Materials Informatics

Explainable artificial intelligence (XAI) is increasingly important in materials informatics because black-box ML models may provide accurate predictions without revealing the physical factors responsible for those predictions. In studies of aluminum alloys, Shah et al. used explainable ML to interpret tensile-property predictions for laser powder bed fusion-processed aluminum alloys Shah et al. (2025). Choi et al. incorporated model interpretation into a data-driven inverse-design framework for HPDC aluminum alloys Choi et al. (2026), and Pichlmann et al. emphasized that physics-based features can improve the interpretability of yield-strength prediction Pichlmann et al. (2025). In broader materials and engineering applications, Alaca et al. used SHAP, LIME, and partial dependence plots to interpret the effect of alloying elements on titanium-alloy density Alaca et al. (2025). Kim and Lee used XAI to identify average valence electron concentration, electronegativity difference, and mixing entropy as important descriptors for alloy electrical resistivity Kim and Lee (2025). Banik et al. used SHAP-derived feature importance to generate task-specific descriptor rankings for elastic-property prediction Banik et al. (2024).

However, XAI outcomes can be affected by feature dependence and multicollinearity. Salih et al. showed that correlated features may affect the robustness of XAI-based rankings and proposed a correction strategy for dependent features Salih et al. (2024). Nazat et al. used SHAP and LIME for anomaly-detection explanation and XAI-based feature selection Nazat et al. (2024), whereas Shakil et al. combined feature selection and SHAP-based explanations in a medical diagnostic framework Shakil et al. (2024). These studies reinforce the need to interpret ML models carefully, especially when input variables are correlated. Therefore, in the present study, SHAP is used to examine target-specific model attributions associated with composition and processing descriptors for YS, UTS, and El, rather than to infer causal material effects.

2.4. Uncertainty-Aware and Reliability-Aware Prediction

Most ML-based aluminum alloy studies focus mainly on point prediction, where the model outputs a single estimated value for YS, UTS, El, hardness, or tensile strength. Although such predictions are useful for accuracy comparison, they are insufficient for candidate screening because they do not indicate prediction confidence. This limitation is important for literature-curated datasets, where heterogeneity may arise from differences in composition reporting, processing conditions, heat-treatment schedules, sample preparation, and mechanical testing protocols. Pichlmann et al. noted the challenges of limited data availability and generalization in literature-extracted aluminum alloy datasets Pichlmann et al. (2025), whereas Soofi et al. emphasized the practical limitations of small datasets in commercial wrought aluminum alloy modeling Soofi et al. (2022).

Uncertainty-aware learning can distinguish reliable predictions from less reliable ones, which is useful when screening candidate compositions or heat-treatment conditions outside well-represented regions of the training data Tavazza et al. (2021). Ensemble uncertainty, Bayesian models, Gaussian-process regression, quantile regression, and conformal prediction are commonly used for uncertainty estimation Zhang et al. (2022). Conformal prediction is particularly attractive because it can provide distribution-free prediction intervals and can be applied on top of trained ML models Varivoda et al. (2023). In the present study, uncertainty estimation is integrated with candidate screening through an ensemble-conformal strategy. Candidate conditions are ranked using predicted YS, UTS, and El while penalizing sample-specific uncertainty, thereby extending conventional point-prediction-based alloy screening toward reliability-aware candidate selection.

2.5. Synthesis of Related Work

The reviewed literature indicates that ML is effective in predicting the mechanical properties of aluminum alloys, particularly when composition, process variables, microstructural descriptors, physics-based features, feature selection, and explainability are incorporated. However, most existing studies still focus on point prediction, use either all descriptors or target-wise feature selection, and provide limited integration of physics-retained descriptor selection, multi-target strength–ductility prediction, explainability, and uncertainty-aware candidate screening within a single framework.

Accordingly, the present study is positioned as a physics-retained, explainable, and uncertainty-aware multi-target ML framework for heat-treatable aluminum alloys. It reduces the descriptor space from 17 to 12 variables while preserving solution treatment temperature, aging temperature, and aging time; predicts YS, UTS, and El in a unified modeling pipeline; interprets target-specific model attributions using SHAP; and ranks candidate composition–process conditions using ensemble-conformal uncertainty.

3. Materials and Methods

3.1. Dataset and Preprocessing

This study used the aged forged aluminum-alloy dataset reported by Chang. The source dataset comprised 524 records containing alloy-composition variables, heat-treatment parameters, and one or more experimentally reported mechanical properties. However, not all records contained complete measurements for the three target properties. After excluding 72 records with at least one missing value in YS, UTS, or El, 452 complete-target records remained and were used for multi-target modeling. These retained records provided both material-composition and processing descriptors relevant to the strength–ductility response of heat-treatable aluminum alloys.

The target properties considered in this study are yield strength (YS), ultimate tensile strength (UTS), and elongation (El). YS and UTS describe strength-related behavior, while El represents ductility. Since strength and ductility are often coupled and may exhibit a trade-off in heat-treatable aluminum alloys, simultaneous prediction of these three properties is important for evaluating the overall mechanical performance of candidate alloy conditions.

The original input space consisted of 17 variables, including 14 alloying-element composition variables and three heat-treatment variables. The composition variables were Si, Fe, Cu, Mn, Mg, Cr, Zn, V, Ti, Zr, Li, Ni, Be, and Sc. The heat-treatment variables were solution treatment temperature ($T_{\mathrm{sol}}$), aging temperature ($T_{\mathrm{age}}$), and aging time ($t_{\mathrm{age}}$). The target variables were YS, UTS, and El. A summary of the input and target variables is presented in Table 1.

Accordingly, the prediction task is formulated as a multi-output regression problem that maps the input feature vector $\mathbf{X}=[\mathrm{Si},\mathrm{Fe},\mathrm{Cu},\mathrm{Mn},\mathrm{Mg},\mathrm{Cr},\mathrm{Zn},\mathrm{V},\mathrm{Ti},\mathrm{Zr},\mathrm{Li},\mathrm{Ni},\mathrm{Be},\mathrm{Sc},T_{\mathrm{sol}},T_{\mathrm{age}},t_{\mathrm{age}}]$ to the target vector $\mathbf{Y}=[YS,UTS,El]$, where YS and UTS are expressed in MPa and El is expressed as percentage elongation. In contrast to single-target regression, this work treats YS, UTS, and El as related responses and evaluates them within a unified multi-target prediction framework.

Table 2 reports the descriptive statistics, while Figure 1 shows the corresponding empirical distributions. The variables are reported using the units listed in Table 2. Several minor-element variables were strongly zero-inflated, and aging time had a pronounced right-skewed distribution, with a maximum of 1440 h compared with a median of 16 h. The target properties also covered broad and nonuniform ranges. These characteristics should be considered when interpreting feature rankings, model explanations, and generalization performance.

The retained dataset spans multiple alloying elements, heat-treatment conditions, and both strength- and ductility-related targets. This supports an internal evaluation of composition–process–property modeling and strength–ductility trade-offs. However, the number of independent compositions is limited, and the dataset should not be considered representative of the complete design space of heat-treatable aluminum alloys.

The 17 input variables contained no missing values. In the source data, YS, UTS, and El contained 49, 15, and 37 missing values, respectively; because some records lacked more than one target, 72 records had at least one missing target and were excluded. The remaining 452 complete-target records contained eight rows belonging to four exact-duplicate pairs, corresponding to four surplus duplicate rows. The locked row-wise analysis retained these records. One exact-duplicate pair crossed the proper-training and calibration partitions, and one crossed the proper-training and final-test partitions.

Among the 452 retained records, 50 rows shared an identical 17-variable predictor condition, giving 420 unique composition–processing conditions. The source file did not include a row-level literature-source identifier or a verified alloy-family label. Consequently, source-aware or alloy-family-aware validation could not be performed without inferring unavailable metadata. All input variables were numerical, and no categorical encoding was required.

Using a master random seed of 42, the 452 complete-target records were assigned to three mutually exclusive partitions. The three partitions were mutually exclusive at the row-index level but were not fully independent at the composition–condition level because exact duplicates and repeated compositions crossed partition boundaries, as documented below. An initial 80:20 split produced 361 development records and 91 final-test records. The 361 development records were then divided into 270 proper-training records and 91 calibration records. Relative to the complete-target dataset, the resulting proportions were therefore 59.73% proper training, 20.13% calibration, and 20.13% final testing. No separate permanent validation partition was used; validation was performed internally through cross-validation within the proper-training partition.

Feature ranking, model-family selection, and hyperparameter optimization were restricted to the 270 proper-training records. The evaluation used five-fold outer cross-validation repeated five times, producing 25 outer-validation evaluations for each prespecified feature count. Within each outer-training fold, a shuffled five-fold inner cross-validation was used to select the hybrid-score weight $\alpha$, model family, and model hyperparameters for the corresponding fixed feature count. The final feature count, $k=12$, was subsequently chosen because it produced the highest mean macro-$R^2$ across the 25 outer-validation evaluations. Thus, feature-count selection was not enclosed within an additional inner selection layer, and the reported comparison involving the selected $k=12$ result is interpreted as a post-selection exploratory comparison rather than an unbiased fully nested estimate. After selecting $k=12$, a final five-fold inner search using all 270 proper-training records locked the model family, $\alpha$, and model hyperparameters. The 91 calibration records were accessed only after this lock for conformal calibration, and the 91 final-test records were accessed only for final point-prediction and interval evaluation.

All data-dependent operations were restricted to the corresponding training fold. This included target standardization for embedded ranking, Pearson-correlation calculation, bootstrap ExtraTrees importance estimation, feature ranking, feature-subset construction, and hyperparameter optimization. When standardization was required by a scale-sensitive benchmark model, the scaler was fitted only on the relevant training fold and then applied to its validation fold. Tree-based models were trained on their original input scales. Consequently, neither the calibration partition nor the final-test partition influenced feature or model selection.

The effective model-development sample consisted of 270 proper-training records and must be regarded as modest for a multi-target materials-informatics problem. The dataset is sufficient for an internally evaluated comparison of tree-ensemble configurations, but it is not sufficient to establish generalization to compositionally novel alloy families or industrial processing conditions. The candidate model families were alternative configurations evaluated within the same model-selection procedure; they did not increase the amount of independent training information available to the final model.Overfitting risk was reduced by restricting all feature-ranking and tuning operations to the relevant training folds, using five-fold inner validation, evaluating alternative tree-depth and minimum-leaf settings during tuning, limiting the randomized-search budget, using ensemble averaging, evaluating variability across 25 outer splits, and retaining separate calibration and final-test record partitions. Nevertheless, the post-selection choice of $k=12$, the presence of repeated composition groups, and the reduction in performance under composition-group-disjoint validation indicate that residual selection optimism and composition dependence cannot be excluded. External validation is therefore required before broader generalization claims can be made.

Because the source dataset did not provide row-level source or alloy-family identifiers, an additional composition-group-disjoint sensitivity analysis was performed. A composition-group identifier was constructed from the 14 alloying-element variables after rounding each reported composition to three decimal places. This produced 138 composition groups across the 452 retained records; 66 groups contained multiple records, 380 records belonged to repeated-composition groups, and the largest group contained 37 records. The permanent row-wise split shared 31 composition groups between proper training and calibration, 36 between proper training and final testing, and 19 between calibration and final testing.

To quantify the sensitivity of performance to this dependence, the locked 12-feature ExtraTrees specification was evaluated using 25 composition-group-disjoint splits within the 270-record proper-training partition. Each split held out approximately 20% of the 106 proper-training composition groups using GroupShuffleSplit; all records belonging to the same rounded-composition group were assigned to either training or validation, never both. The selected features and locked ExtraTrees hyperparameters were kept fixed so that the analysis isolated the effect of row-wise versus composition-disjoint partitioning. This was a post-selection generalization diagnostic and was not used to reselect the feature count, model family, or hyperparameters.

3.2. Matched Full-Descriptor and Physics-Retained Modeling Framework

The formal comparison used two matched pipelines. The proposed pipeline used the selected 12-feature physics-retained subset, whereas the comparator used all 17 available descriptors. Both pipelines were evaluated using identical outer-training and outer-validation indices, the same inner cross-validation procedure, the same model families, the same hyperparameter ranges, and the same search budget. This matched design isolated the effect of descriptor reduction while avoiding comparisons between independently generated splits or differently tuned models.

Four model implementations were considered in both pipelines. The first was a native multi-output ExtraTrees regressor, implemented by fitting ExtraTreesRegressor directly to the three-column response matrix. This implementation learns a common tree-partition structure and produces vector-valued leaf predictions for YS, UTS, and El. The second was a native multi-output random-forest regressor implemented in the same manner using RandomForestRegressor.

The third implementation was an independent target-wise ExtraTrees model constructed using MultiOutputRegressor. This wrapper cloned and fitted one separate ExtraTrees regressor for each of YS, UTS, and El. It therefore served as both the wrapper-based and independent target-wise comparator and did not directly model dependence among the three responses. Under the common base-estimator configuration used here, separately coding three target-wise ExtraTrees regressors would reproduce the same independent formulation.

The fourth implementation was an ExtraTrees regressor chain. Unlike the independent wrapper, the chain allowed preceding target predictions to be used as additional inputs for subsequent targets. The chain order was fixed as YS–UTS–El. Thus, only the regressor-chain implementation explicitly incorporated preceding target predictions, whereas the native tree ensembles used shared tree structures without explicitly parameterizing conditional target dependence.

The benchmark set deliberately emphasized tree-based ensemble models because the proper-training partition contained only 270 observations and all predictors were structured tabular variables. Native random forest and ExtraTrees provided two nonlinear ensemble formulations, while the independent target-wise implementation and regressor chain evaluated alternative multi-output strategies. A deep neural-network benchmark was not included because architecture and regularization selection with this sample size would impose a substantial model-selection burden relative to the number of independent observations. Moreover, recent deep-learning studies discussed in the related literature frequently use microstructural images, data-augmentation procedures, or different alloy systems and therefore cannot be compared directly with the present tabular dataset. Comparisons with gradient-boosted ensembles and appropriately regularized neural models using larger external datasets are recommended for future work.

Within every outer-training partition, the model implementation and its hyperparameters were selected through five-fold inner cross-validation. For the 12-feature pipeline, feature ranking and $\alpha$ selection were also repeated using only the corresponding outer-training records. The 17-feature comparator did not remove any descriptors. All four implementations were evaluated using identical inner folds, feature budgets, hyperparameter ranges, and search budgets. The structures of the two matched pipelines are summarized in Table 3.

Let ${\mathbf{x}}_i\in {\mathbb{R}}^d$ denote the input feature vector of alloy condition i, where $d=12$ for the proposed pipeline and $d=17$ for the full-descriptor comparator. The corresponding multi-target response was defined as ${\mathbf{y}}_i=[Y{S}_i,UT{S}_i,E{l}_i]$. The regression model learned the mapping

$${\widehat{\mathbf{y}}}_i=f({\mathbf{x}}_i)=[{\widehat{YS}}_i,{\widehat{UTS}}_i,{\widehat{El}}_i].$$

(1)

The native multi-output ExtraTrees and random-forest implementations fitted the three-column response matrix directly and produced all three predictions from a shared ensemble structure. In contrast, the independent MultiOutputRegressor implementation fitted three separate estimators without exchanging target information.

For the regressor-chain implementation, the target order was fixed as $\mathrm{YS}\to \mathrm{UTS}\to \mathrm{El}$. Its inference-stage predictions were therefore generated as

$$\begin{array}{cc}\hfill {\widehat{YS}}_i& =f_1\left({\mathbf{x}}_i\right),\hfill \\ \hfill {\widehat{UTS}}_i& =f_2\left({\mathbf{x}}_i,{\widehat{YS}}_i\right),\hfill \\ \hfill {\widehat{El}}_i& =f_3\left({\mathbf{x}}_i,{\widehat{YS}}_i,{\widehat{UTS}}_i\right).\hfill \end{array}$$

(2)

During chain fitting, the measured preceding targets were appended to the training inputs, whereas during inference they were replaced by the corresponding preceding predictions. The YS–UTS–El order was fixed throughout the reported comparison and was not changed using calibration or final-test performance.

The literature-derived methods discussed in Section 2 were adapted and implemented in-house to provide methodological context; they were not direct reproductions of the original studies. The matched numerical comparison was restricted to the 12- and 17-feature pipelines described above; because the 12-feature count was selected from the same outer-validation results, this comparison is interpreted as exploratory rather than confirmatory.

3.3. Physics-Retained Hybrid Feature Selection

Feature selection was performed to obtain a compact and physically meaningful input subset for strength–ductility prediction. Because the original dataset contained only 17 input variables, the objective of feature selection was not merely to reduce dimensionality but also to remove weak or redundant descriptors while retaining variables that are metallurgically important for heat-treatable aluminum alloys. Therefore, a physics-retained hybrid feature-selection strategy was adopted.

The proposed strategy combines three complementary criteria: statistical relevance, embedded model-based importance, and physics-guided retention. Statistical relevance was used to estimate the direct association between each input variable and the target properties. Embedded importance was obtained from tree-based ensemble models, which can capture nonlinear feature–property relationships and interactions among alloying elements and heat-treatment variables. Physics-guided retention was then applied to preserve heat-treatment descriptors that are essential for interpreting the mechanical response of precipitation-strengthened aluminum alloys. The role of each criterion is summarized in Table 4.

To formalize the statistical relevance criterion, the relationship between each input feature $x_j$ and each target property $y_k$ was first estimated. Because the target vector contained three mechanical properties, $y_k\in \{YS,UTS,El\}$, the multi-target relevance score of feature $x_j$ was computed using Eq. (3), where $\rho (x_j,y_k)$ denotes the Pearson correlation coefficient between feature $x_j$ and target $y_k$. The absolute value was used because both positive and negative correlations can be meaningful in strength–ductility modeling. Thus, $R_j$ represents the average linear association of a feature with YS, UTS, and El.

$$R_j={\displaystyle \frac{1}{3}}{\displaystyle \sum _{k=1}^3}\left|\rho (x_j,y_k)\right|$$

(3)

To account for nonlinear predictive contribution, embedded importance was estimated separately inside each training fold. The three target variables were first standardized using the means and standard deviations calculated from that training fold. For each of $B=5$ bootstrap samples, an ExtraTrees regressor containing 200 trees and using $\max \text{\_}\mathrm{features}=\sqrt{p}$ was fitted to the standardized multi-target response. Let $I_j^{(b)}$ denote the impurity-based importance of feature $x_j$ from bootstrap run b. The mean embedded importance was calculated as

$$E_j={\displaystyle \frac{1}{B}}{\displaystyle \sum _{b=1}^B}{I}_j^{(b)},B=5.$$

(4)

Equation (4) represents mean embedded importance and not feature-selection stability. Genuine stability was evaluated separately from the feature subsets selected across the 25 outer cross-validation evaluations.

Before combining the Pearson relevance and embedded importance, each criterion was independently min–max normalized across the 17 input variables:

$${\tilde{R}}_j={\displaystyle \frac{R_j-{\min }_{\ell }{R}_{\ell }}{{\max }_{\ell }{R}_{\ell }-{\min }_{\ell }{R}_{\ell }}},{\tilde{E}}_j={\displaystyle \frac{E_j-{\min }_{\ell }{E}_{\ell }}{{\max }_{\ell }{E}_{\ell }-{\min }_{\ell }{E}_{\ell }}}.$$

(5)

When the maximum and minimum of a criterion were equal, its normalized values were set to zero. The hybrid feature score was then calculated as

$$H_j(\alpha )=\alpha {\tilde{R}}_j+(1-\alpha ){\tilde{E}}_j.$$

(6)

The candidate values $\alpha \in \{0,0.25,0.50,0.75,1.00\}$ were evaluated only through inner cross-validation. Here, $\alpha =0$ represents embedded-importance-only ranking, whereas $\alpha =1$ represents Pearson-relevance-only ranking. The final proper-training search selected $\alpha =0.25$, thereby assigning a weight of 0.25 to normalized Pearson relevance and 0.75 to normalized mean embedded importance.

Sensitivity to $\alpha$ was summarized using the existing final $k=12$ proper-training search. For each candidate value of $\alpha$, the evaluated configuration having the lowest mean five-fold inner-validation macro-NRMSE was identified. The frequency with which each value of $\alpha$ was selected across the 25 outer-training searches at $k=12$ was also recorded. Neither the calibration partition nor the final-test partition was accessed for this analysis. Because RandomizedSearchCV jointly sampled $\alpha$ and the model hyperparameters, the analysis represents a search-based sensitivity diagnostic rather than an isolated one-factor experiment.

Physics-guided retention was applied by including ${\mathcal{F}}_{\mathrm{phys}}=\{T_{\mathrm{sol}},T_{\mathrm{age}},t_{\mathrm{age}}\}$ in every evaluated feature subset. For a total feature budget k, the remaining $k-3$ variables were selected from the non-mandatory descriptors in descending order of $H_j(\alpha )$:

$${\mathcal{F}}_k={\mathcal{F}}_{\mathrm{phys}}\cup {\mathrm{Top}}_{k-3}\left[H_j(\alpha ):x_j\notin {\mathcal{F}}_{\mathrm{phys}}\right].$$

(7)

The three heat-treatment descriptors were structurally retained in every evaluated subset irrespective of their hybrid-score ranks. Consequently, their inclusion frequency cannot be interpreted as empirical evidence of feature importance or selection stability. Their subsequent SHAP attributions indicate only how the locked model used these prespecified variables.

Feature budgets from $k=3$ to $k=17$ were evaluated independently. Exact hybrid-score ties were resolved lexically after inserting the three mandatory heat-treatment descriptors.

The complete procedure of the proposed physics-retained hybrid feature-selection strategy is summarized in Algorithm 1. The algorithm combines statistical relevance, mean embedded importance, and physics-guided retention to obtain a compact feature subset while preserving metallurgically important heat-treatment descriptors.

Algorithm 1 Training-fold-contained physics-retained feature selection

Require: Training-fold features ${\mathbf{X}}_{\mathrm{tr}}$, targets ${\mathbf{Y}}_{\mathrm{tr}}$, feature budget k, hybrid weight $\alpha$, and ${\mathcal{F}}_{\mathrm{phys}}=\{T_{\mathrm{sol}},T_{\mathrm{age}},t_{\mathrm{age}}\}$ Ensure: Selected feature subset ${\mathcal{F}}_k$ 1: Standardize each target using only ${\mathbf{Y}}_{\mathrm{tr}}$ 2: for each input feature $x_j$ do 3:     Compute $R_j$ as the mean absolute Pearson correlation with YS, UTS, and El 4: end for 5: for $b=1$ to 5 do 6:     Draw a bootstrap sample from the current training fold 7:     Fit a 200-tree ExtraTrees regressor 8:     Extract the embedded importance $I_j^{(b)}$ for every feature 9: end for 10: for each input feature $x_j$ do 11:     Compute the mean embedded importance $E_j$ using Eq. (4) 12:     Compute ${\tilde{R}}_j$ and ${\tilde{E}}_j$ using Eq. (5) 13:     Compute $H_j(\alpha )$ using Eq. (6) 14: end for 15: Initialize ${\mathcal{F}}_k\leftarrow {\mathcal{F}}_{\mathrm{phys}}$ 16: Add the highest-ranked non-mandatory variables until $|{\mathcal{F}}_k|=k$ 17: return ${\mathcal{F}}_k$

The final configuration used $k=12$ and $\alpha =0.25$. The selected descriptors were $T_{\mathrm{sol}}$, $T_{\mathrm{age}}$, $t_{\mathrm{age}}$, Zn, Mg, Cu, Ti, Fe, Zr, Si, Sc, and Cr. The 12-feature count was selected because it produced the highest observed mean macro-$R^2$ across the 25 outer-validation evaluations performed within the 270-record proper-training partition. Neither the calibration partition nor the final-test partition was used for this selection. However, because the same outer-validation results were used both to select $k=12$ and to summarize its cross-validation performance, the resulting comparison is post-selection and exploratory. A secondary one-standard-error analysis identified $k=8$ as a more parsimonious alternative, but the reported model retained $k=12$ according to the maximum-observed-mean-$R^2$ criterion.

3.4. Hyperparameter Optimization

Model and feature-selection hyperparameters were optimized using RandomizedSearchCV. Four physics-retained model families were considered: native multi-output ExtraTrees, native multi-output random forest, independent target-wise ExtraTrees implemented through a multi-output wrapper, and an ExtraTrees regressor chain using the target order YS–UTS–El.

For each fixed feature count and model family, 24 randomly sampled configurations were evaluated using five-fold inner cross-validation. Thus, the final fixed-$k=12$ lock evaluated $4\times 24=96$ candidate configurations and generated $96\times 5=480$ inner-fold evaluation records. Feature counts from 3 to 17 were evaluated independently rather than by cloning one fitted configuration across different subset sizes.

For target k, the normalized validation error was calculated using the target standard deviation in the corresponding inner-validation fold:

$${\mathrm{NRMSE}}_{k,\mathrm{val}}={\displaystyle \frac{\sqrt{\frac{1}{n_{\mathrm{val}}}{\sum }_{i=1}^{n_{\mathrm{val}}}{\left(y_{ik}-{\widehat{y}}_{ik}\right)}^2}}{s_{k,\mathrm{val}}}},$$

(8)

where $s_{k,\mathrm{val}}$ is the sample standard deviation of target k in the corresponding validation fold. Validation targets were used only to calculate validation metrics, as in conventional cross-validation; they were not used for model fitting or feature ranking. For final-test reporting, the corresponding final-test target standard deviation was used. This predictive-metric normalization is separate from the uncertainty penalty in Eq. (33), which uses the proper-training target standard deviations. The search objective was the mean target-wise NRMSE:

$$\mathcal{L}(\theta )={\displaystyle \frac{1}{3}}\left({\mathrm{NRMSE}}_{\mathrm{YS}}+{\mathrm{NRMSE}}_{\mathrm{UTS}}+{\mathrm{NRMSE}}_{\mathrm{El}}\right).$$

(9)

The configuration with the lowest mean inner-validation NRMSE was selected. Table 5 reports the complete search space.

The final $k=12$ search selected the native multi-output ExtraTrees model with 800 trees, unrestricted maximum depth, $\max \text{\_}\mathrm{features}=0.7$, $\min \text{\_}\mathrm{samples}\text{\_}\mathrm{leaf}=1$, and $\alpha =0.25$. All 270 proper-training records were used to refit this locked configuration before conformal calibration.

The master seed was 42 and was used for the permanent data split and repeated outer cross-validation. Fold-specific ranking, estimator, and search seeds were generated deterministically from this master seed and retained with the corresponding index hashes. The final five-fold inner split used seed 70042, the final randomized search used seed 72043, and the 20 bootstrap ensemble members used seeds 60043–60062. Calibration and final-test records were unavailable to all hyperparameter-selection operations.

3.5. Computational Environment and Reproducibility

All analyses were executed in Python 3.8.20 on a Linux x86_64 workstation equipped with AMD Ryzen 7 5800X, 16 logical CPU cores, and 128 GB of system memory. The scikit-learn tree ensembles and SHAP calculations were executed on the CPU; no GPU acceleration was used. The principal software versions were NumPy 1.24.3, pandas 2.0.3, SciPy 1.10.1, scikit-learn 1.3.2, SHAP 0.44.1, Matplotlib 3.7.5, joblib 1.4.2, and openpyxl 3.1.5. The master random seed was 42, while deterministic derived seeds were used for the inner folds, randomized search, bootstrap ensemble, and SHAP-stability analysis as reported in the preceding sections.

3.6. Explainability and Uncertainty Estimation

Explainability analysis was performed to characterize the contributions of alloying elements and heat-treatment variables to the model predictions. Because prediction accuracy alone is insufficient for interpreting model behavior, SHAP analysis was used to examine how each selected descriptor contributed to the predictions of YS, UTS, and El. This target-wise analysis allowed differences in the model’s use of the same descriptor across strength and ductility predictions to be examined.

SHAP is based on the concept of Shapley values from cooperative game theory. In the present context, each input feature is treated as a contributor to the model prediction, and its Shapley value represents its marginal contribution to the predicted target property. For a trained model f and a feature set $\mathcal{M}$, the SHAP value of feature j for a given prediction is defined in Eq. (10), where $\mathcal{S}$ is a subset of features excluding feature j, $\mathcal{M}$ is the complete set of input features, and $f(\mathcal{S})$ denotes the model prediction using the feature subset $\mathcal{S}$. A positive SHAP value indicates that the feature increases the predicted target value, whereas a negative SHAP value indicates that the feature decreases the prediction.

$${\varphi }_j=\sum _{\mathcal{S}\subseteq \mathcal{M}\setminus \{j\}}{\displaystyle \frac{|\mathcal{S}|!(|\mathcal{M}|-|\mathcal{S}|-1)!}{|\mathcal{M}|!}}\left[f(\mathcal{S}\cup \{j\})-f(\mathcal{S})\right]$$

(10)

For global interpretation, the mean absolute SHAP value was computed for each selected feature, as shown in Eq. (11), where ${\varphi }_{ij}$ is the SHAP value of feature j for the ith sample and n is the number of samples. A higher value of $G_j$ indicates a larger average absolute model attribution for feature j.

$$G_j={\displaystyle \frac{1}{n}}{\displaystyle \sum _{i=1}^n}|{\varphi }_{ij}|$$

(11)

SHAP values were calculated directly for the locked 12-feature native multi-output ExtraTrees model using shap.TreeExplainer with feature_perturbation = ’tree_path_dependent’ and check_additivity = False. No external background matrix was supplied; the explainer used the internal tree-path distribution learned by the fitted ensemble. The evaluation data comprised the 91 row-wise held-out final-test records, and the resulting attributions were separated target-wise for YS, UTS, and El.

The tree-path-dependent formulation does not causally disentangle correlated composition and heat-treatment descriptors. Correlated variables may share or redistribute attribution according to the fitted tree structure. Therefore, the reported SHAP values describe how the locked model generated its predictions within the observed data distribution and do not establish isolated physical or causal effects.

Sensitivity of the global SHAP ranking to the finite evaluation subset was assessed using 1000 bootstrap resamples of the 91 evaluation records with seed 81042. For each resample, mean absolute SHAP values and feature ranks were recalculated, and the resulting ranking was compared with the full-evaluation ranking using Spearman correlation.

Target-wise SHAP analysis was performed separately for YS, UTS, and El. This enabled the identification of features that mainly contribute to strength prediction and those that affect ductility. For example, some alloying elements or heat-treatment variables may increase YS and UTS while reducing El, reflecting the strength–ductility trade-off. Therefore, SHAP analysis was used not only to rank important variables but also to interpret the direction and magnitude of their effects on each target property.

In addition to explainability, uncertainty estimation was incorporated to evaluate the reliability of the predicted mechanical properties. Of the 524 records in the source dataset, 452 records contained complete values for YS, UTS, and El and were used for multi-target modeling. Using random seed 42, these records were divided into 270 proper-training records, 91 calibration records, and 91 final-test records. Feature selection, hyperparameter selection, and model fitting were restricted to the proper-training partition. The calibration partition was used exclusively for conformal calibration, while the final-test partition was used only for final evaluation. The three partitions were mutually exclusive at the row-index level; the exact-duplicate and repeated-composition dependence across partitions was audited as described in Section 3.1.

An ensemble of $M=20$ ExtraTrees models was constructed by training each ensemble member on a bootstrap sample drawn from the 270 proper-training records. For sample i and target k, the ensemble mean prediction was calculated as

$${\widehat{y}}_{ik}^{\mathrm{ens}}={\displaystyle \frac{1}{M}}{\displaystyle \sum _{m=1}^M}{\widehat{y}}_{ik}^{(m)},$$

(12)

where ${\widehat{y}}_{ik}^{(m)}$ denotes the prediction produced by ensemble member m. The sample-specific ensemble standard deviation was calculated using

$${\sigma }_{ik}^{\mathrm{ens}}=\sqrt{{\displaystyle \frac{1}{M}}{\displaystyle \sum _{m=1}^M}{\left({\widehat{y}}_{ik}^{(m)}-{\widehat{y}}_{ik}^{\mathrm{ens}}\right)}^2}.$$

(13)

The ensemble standard deviation measures the disagreement among the bootstrap models and allows the uncertainty interval to vary across samples. Locally normalized nonconformity scores were calculated on the independent calibration partition as

$$a_{ik}^{\mathrm{cal}}={\displaystyle \frac{\left|y_{ik}^{\mathrm{cal}}-{\widehat{y}}_{ik}^{\mathrm{ens},\mathrm{cal}}\right|}{{\sigma }_{ik}^{\mathrm{ens},\mathrm{cal}}+\epsilon }},\epsilon ={10}^{-8},$$

(14)

where $y_{ik}^{\mathrm{cal}}$ is the measured value and ${\widehat{y}}_{ik}^{\mathrm{ens},\mathrm{cal}}$ is the ensemble mean prediction for calibration sample i and target k.

For a desired marginal coverage level of $1-\delta$, the finite-sample conformal rank was determined as

$$r=\min \left\{n_{\mathrm{cal}},⌈(n_{\mathrm{cal}}+1)(1-\delta )⌉\right\},q_k=a_{(r)k}^{\mathrm{cal}},$$

(15)

where $a_{(r)k}^{\mathrm{cal}}$ denotes the rth ordered calibration score for target k. In the present study, $n_{\mathrm{cal}}=91$, $\delta =0.10$, and therefore

$$r=⌈(91+1)(0.90)⌉=83.$$

(16)

The resulting normalized conformal quantiles were 4.8429 for YS, 4.4340 for UTS, and 5.6609 for El. The sample-specific prediction interval was then defined as

$$I_{ik}=\left[{\widehat{y}}_{ik}^{\mathrm{ens}}-q_k\left({\sigma }_{ik}^{\mathrm{ens}}+\epsilon \right),{\widehat{y}}_{ik}^{\mathrm{ens}}+q_k\left({\sigma }_{ik}^{\mathrm{ens}}+\epsilon \right)\right].$$

(17)

Accordingly, the full interval width for sample i and target k was

$$W_{ik}=2q_k\left({\sigma }_{ik}^{\mathrm{ens}}+\epsilon \right).$$

(18)

Although $q_k$ is fixed for a given target, ${\sigma }_{ik}^{\mathrm{ens}}$ varies across samples. Therefore, the resulting conformal intervals are sample-specific rather than constant-width intervals.

For comparison, the same bootstrap ensemble was also evaluated without conformal calibration. A central 90% Gaussian ensemble-disagreement interval was constructed as

$$I_{ik}^{\mathrm{ens}}=\left[{\widehat{y}}_{ik}^{\mathrm{ens}}-z_{0.95}{\sigma }_{ik}^{\mathrm{ens}},{\widehat{y}}_{ik}^{\mathrm{ens}}+z_{0.95}{\sigma }_{ik}^{\mathrm{ens}}\right],z_{0.95}=1.6449.$$

(19)

This comparator used only the dispersion among the same 20 bootstrap-model predictions and did not use the calibration partition. It therefore represents an uncalibrated ensemble-based uncertainty estimate rather than a distribution-free prediction interval. Comparing it with the ensemble-conformal interval isolates the effect of conformal calibration while keeping the underlying ensemble predictions unchanged.

The target-wise marginal prediction interval coverage probability was calculated as

$${\mathrm{PICP}}_k={\displaystyle \frac{1}{n_{\mathrm{test}}}}{\displaystyle \sum _{i=1}^{n_{\mathrm{test}}}}\mathbb{I}\left(y_{ik}\in I_{ik}\right),$$

(20)

where $\mathbb{I}(·)$ equals 1 when the measured value lies within the corresponding interval and 0 otherwise. Because candidate acceptance depends on all three properties, simultaneous coverage was calculated separately as

$${\mathrm{PICP}}_{\mathrm{sim}}={\displaystyle \frac{1}{n_{\mathrm{test}}}}{\displaystyle \sum _{i=1}^{n_{\mathrm{test}}}}\mathbb{I}\left[\bigcap _{k\in \{\mathrm{YS},\mathrm{UTS},\mathrm{El}\}}\left(y_{ik}\in I_{ik}\right)\right].$$

(21)

The mean prediction interval width for target k was calculated as

$${\mathrm{MPIW}}_k={\displaystyle \frac{1}{n_{\mathrm{test}}}}{\displaystyle \sum _{i=1}^{n_{\mathrm{test}}}}{W}_{ik}.$$

(22)

The nominal 90% coverage level applies separately to the marginal interval for each target and does not imply 90% simultaneous coverage across YS, UTS, and El. Coverage and interval width were therefore interpreted jointly to assess the calibration–sharpness trade-off.

3.7. Evaluation Metrics and Reliability-Aware Candidate Screening

Predictive metrics were calculated separately for YS, UTS, and El. Because YS and UTS errors are expressed in MPa while El errors are expressed in percentage points, raw MAE and RMSE values were reported only target-wise. Raw errors were not averaged across the three targets. Cross-target aggregate comparisons instead used macro-$R^2$, macro-NRMSE, and macro-NMAE.

The coefficient of determination, $R^2$, was used to measure how well the model explains the variance in the target property. For the kth target, $R^2$ was computed using Eq. (23), where $y_{ik}$ and ${\widehat{y}}_{ik}$ are the measured and predicted values of the kth target for the ith sample, respectively, and ${\overline{y}}_k$ is the mean measured value of the kth target.

$$R_k^2=1-{\displaystyle \frac{{\sum }_{i=1}^n{(y_{ik}-{\widehat{y}}_{ik})}^2}{{\sum }_{i=1}^n{(y_{ik}-{\overline{y}}_k)}^2}}$$

(23)

The mean absolute error (MAE) was used to measure the average absolute deviation between measured and predicted values, as shown in Eq. (24).

$${\mathrm{MAE}}_k={\displaystyle \frac{1}{n}}{\displaystyle \sum _{i=1}^n}|y_{ik}-{\widehat{y}}_{ik}|$$

(24)

The root mean square error (RMSE) was used to penalize larger prediction errors more strongly and was calculated using Eq. (25).

$${\mathrm{RMSE}}_k=\sqrt{{\displaystyle \frac{1}{n}}\sum _{i=1}^n{(y_{ik}-{\widehat{y}}_{ik})}^2}$$

(25)

Macro-$R^2$ was calculated as the unweighted mean of the three target-wise $R^2$ values:

$$R_{\mathrm{macro}}^2={\displaystyle \frac{1}{3}}\left(R_{\mathrm{YS}}^2+R_{\mathrm{UTS}}^2+R_{\mathrm{El}}^2\right).$$

(26)

The normalized mean absolute error for target k in evaluation partition $\mathcal{E}$ was defined as

$${\mathrm{NMAE}}_{k,\mathcal{E}}={\displaystyle \frac{{\mathrm{MAE}}_{k,\mathcal{E}}}{s_{k,\mathcal{E}}}},$$

(27)

where $s_{k,\mathcal{E}}$ is the sample standard deviation of target k in the corresponding outer-validation or final-test partition. Macro-NRMSE and macro-NMAE were calculated by averaging their three target-wise normalized values.

The 12- and 17-feature pipelines were statistically compared using their matched results from the same 25 outer evaluations. For metric Q, the paired difference for outer evaluation m was

$$d_m=Q_m^{(12)}-Q_m^{(17)},m=1,\dots ,25.$$

(28)

Because the repeated cross-validation results are not independent, the standard error of the mean paired difference was corrected as

$${\mathrm{SE}}_{\mathrm{corr}}=\sqrt{{\displaystyle \frac{1}{RK}}+{\displaystyle \frac{n_{\mathrm{val}}}{n_{\mathrm{tr}}}}}{s}_d,$$

(29)

where $R=5$ repetitions, $K=5$ folds, $RK=25$ paired evaluations, $n_{\mathrm{val}}/n_{\mathrm{tr}}=54/216=0.25$, and $s_d$ is the sample standard deviation of the 25 paired differences. The corrected statistic was

$$t_{\mathrm{corr}}={\displaystyle \frac{\overline{d}}{{\mathrm{SE}}_{\mathrm{corr}}}},$$

(30)

with $RK-1=24$ degrees of freedom. Two-sided p-values and corrected 95% confidence intervals were reported. For macro-$R^2$, a positive difference favors the 12-feature pipeline; for macro-NRMSE and macro-NMAE, a negative difference favors it.

The final stage of the framework was a retrospective reliability-aware screening of row-wise held-out composition–processing conditions. No synthetic compositions, new alloy candidates, or optimized heat-treatment conditions were generated. Instead, the input conditions contained in the 91 untouched final-test records were screened without using their measured YS, UTS, or El values. Conditions were defined using exact agreement across the 14 reported elemental variables and the three heat-treatment variables. This produced 87 unique held-out conditions because four records duplicated composition–processing conditions represented elsewhere in the final-test partition. All original record indices were retained for traceability.

For each target $k\in \{\mathrm{YS},\mathrm{UTS},\mathrm{El}\}$, the ensemble predictions were min–max normalized within the 91-record held-out screening set as

$${\tilde{y}}_{ik}={\displaystyle \frac{{\widehat{y}}_{ik}^{\mathrm{ens}}-{\min }_{j\in \mathcal{C}}{\widehat{y}}_{jk}^{\mathrm{ens}}}{{\max }_{j\in \mathcal{C}}{\widehat{y}}_{jk}^{\mathrm{ens}}-{\min }_{j\in \mathcal{C}}{\widehat{y}}_{jk}^{\mathrm{ens}}+\epsilon }},$$

(31)

where $\mathcal{C}$ denotes the held-out screening records and $\epsilon ={10}^{-8}$ prevents numerical division by zero. Only predicted values were used in this normalization; the measured final-test target values were not used.

The normalized strength and ductility components were defined as

$$S_i^{\mathrm{str}}=0.5{\tilde{y}}_{i,\mathrm{YS}}+0.5{\tilde{y}}_{i,\mathrm{UTS}},S_i^{\mathrm{duct}}={\tilde{y}}_{i,\mathrm{El}}.$$

(32)

The sample-specific uncertainty penalty was calculated using the full conformal interval widths from Eq. (18):

$$P_i^{\mathrm{unc}}={\displaystyle \frac{1}{3}}\sum _{k\in \{\mathrm{YS},\mathrm{UTS},\mathrm{El}\}}{\displaystyle \frac{W_{ik}}{s_{k,\mathrm{tr}}}},$$

(33)

where $s_{k,\mathrm{tr}}$ is the standard deviation of target k calculated using only the 270 proper-training records. The corresponding values were 120.4265 MPa for YS, 103.8424 MPa for UTS, and 7.2282 percentage points for El.

The baseline reliability-aware score was then calculated as

$$\begin{array}{cc}\hfill {\mathcal{S}}_i& =0.4S_i^{\mathrm{str}}+0.4S_i^{\mathrm{duct}}-0.2P_i^{\mathrm{unc}}\hfill \\ & =0.2{\tilde{y}}_{i,\mathrm{YS}}+0.2{\tilde{y}}_{i,\mathrm{UTS}}+0.4{\tilde{y}}_{i,\mathrm{El}}-0.2P_i^{\mathrm{unc}}.\hfill \end{array}$$

(34)

These baseline weights assign equal aggregate emphasis to strength and ductility. Within the strength component, YS and UTS receive equal weights. The uncertainty weight of 0.2 penalizes conditions with wide prediction intervals. These weights were treated as decision preferences rather than statistically optimized parameters; therefore, their influence was examined through sensitivity analysis.

No minimum-property threshold was imposed because the dataset did not provide an application-specific acceptance standard for YS, UTS, or El. Introducing arbitrary minimum values would have produced an unsupported binary acceptance decision. The conditions were therefore ranked continuously, while composition feasibility and applicability-domain membership were evaluated separately.

Elemental feasibility was evaluated using all 14 reported composition variables. Because Al was not explicitly included in the source dataset, its balance was calculated as

$$x_{i,\mathrm{Al}}^{\mathrm{bal}}=100-\sum _{e\in \mathcal{E}}{x}_{ie},$$

(35)

where $\mathcal{E}=\{\mathrm{Si},\mathrm{Fe},\mathrm{Cu},\mathrm{Mn},\mathrm{Mg},\mathrm{Cr},\mathrm{Zn},\mathrm{V},\mathrm{Ti},\mathrm{Zr},\mathrm{Li},\mathrm{Ni},\mathrm{Be},\mathrm{Sc}\}$. A condition was considered composition-feasible when all reported elemental values were present and nonnegative, their sum did not exceed 100 wt.%, and the calculated Al balance was nonnegative. The calculated balance represents only the remainder after subtracting the reported elements; it does not constitute chemical certification, thermodynamic phase-stability verification, or experimental manufacturability validation.

Extrapolation was controlled using proper-training feature ranges and a distance-based applicability domain. The observed proper-training ranges used in the audit are summarized in Table 6.

A StandardScaler was fitted exclusively to the selected feature values of the 270 proper-training records. For each proper-training record, the mean Euclidean distance to its five nearest proper-training neighbours was calculated in the standardized selected-feature space using leave-one-out evaluation. The applicability-domain threshold was the 95th percentile of these proper-training distances:

$${\tau }_{\mathrm{AD}}=Q_{0.95}\left(d_{j,\mathrm{tr}}^{(5)}\right)=2.6849.$$

(36)

For held-out condition i, the five-nearest-neighbour distance was

$$d_i^{(5)}={\displaystyle \frac{1}{5}}{\displaystyle \sum _{\ell =1}^5}{∥{\mathbf{z}}_i-{\mathbf{z}}_{N_{\ell }(i)}∥}_2,$$

(37)

where ${\mathbf{z}}_i$ is the standardized selected-feature vector and $N_{\ell }(i)$ is its ℓth nearest proper-training neighbour. A condition was considered inside the applicability domain only when $d_i^{(5)}\le {\tau }_{\mathrm{AD}}$ and all selected feature values remained within their proper-training ranges. Recommendation eligibility required both composition feasibility and applicability-domain membership. Ranking ties were resolved deterministically using ascending source-record index.

4. Results and Discussion

4.1. Matched Full-Descriptor Comparison

The predictive effect of descriptor reduction was evaluated by comparing the 12-feature physics-retained pipeline with the matched 17-feature pipeline across the same 25 outer evaluations. Both pipelines used identical outer splits, inner-validation folds, eligible model families, hyperparameter ranges, and search budgets. Therefore, every reported difference was paired by repetition and outer fold. Values are reported as mean ± standard deviation across the 25 evaluations, and paired differences are defined as the 12-feature result minus the 17-feature result.

Table 7. Matched repeated-cross-validation performance of the 12- and 17-feature pipelines.

Metric	12 features	17 features	Paired difference	Corrected 95% CI	p
Macro-$R^2$	$0.8070\pm 0.0428$	$0.8002\pm 0.0474$	$+0.0068$	$[-0.0198,0.0335]$	0.602
Macro-NRMSE	$0.4243\pm 0.0472$	$0.4314\pm 0.0522$	$-0.0071$	$[-0.0339,0.0197]$	0.590
Macro-NMAE	$0.2926\pm 0.0335$	$0.2935\pm 0.0400$	$-0.0009$	$[-0.0223,0.0204]$	0.928

Table 7. Matched repeated-cross-validation performance of the 12- and 17-feature pipelines.

Metric	12 features	17 features	Paired difference	Corrected 95% CI	p
Macro-$R^2$	$0.8070\pm 0.0428$	$0.8002\pm 0.0474$	$+0.0068$	$[-0.0198,0.0335]$	0.602
Macro-NRMSE	$0.4243\pm 0.0472$	$0.4314\pm 0.0522$	$-0.0071$	$[-0.0339,0.0197]$	0.590
Macro-NMAE	$0.2926\pm 0.0335$	$0.2935\pm 0.0400$	$-0.0009$	$[-0.0223,0.0204]$	0.928

The 12-feature pipeline produced numerically higher macro-$R^2$ and lower normalized errors than the 17-feature comparator. The exploratory paired confidence intervals included zero, and the corresponding post-selection p-values exceeded 0.05. Because $k=12$ was selected from the same outer-validation results, these statistics should not be interpreted as an unbiased confirmatory hypothesis test. The results indicate similar observed mean performance while reducing the descriptor count from 17 to 12, corresponding to a 29.4% reduction; they establish neither predictive superiority nor formal statistical equivalence.

The practical advantage of the selected pipeline is consequently descriptor compactness and retention of the three metallurgically important heat-treatment variables, rather than a demonstrated accuracy advantage over the full-descriptor pipeline.

4.2. Comparison of Multi-Output Implementations

The four multi-output implementations were compared within the model-selection procedure. For the final fixed-$k=12$ search, each implementation was assigned the same budget of 24 randomly sampled configurations and evaluated using the same five inner-validation folds. Table 8 reports the best configuration obtained within each implementation. Inner NRMSE values are reported as mean ± standard deviation across the five folds, with lower values indicating better performance. Outer selections indicate how often each implementation was selected at $k=12$ across the 25 outer-training searches. Because these values were used for model selection, they are presented as inner-validation evidence and do not replace the outer-validation and untouched final-test results reported elsewhere.

The best family-specific configuration was identical in its principal tree and feature-selection settings across the four implementations: 800 trees, unrestricted maximum depth, $\max \text{\_}\mathrm{features}=0.7$, $\min \text{\_}\mathrm{samples}\text{\_}\mathrm{leaf}=1$, and $\alpha =0.25$. The native multi-output ExtraTrees implementation achieved the lowest final inner-validation NRMSE of 0.4074 and was selected in 16 of the 25 outer-training searches. It was therefore retained as the final model.

The independent target-wise wrapper produced a closely comparable inner-validation NRMSE of 0.4084. This confirms that the wrapper provided competitive predictions, but it should not be interpreted as modeling target dependence because its three regressors were fitted independently. The YS–UTS–El regressor chain achieved an inner-validation NRMSE of 0.4155 and was selected in two outer-training searches. Therefore, introducing sequential target predictions through this chain order did not improve the model-selection objective relative to native multi-output ExtraTrees for the present dataset. These results do not establish the absence of physical dependence among YS, UTS, and El; they show only that the evaluated chain formulation did not provide a predictive advantage under the present data and search protocol.

4.3. Effect and Stability of Physics-Retained Feature Selection

Feature-count performance was evaluated using five-fold outer cross-validation repeated five times within the proper-training partition. For every outer split and feature count, feature ranking, $\alpha$ selection, and model tuning were repeated using only the corresponding outer-training records. Table 9 presents the resulting performance as mean ± standard deviation across the 25 outer evaluations.

The highest observed mean macro-$R^2$ was obtained with 12 features. Its mean value of 0.8070 was numerically higher than the 0.8002 obtained using all 17 descriptors, while its mean macro-NRMSE was numerically lower. Accordingly, $k=12$ was selected using the maximum-observed-mean-$R^2$ rule. Because the same outer-validation results supported both feature-count selection and the subsequent 12-versus-17-feature comparison, the paired statistics in Table 7 are exploratory. The selection of $k=12$ should therefore be interpreted as choosing the highest observed mean together with a 29.4% descriptor reduction, not as evidence of superior or statistically equivalent predictive performance.

The standard error of the $k=12$ result was 0.0086, giving a one-standard-error threshold of 0.7985. Feature counts 8, 9, 10, 11, 12, 14, 15, 16, and 17 satisfied this secondary criterion, with $k=8$ being the smallest qualifying subset. Thus, $k=8$ represents a defensible parsimonious alternative, whereas the reported $k=12$ model follows the primary maximum-mean-$R^2$ selection rule.

Table 10 summarizes the search-based sensitivity to $\alpha$. For each value, the inner NRMSE is reported as mean ± standard deviation across five proper-training inner folds for the best evaluated configuration. The outer-selection frequency indicates how often each value was selected at $k=12$ across the 25 outer-training searches. Lower NRMSE indicates better performance.

The lowest final proper-training inner-validation NRMSE was obtained at $\alpha =0.25$, which motivated its use in the locked configuration. However, the differences among $\alpha =0.25$, 0.50, and 0.75 were small, and the outer-fold selections were distributed across these intermediate values. The results therefore support combining statistical relevance and embedded importance but do not indicate a sharply determined universal optimum. The selected value $\alpha =0.25$ should be interpreted as specific to the present proper-training data and search protocol.

The final 12-feature set contained $T_{\mathrm{sol}}$, $T_{\mathrm{age}}$, $t_{\mathrm{age}}$, Zn, Mg, Cu, Ti, Fe, Zr, Si, Sc, and Cr. This subset reduced the descriptor count by 29.4% while retaining all three metallurgically important heat-treatment variables.

Genuine feature-selection stability was evaluated across the 25 outer-fold subsets. For outer selections ${\mathcal{F}}_a$ and ${\mathcal{F}}_b$, pairwise Jaccard similarity was calculated as

$$J_{ab}={\displaystyle \frac{|{\mathcal{F}}_a\cap {\mathcal{F}}_b|}{|{\mathcal{F}}_a\cup {\mathcal{F}}_b|}}.$$

(38)

All $\left({\displaystyle \genfrac{}{}{0pt}{}{25}{2}}\right)=300$ outer-fold pairs were evaluated. Nogueira stability was also calculated from the binary selection matrix $\mathbf{Z}\in {\{0,1\}}^{M\times p}$:

$$\Phi =1-{\displaystyle \frac{\frac{1}{p}{\sum }_{j=1}^p{s}_j^2}{\frac{\overline{k}}{p}\left(1-\frac{\overline{k}}{p}\right)}},s_j^2={\displaystyle \frac{M}{M-1}}{\widehat{p}}_j\left(1-{\widehat{p}}_j\right),$$

(39)

where $M=25$, ${\widehat{p}}_j$ is the observed selection frequency of feature j, and $\overline{k}$ is the mean selected-subset size. The factor $M/(M-1)$ provides the finite-sample variance correction. The resulting Jaccard and Nogueira stability values are reported in Table 11.

The stability values remained high after excluding the three always-retained heat-treatment variables. Therefore, the observed stability did not arise solely from mandatory physics-guided retention; the data-driven non-mandatory feature selections were also reasonably consistent across the outer folds.

4.4. Sensitivity to Composition-Group Dependence

Table 12 compares two post-selection diagnostics using the same locked 12-feature ExtraTrees specification: repeated row-wise cross-validation and composition-group-disjoint validation. In both diagnostics, the selected features and hyperparameters were fixed before the splits were evaluated. Values are reported as mean ± standard deviation across 25 evaluations.

The row-wise diagnostic value differs slightly from the matched-pipeline result in Table 7 because the present diagnostic used one fixed feature set and hyperparameter specification, whereas the matched-pipeline evaluation repeated training-fold-contained feature ranking and model tuning.

The substantial performance reduction under composition-group-disjoint validation indicates that the random row-wise protocol benefited from similarities among conditions derived from the same or closely repeated alloy compositions. Therefore, the higher row-wise cross-validation and final-test results should be interpreted as interpolation performance within the represented composition domain and may be optimistic for previously unseen alloy compositions. The grouped result provides a more conservative estimate of composition-level generalization.

The grouped diagnostic used the already locked feature set and hyperparameters and therefore does not replace a fully nested composition-grouped model-selection experiment. A complete grouped rerun could select a different feature count, model family, or hyperparameter configuration. Nevertheless, the diagnostic directly demonstrates the effect of composition dependence and prevents the row-wise results from being interpreted as evidence of equivalent performance on compositionally novel alloys.

4.5. Target-Wise Prediction Performance

After evaluating the average performance across all targets, the prediction capability of the proposed framework was further analyzed separately for YS, UTS, and El. Target-wise analysis is important because the three mechanical properties do not have the same prediction difficulty. YS and UTS are strength-related properties and are generally more strongly coupled, whereas El represents ductility and is often more sensitive to microstructural variation, alloy processing history, and strength–ductility trade-off behavior.

Table 13 presents the point-prediction performance of the locked 12-feature ExtraTrees model on the untouched 91-record final-test partition. This evaluation was performed only after feature-count selection, model-family selection, and hyperparameter optimization were completed using the proper-training partition. MAE and RMSE are reported in MPa for YS and UTS and in percentage points for El, whereas NMAE and NRMSE are dimensionless.

The locked model achieved final-test $R^2$ values of 0.857, 0.842, and 0.826 for YS, UTS, and El, respectively. El had the largest normalized errors, with NMAE of 0.275 and NRMSE of 0.415, indicating greater prediction difficulty for ductility than for the two strength properties on the present final-test partition. Raw MAE and RMSE values were not averaged because the targets have different physical units and numerical scales.

Figure 2 shows the measured and predicted values obtained using the physics-retained ExtraTrees model. In each panel, the dashed line represents the ideal relation between measured and predicted values. The YS and UTS predictions are generally distributed near this line, indicating agreement within the present row-wise final-test partition. This alignment is consistent with the target-wise $R^2$ values reported in Table 13.

For El, the predicted values also follow the overall increasing trend of the measured values. However, a larger relative scatter is observed compared with YS and UTS, particularly at higher elongation values. This behavior is expected because elongation is more sensitive to local microstructural features, defect populations, precipitate distributions, and fracture behavior, which may not be fully captured by composition and heat-treatment variables alone. The higher relative error for El may therefore reflect both the greater sensitivity of ductility to unobserved microstructural variables and the limitations of the present composition–processing dataset.

Overall, the target-wise results and actual-versus-predicted plots indicate that the locked model captured substantial variation in all three mechanical properties within the represented row-wise final-test domain. The unified multi-target formulation permits the strength–ductility balance to be evaluated within the represented row-wise data domain because YS, UTS, and El can be considered simultaneously rather than as isolated properties.

4.6. Explainability Analysis

SHAP analysis was performed to interpret the predictions of the proposed physics-retained ExtraTrees model. Since the model predicts YS, UTS, and El simultaneously, target-wise SHAP analysis was conducted to identify the dominant descriptors for each mechanical property. The mean absolute SHAP value was used to quantify the overall magnitude of each feature’s model attribution, and target-wise bar plots were used to visualize the resulting global importance rankings.

Table 14 reports the five largest model attributions for each target. Zn received the largest mean absolute attribution for both strength properties. For El, the largest attributions were assigned to aging temperature, Zn, Cu, aging time, and Zr.

To assess the sampling stability of the global attribution rankings, 1000 bootstrap resamples of the 91-record evaluation subset were analyzed. The resulting target-wise Spearman rank correlations are reported in Table 15.

The bootstrap results indicate that the global SHAP rankings were stable with respect to resampling of the row-wise evaluation subset. This analysis evaluates sampling variability of the reported model explanations; it does not eliminate the dependence-related limitations caused by correlated input variables.

Figure 3 presents the target-wise global attribution rankings across the 91 final-test records. Longer bars indicate larger average absolute model attributions. These rankings explain the behavior of the locked model and should not be interpreted as causal material effects. Zn received the largest mean absolute attribution for both YS and UTS, while Cu and the heat-treatment descriptors also received substantial strength-related attributions. The attributions assigned to Zn, Mg, and Cu are qualitatively compatible with precipitation-strengthening mechanisms in heat-treatable aluminum-alloy systems, in which alloy chemistry influences solute availability and precipitation response. Aging temperature and duration regulate solute clustering, precipitate nucleation, growth, and coarsening and thereby influence the transition among under-aged, peak-aged, and over-aged conditions and the associated strength–ductility balance Girgis et al. (2023); Wang et al. (2007); Zhang et al. (2024). However, because the dataset combines different alloy chemistries and lacks direct descriptors of precipitate phase, morphology, coherency, number density, and spatial distribution, the SHAP patterns cannot identify a specific precipitation sequence or confirm a causal strengthening mechanism.

For El, aging temperature received the largest mean absolute attribution, followed by Zn and Cu. Higher $T_{\mathrm{age}}$ values were generally associated with negative SHAP values, whereas lower values more frequently produced positive contributions. This pattern is compatible with the established role of aging temperature and time in controlling precipitation kinetics and the transition among under-aged, peak-aged, and over-aged conditions Zhang et al. (2024). It may also reflect the conventional strength–ductility trade-off, but SHAP cannot distinguish the direct effect of aging temperature from correlated alloy chemistry, treatment duration, precipitate state, or unmeasured microstructural factors.

The three heat-treatment descriptors were retained a priori because of their established metallurgical relevance, not because of their data-driven rankings or subsequent SHAP values. Their appearance among influential descriptors therefore cannot independently validate or justify their mandatory retention. The SHAP analysis should be interpreted solely as a post hoc explanation of how the locked model used the prespecified inputs, rather than as causal metallurgical confirmation.

4.7. Uncertainty-Aware Reliability and Candidate Screening

After evaluating prediction accuracy and explainability, sample-specific ensemble-conformal intervals were evaluated on the 91 untouched final-test records. The intervals were generated at a nominal target-wise marginal coverage level of 90%. Table 16 reports the target-wise coverage, normalized conformal quantiles, and ensemble disagreement.

Table 17 compares the uncalibrated bootstrap-ensemble intervals with the ensemble-conformal intervals on the same 91 final-test observations. Coverage uncertainty is reported using 95% exact binomial confidence intervals. Mean interval widths are expressed in MPa for YS and UTS and in percentage points for El.

The uncalibrated bootstrap-ensemble intervals were substantially narrower but severely under-covered all three targets. Conformal calibration increased interval widths and improved marginal coverage from 49.45% to 89.01% for YS, from 51.65% to 80.22% for UTS, and from 57.14% to 91.21% for El. Simultaneous coverage increased from 29.67% to 73.63%. Thus, the comparison demonstrates that ensemble disagreement alone did not adequately represent total predictive error. Nevertheless, UTS remained below the nominal 90% level even after conformal calibration, indicating that the ensemble-conformal method should not be described as uniformly calibrated across all three targets.

The marginal coverage was close to the nominal level for YS and El, whereas the UTS interval showed under-coverage on the present final-test partition. Simultaneous coverage, requiring all three measured properties to fall within their respective intervals for the same test observation, was 67 of 91 observations, corresponding to 73.63%. This value was lower than the target-wise coverages because the intervals were calibrated separately for each property and were not designed to guarantee 90% simultaneous three-property coverage.

Figure 4 presents the measured values, ensemble mean predictions, and sample-specific conformal intervals for the 91 final-test observations. Black points represent measured values, while blue points and vertical error bars represent ensemble mean predictions and their 90% target-wise marginal prediction intervals. The observations were sorted by their measured target values only to improve visual readability and were not connected because they do not form a continuous sequence. YS and UTS are expressed in MPa, whereas El is expressed in percent. The pronounced variation in error-bar lengths demonstrates that the conformal interval width changes according to sample-specific ensemble disagreement.

The interval-width statistics are presented in Table 18. The substantial differences between the minimum and maximum widths confirm that the intervals are sample-specific. Wider intervals indicate greater disagreement among the bootstrap ensemble members and therefore a larger uncertainty penalty.

The 91 final-test records represented 87 unique composition–processing conditions. All 87 conditions passed the elemental feasibility check, while 81 conditions were within the proper-training applicability domain and were therefore eligible for recommendation. Six conditions were excluded from recommendation because their standardized five-nearest-neighbour distances or selected feature values were outside the proper-training domain. The measured final-test target values were not used for condition construction, scoring, feasibility assessment, applicability-domain estimation, sensitivity analysis, or ranking.

The ten highest-ranked eligible conditions involved three distinct compositions, denoted C1–C3 in Table 19. All omitted elemental values among Si, Fe, Cu, Mn, Mg, Cr, Zn, V, Ti, Zr, Li, Ni, Be, and Sc were zero. The elemental compositions and calculated Al balances are expressed in wt.%.

Table 20 reports the ranked heat-treatment conditions, predicted properties, conformal bounds, scores, and applicability-domain information. YS and UTS predictions and bounds are expressed in MPa, El predictions and bounds in percent, treatment temperatures in °C, and aging time in h. The heat-treatment variables retain their original source-data scale. Values are rounded for presentation, whereas ranking was performed using the unrounded values.

Here, $S_i^{\mathrm{prop}}=0.4S_i^{\mathrm{str}}+0.4S_i^{\mathrm{duct}}$, while ${\mathcal{S}}_i=S_i^{\mathrm{prop}}-0.2P_i^{\mathrm{unc}}$. All ten conditions in Table 20 were composition-feasible, remained within the audited proper-training feature ranges, and were inside the applicability domain with $d_i^{(5)}<2.6849$. HTC-026 was ranked first. Although its property-only rank was eighth, incorporating the sample-specific uncertainty penalty moved it to reliability-aware rank one, demonstrating that the sample-specific uncertainty term affected the ranking.

Sensitivity to the target weights and uncertainty penalty was evaluated using six prespecified scenarios without model refitting or conformal recalibration. The results are reported in Table 21.

HTC-026 remained the leading condition under equal property emphasis, YS-priority, UTS-priority, and high-uncertainty-penalty scenarios. Reducing the uncertainty penalty from 0.2 to 0.1 changed the leading condition to HTC-041 and reduced the top-ten overlap to 60%. Thus, the leading recommendation was stable under changes in the relative property emphasis but was sensitive to reducing the importance assigned to uncertainty. The ranking should therefore be interpreted as decision-weight dependent rather than as a unique universal optimum.

Overall, the analysis provides a transparent retrospective ranking of held-out conditions rather than evidence of newly generated or experimentally validated alloys. The composition-balance and applicability-domain checks reduce unsupported extrapolation, but they do not establish industrial manufacturability, thermodynamic stability, or performance outside the represented dataset. In addition, the UTS conformal intervals showed under-coverage on the present final-test partition; consequently, the ranked conditions should be treated as candidates for further experimental evaluation rather than direct design recommendations.

From an industrial perspective, the present workflow may be used as a preliminary prioritization stage for reducing the number of previously reported composition–processing conditions requiring detailed laboratory investigation. It does not replace alloy synthesis, thermodynamic or kinetic calculations, process-scale assessment, metallographic characterization, or mechanical testing. In particular, the applicability-domain check measures statistical proximity to the proper-training data and does not establish industrial scalability, forging suitability, quench sensitivity, heat-treatment reproducibility, cost effectiveness, or compliance with application-specific material standards. Predictions outside the audited composition and processing domain should not be used for recommendation without new data and external validation.

5. Conclusions

This study proposed a physics-retained, explainable, and uncertainty-aware multi-target machine-learning framework for predicting the strength–ductility response of aged and heat-treatable aluminum alloys. Using the aged forged aluminum alloy dataset reported by Chang, the prediction task was formulated as a multi-output regression problem to simultaneously estimate YS, UTS, and El from alloy composition and heat-treatment descriptors. The proposed physics-retained hybrid feature-selection strategy reduced the original input space from 17 variables to 12 variables while explicitly retaining the heat-treatment descriptors $T_{\mathrm{sol}}$, $T_{\mathrm{age}}$, and $t_{\mathrm{age}}$. The selected subset included $T_{\mathrm{sol}}$, $T_{\mathrm{age}}$, $t_{\mathrm{age}}$, Zn, Mg, Cu, Ti, Fe, Zr, Si, Sc, and Cr.

Across the 25 outer evaluations, the 12-feature pipeline achieved a macro-$R^2$ of $0.8070\pm 0.0428$, macro-NRMSE of $0.4243\pm 0.0472$, and macro-NMAE of $0.2926\pm 0.0335$. The corresponding 17-feature results were $0.8002\pm 0.0474$, $0.4314\pm 0.0522$, and $0.2935\pm 0.0400$, respectively. The exploratory post-selection paired comparison yielded a mean macro-$R^2$ difference of 0.0068, with a corrected 95% interval of $[-0.0198,0.0335]$ and $p=0.602$. Because $k=12$ was chosen from the same outer-validation summaries, these statistics are not treated as confirmatory evidence. The principal observed benefit of the selected pipeline was therefore descriptor compactness while retaining the three metallurgically important heat-treatment variables; superiority or formal statistical equivalence to the full-descriptor pipeline was not established. On the untouched final-test partition, the locked 12-feature ExtraTrees model achieved target-wise $R^2$ values of 0.857, 0.842, and 0.826 for YS, UTS, and El, respectively. SHAP analysis assigned the largest model attribution to Zn for both YS and UTS and to $T_{\mathrm{age}}$ for El, followed by Zn and Cu; these attributions describe model behavior and are not interpreted as causal metallurgical evidence. At a nominal target-wise coverage level of 90%, the ensemble-conformal analysis achieved marginal coverages of 89.01%, 80.22%, and 91.21% for YS, UTS, and El, respectively, while simultaneous three-property coverage was 73.63%. These results indicate near-nominal coverage for YS and El but under-coverage for UTS on the present internal final-test partition. Retrospective reliability-aware screening identified 81 of 87 unique held-out composition–processing conditions as composition-feasible and within the proper-training applicability domain, with HTC-026 ranked first under the baseline decision weights; however, these conditions require experimental verification and should not be interpreted as newly designed or externally validated alloys.

Overall, the proposed framework provides a compact and interpretable internally evaluated workflow for aluminum-alloy strength–ductility prediction. Its point predictions, conformal intervals, SHAP explanations, and retrospective screening results may support preliminary model-based prioritization within the represented composition domain, but they do not establish external validity, experimental performance, manufacturability, or practical alloy-design reliability.

An important limitation is that the source dataset contains alloy composition and heat-treatment conditions but does not provide direct microstructural descriptors such as grain size, precipitate phase, precipitate morphology, size distribution, number density, coherency, texture, porosity, or defect population. It also lacks detailed processing-history variables such as quenching rate, deformation history, forging schedule, specimen orientation, and testing protocol. These unobserved variables may explain part of the remaining prediction error, particularly for elongation, and limit the physical and causal interpretation of the learned composition–process relationships.

Future work should perform fully nested composition-group- or source-disjoint model selection and validate the locked workflow using external alloy datasets. Incorporating grain-scale and precipitate-scale descriptors, CALPHAD-derived variables, quenching and deformation history, and testing-condition information may improve both predictive performance and metallurgical interpretation. Additional uncertainty approaches, including Gaussian-process regression and quantile-regression models, should be evaluated across multiple external datasets. Most importantly, the highest-ranked conditions should be synthesized and experimentally tested before the workflow is considered for practical alloy-design or industrial decision support.