Longitudinal Hardness Prediction of a Weld Bead Using Linear and Nonlinear Mixed‐Effects Models and Artificial Neural Networks

1. Introduction

In the initial stage of electric arc welding, inherent to the process, there will always be a regime called the transient regime, characterized by the continuous variation of energy balance in the welded joint, consequently resulting in the variation of cooling rate values in the weld bead and heat-affected zone. After the transient regime, with enough time for energy dissipation, the quasi-stationary regime can be reached. If we assume that an observer is positioned at the heat source (electric arc), the quasi-stationary regime will be reached from the moment the observer does not notice variation in the distribution of temperatures around him [1], consequently, the temperature gradient as well as the cooling rates are constant throughout this regime.

In 2018, in a previous study by one of the present authors [2], the effect of varying cooling rates during the transient regime on the microstructure and on the mechanical property (hardness) along the weld bead was investigated, using steel ASTM A36 as base material and ER70S6 filler metal. Twelve specimens were made, each of the specimen being welded with different values of electric arc power ($p$), welding speed ($v$) and plate thickness ($d$), these being the process variables or factors of the design of experiments (DOE). In addition to these three DOE variables, the values of the carbon equivalent of the steel ($C$) were stored, using the formulation proposed by the International Institute of Welding (IIW), as well as the values of the relative thickness ($\tau$), a dimensionless variable used for classify the heat flow in the material, i.e., classify if the heat flow is three-dimensional, downward as well as laterally from the weld, or essentially lateral. The value of $\tau$ can be calculated from the values of $p$, $v$, $d$ and thermophysical properties of the welded metallic alloy. As a result of the study, it was found that the variation of Vickers hardness values ($H$) along the weld bead ($l$), during the transient regime, could be approximated to an asymptotic exponential function (AEF), being the asymptote of function referring to the hardness value of the quasi-stationary regime. Finally, twelve non-linear regressions were performed, one for each specimen, using the Levenberg-Marquardt method, obtaining the coefficients of each of the twelve AEF. As a main result, it was found that the hardness value at the beginning of the bead was approximately 20% above the quasi-stationary hardness value, thus being considered a sharp increase. It is noteworthy that, prior to the study by Cruz Neto [2], an increase in hardness at the beginning of the bead in this range of values had not been documented in the literature. A year later, in 2019, Arruda et al. [3] found the same trend, even comparing two different steels in the study, ASTM A36 and SAE 1045.

However, there was a common limitation to the aforementioned studies[2,3], the authors were able to adjust an AEF for each specimen, but they were unable to relate the AEF to the process variables. Consequently, it was not possible to evaluate the influence of each process variable on the variation of hardness values $H$ along the weld bead ($l$). In case at least one process variable has a significant influence on the values of $H$, and if there is a hierarchical structure in the data, multilevel modeling could be applied. Multilevel models allow the identification and analysis of individual heterogeneities and between groups to which these individuals belong, making it possible to specify random components at each level of analysis [4]. Therefore, this paper proposes to investigate the applicability of multilevel modeling to the longitudinal hardness of weld bead and, if applicable, compare the multilevel model obtained with an artificial neural network (ANN). It is worth mentioning that we can found a several welding studies applying ANN. For example, Li et al. [5] used ANN to predict the bead geometry with changing of welding speed, Hu et al. [6] realize a prediction of resistance spot welding quality, De Filippis et al. [7] apply ANN to predict Vickers microhardness and tensile strength of friction stir welding butt joints, and we can easily find other studies applying ANN for other purposes in the welding field. On the other hand, the use of multilevel modeling in welding studies is still very scarce, with more studies being found in the health context, studies that mainly sought to investigate the effects of welding fumes on the health of welders [8,9,10].

2. Materials and Methods

Aiming for a brief contextualization, Figure 1 presents a summary of the experimental methodology adopted by Cruz Neto [2], from the fixation and alignment of the specimens for the welding process to the adjustment of the AEF to the hardness profile, H versus l, for each specimen.

Table 1 presents the $p$, $v$, $d$, $C$, $\tau$, $l$ and $H$ values. The first column presents the identification of the specimen, with a total of 12. Each specimen is the result of a DOE testing, where the variables of process, consequently, in Table 1, the values of $p$, $v$, $d$, $C$ and $\tau$ appear repeated for the same specimen, varying only the values of $l$ and $H$ (hardness profile). The $p$, $v$, $d$, $C$ and $\tau$ variables are candidates for level in multilevel modeling, since, in principle, hardness values could be nested, or grouped, according to any of these variables. Initially, specimen was used as a level and, although it was possible to adjust multilevel models, having specimen as a level represented little physically. Finally, we chose to use $\tau$ as a level, since its value is associated with the heat flow in the material, having a physical relationship with the cooling rate and hardness values. The data set has a total of 871 observations, a small part of which is presented in Table 1. To adjust the models, both the multilevel models and the ANN models, a randomly generated sample containing 70% of the observations was used, the remaining 30% was separated to evaluate and compare the different models obtained.

Before fitting the multilevel models, a single nonlinear regression model (NLS, nonlinear squares) - model_1 was fitted, as shown in Figure 2, using the Gauss-Newton method, the standard method when using the package “nls”. model_1 will serve as a reference for evaluating the adjusted multilevel models later.

After fitting a single non-linear model, the applicability of multilevel modeling was investigated. To adjust the models, the package “nlme” of the R language was used. First, a linear multilevel model was adjusted, using the function lme, requiring the linearization of an asymptotic exponential model for this purpose, Equation (1).

$$X=\mathrm{e}\mathrm{x}\mathrm{p}(-c3\times l)$$

(1)

where $l$ is the length of the weld bead, it can be noted that for the linearization the $c3$ coefficient obtained in the adjustment of model_1 was used. Therefore, we have the linearized multilevel model (LME, linear mixed-effects) as shown in Equation (2):

$$H_{ij}={\beta }_{0j}+{\beta }_{1j}{\mathrm{X}}_{ij}+{\epsilon }_{ij}$$

(2.1)

$${\beta }_{0j}={\gamma }_{00}+{\gamma }_{01}{\mathrm{C}}_j+{\nu }_{0j}$$

(2.2)

$${\beta }_{1j}={\gamma }_{10}+{\gamma }_{11}{\mathrm{C}}_j+{\nu }_{1j}$$

(2.3)

where ${\beta }_{0j}$ and ${\beta }_{1j}$ are the LME coefficients and ${\epsilon }_{ij}$ the respective error terms, with observations i nested in j groups. It is worth nothing that, in this study, the observations will be grouped according to the values of $\tau$, as described in the previous section. The coefficients ${\gamma }_{00}$ and ${\gamma }_{10}$ represent the expected value of H for a given observation i belonging to group j, respectively, when $X=C=0$ (general intercept) and where there is a unit change in the value of $X$, ceteris paribus. And ${\gamma }_{01}$ e ${\gamma }_{11}$ represent the change in the expected value of $H$ for a given observation i belonging to group j when there is unit variation, respectively, in the value of $C$ and in the value of the product $C\times X$, ceteris paribus. In addition, ${\nu }_{0j}$ and ${\nu }_{1j}$ represent the error terms arising from the presence of observations from different groups in the dataset and indicate, respectively, the existence of randomness in the intercepts and the existence of randomness in the slopes of the models referring to the groups [4]. Combining Equations (2.1)–(2.3), we arrive at Equation (3):

$$H_{ij}=\underset{randomintercepts}{\underbrace{\left({\gamma }_{00}+{\gamma }_{01}{C}_j+{\nu }_{0j}\right)}}+\underset{randomslopes}{\underbrace{\left({\gamma }_{10}+{\gamma }_{11}{C}_j+{\nu }_{1j}\right)}}{X}_{ij}+{\epsilon }_{ij}$$

(3)

Rearranging the terms of Equation (3), the components of fixed effects and components of random effects are evident, as shown in Equation (4):

$$H_{ij}=\underset{fixedeffects}{\underbrace{{\gamma }_{00}+{\gamma }_{01}{C}_j+{\gamma }_{10}{X}_{ij}+{\gamma }_{11}{C}_j{X}_{ij}}}+\underset{randomeffects}{\underbrace{{\nu }_{0j}+{\nu }_{1j}{X}_{ij}}}+{\epsilon }_{ij}$$

(4)

Analyzing Equation (4), if the statistical significance of the variances of the error terms ${\nu }_{0j}$ and ${\nu }_{1j}$ is not verified, that is, if both were statistically equal to zero, the model would be reduced to a multiple linear regression model, as shown in Equation (5):

$$H_{ij}={\gamma }_{00}+{\gamma }_{01}{\mathrm{C}}_j+{\gamma }_{10}{X}_{ij}+{\gamma }_{11}{X}_{ij}{\mathrm{C}}_j+{\epsilon }_{ij}$$

(5)

In this scenario of the non-significance of the error terms ${\nu }_{0j}$ and ${\nu }_{1j}$, Equation (5) could be estimated by traditional methods, such as Ordinary Least Squares (OLS). Therefore, we see the importance of running a test to determine whether the variance of the error terms ${\nu }_{0j}$ and ${\nu }_{1j}$ are significantly different from zero, and the results of these tests are presented in the next section. Figure 3 shows how the LME named model_2 was adjusted in the R script.

After adjusting the LME model, a nonlinear multilevel model (NLME, nonlinear mixed-effects) was adjusted, called model_3 in the R script, using the function nlme, as shown in Figure 4.

Finally, after obtaining the LME and NLME models, several ANN models were trained, always having as input neurons $l$, $C$ and $\tau$ and as output neuron the hardness of the weld bead $H$, as can be seen in Figure 5 an architectural example. However, it is worth mentioning that several architectures were tested, modifying the number of neurons per hidden layer, number of hidden layers, type of activation function and the optimization method. For the training of the ANN, the packages “neuralnet” and “h2o” were used, and the commands in R language can be seen in Figure 6.

3. Results and Discussion

3.1. Non-Linear Regression Model versus Multilevel Models

Figure 7 presents the experimental values of hardness H versus bead length l for all specimens obtained in DOE and the values predicted by the NLS model. After visual analysis, it can be noted that the residuals are quite high and that there is heteroscedasticity.

Figure 8 shows the result of the variance significance test for the error terms ${\nu }_{0j}$ and ${\nu }_{1j}$. In the last column, one can visualize the p-values, noting that all values were significant at a significance level of 0.05. Therefore, it can be concluded that the variance values of the error terms ${\nu }_{0j}$ and ${\nu }_{1j}$ are statistically different from zero, with randomness in the intercepts and randomness in the slopes.

From the variance values, it was possible to calculate the intraclass correlation through $var\left({\nu }_{0j}\right)/(var\left({\nu }_{0j}\right)+var\left({\nu }_{1j}\right)+var(e\left)\right)$, obtaining that the variable $\tau$ is responsible for, approximately, 48% of the $H$ hardness variation. This result is somewhat predictable, due to the nature of the relationship between $\tau$ and $H$ being well known. However, in none of the cited studies [2,3], it was possible to isolate the effect of the variable $\tau$ on the values of $H$, which proved to be a relevant contribution of multilevel modeling.

In Figure 9, one can see how the LME model, compared to the NLS model, fits better to the data by considering the hierarchical structure of the data. Figure 9a shows only the experimental values of $H$ and Figure 9b the same values, but this time, discriminating the values of each specimen in colors. Figure 9c shows the values predicted by the LME for each level and Figure 9d only the predicted values.

Figure 10 shows the values of each of the terms ${\nu }_{0j}$ and ${\nu }_{1j}$ and how they vary according to the values of τ. It is worth mentioning that the LME is a linearized model and, therefore, the intercept terms ${\nu }_{0j}$, have the effect of translating the FEA up and down, impacting the values of maximum hardness at the beginning of the bead and hardness in the quasi-stationary regime. In Figure 10a, it can be observed that as the value of $\tau$ increases, a trend towards an increase in the values of ${\nu }_{0j}$ is visualized. This result presents physical coherence, since the higher the value of $\tau$, the higher the cooling rate values, which results in the formation of phases of greater hardness for the alloy under study, such as, for example, the greater formation of acicular ferrite. In relation to the terms of ${\nu }_{1j}$, greater difficulty was found in the interpretation, since it was not possible to visualize a clear trend of variation in the terms of ${\nu }_{1j}$ as the values of $\tau$ change.

Figure 11 presents the graph of predicted values versus experimental hardness, the desired thing is that the results are as close as possible to the diagonal line, in black, which indicates equivalence between the predicted and experimental values. In this graph, the results obtained by the NLS, LME and NLME models were compared. It can be observed that the NLS model obtained the worst performance and the LME and NLME models presented predicted values closer to the experimental ones.

Among the metrics presented, one of the most relevant is the log-likelihood (LogLik) value, Table 2. It is classified as the best model, the one that presents the highest LogLik value, therefore, by this metric, the best model is the LME, with the value of the NLME being quite close and, again, the NLS is classified as the worst model. However, LogLik is not usually used as a metric for ANN, the next model that will be compared, for this reason, the Sum Squared Error (SSE) and Mean Squared Error (MSE) metrics were also added in Table 2. Since SSE and MSE are error terms, the lowest possible values are desired, with the lowest value obtained for the NLME, with the LME having a very close value. It is worth noting that, unlike the LogLik values, the SSE and MSE values were calculated from the test data set (sample of 30% of the original data).

3.2. Multilevel Models versus Artificial Neural Networks

Considering that the metrics of the LME and NLME models were close, a significance test was performed, based on the LogLik values, to assess whether the models are significantly different. For this purpose, the function lrtest from the package “lmtest” was used, the results are shown in Figure 12. Adopting a 95% confidence interval for the test, since the p-value was below the significance level, it is assumed that the LME and NLME models are statistically different.

Next, an ANN was trained from the training dataset, using the package “neuralnet”. After some attempts, we finally arrived at an ANN, called ANN1, with three hidden layers containing, respectively, 8, 5 and 3 neurons, and a sigmoid logistic activation function. Figure 13 shows ANN1 after training. The thicker the line that connects two neurons, the greater, in module, is the value of the weight associated with this connection; in black, the positive weights are represented and, in gray, the negative weights.

Figure 14 shows the Real vs. Predicted for the LME, NLME and ANN1 models, with a worse performance of ANN1 compared to the multilevel models.

Table 3 presents the LME, NLME and ANN1 metrics, observing the lowest SSE and MSE values for ANN1, confirming the result of the visual analysis of Figure 14.

Due to the results obtained with ANN1, a new ANN was trained, called ANN2, with the same architecture as ANN1, only modifying the activation function to hyperbolic tangent and number of epochs to 50000. To train ANN2, the package “h2o” was used. The same training set as the other models was used. However, for cross-validation, the function h2o.deeplearning also requires a validation dataset, obtained via random sampling, containing 30% of the original training set. Consequently, the new training set contains 70% of the remaining data from the original training set. Figure 15 presents the predicted versus experimental values of $H$, where a better fit of the ANN2 in relation to the multilevel models is observed. Finally, for the calculation of the SSE and MSE of ANN2 (Table 4), the same test set already used for the other models was used. After several rounds of training, finally, an ANN with better metrics was obtained compared to the other models.