Modeling the Bending Strength of MDF Faced, Polyurethane Foam-Cored Sandwich Panels Using Response Surface Methodology (RSM) and Artiﬁcial Neural Network (ANN)

: The present study evaluates and compares predictions on the performance and the ap-proaches of the response surface methodology (RSM) and the artiﬁcial neural network (ANN) so to model the bending strength of the polyurethane foam-cored sandwich panel. The effect of the independent variables (formaldehyde to urea molar ratio (MR), sandwich panel thickness (PT) and the oxidized protein to melamine-urea-formaldehyde synthesized resin weight ratio (WR)) was examined based on the bending strength by the central composite design of the RSM and the multilayer perceptron of the ANN. The models were statistically compared based on the training and validation data sets via the determination coefﬁcient ( R 2 ), the root mean squares error (RMSE), the absolute average deviation (AAD) and the mean absolute percentage error (MAPE). The R 2 calculated for the ANN and the RSM models was 0.9969 and 0.9960, respectively. The models offered good predictions; however, the ANN model was more precise than the RSM model, thus proving that the ANN and the RSM models are valuable instruments to model and optimize the bending properties of the sandwich panel.


Introduction
A typical sandwich panel is composed of face sheets and a lightweight core with adhesives used to connect them. The common patterns of the sandwich cores are foam, truss, honeycombs or corrugated patterns. Building foams with different components, such as foams with open and closed cells, can be used as the core of a sandwich panel [1,2]. The core layer offers a high load-bearing capacity coupled by a rather low weight [3] so that in many applications the structure's stiffness is very important. It has been observed that the increase in the core thickness by 4 times can significantly affect the stiffness and the bending strength by up to 37 times, while the increase in the panel's weight can be slight and up to 6% [4]. All transverse forces that cause normal stresses in the core layer are frequently small, and a slight decrease in the core thickness should result in a further reduction in the flexural stiffness. However, the Young's modulus perpendicular to the forces should also be at a high level. Generally, the core is mainly exposed to the shear, and core shear strains cause both global deformation and core shear stress. In this case, the increase in the core thickness reduces the critical force and maximum shear stress [3].
In a sandwich panel, the face sheets must be connected well to the core so that during the load transfer between the components the compressive load is borne by one sheet and the tensile loads are borne by another sheet. However, the shear loads are stabilized by the core [5]. Therefore, the strong connection between the core and the face sheets is very important so that the structure can successfully tolerate the tensile and shear stresses [6]. However, due to construction defects or loading in the service, this connection can result in a regional debonding between the face sheets and the core, resulting in face wrinkling, buckling and shear core failure [7]. Therefore, the sandwich panel properties depend not only on the thickness but on the glue line strength between the sheets and the core.
Many types of paper and wood composite-based materials, such as medium density fiber board, can be used as the face sheet in the production of sandwich panels [8]. These natural fiber productions consume less energy than others in sensitized fiber production and result in both lower cost and reduced air pollution [9]. Moreover, with respect to the high performance of the MDF and other wood fiber and particle-based composites in the production of sandwich panels used as flooring, roofing or panelized walls, one can expect to save building time and labor costs while also contributing to the acceleration of forest restoration efforts [10]. On the other hand, in order to improve the mechanical properties and the efficiency of the MDF, researchers have shown that using glass fiber reinforced laminates as a face sheet can be used to produce MDF-cored sandwich panels [11]. Hussain et al. [12] have studied the application of hollow MDFs as a core layer in the production of sandwich panels with face sheets made of glass fiber reinforced polymer and showed that the MDF honeycomb structures were superior due to their capacity to sustain against the load after failure, although the mechanical characteristics of the solid MDF were similar when compared to the MDF honeycomb panels. Hassan et al. [13] have evaluated the separation strength between the glass fiber/epoxy laminates as face sheets with the MDF core layer and showed that the strength at the interfaced layer was stronger than the debonding of the core layer based on UF resin. In the production of the sandwich panel and joining of the face sheets to the core layer, the application of a safe adhesive with a good quality to create a good bonding line is necessary. With respect to the renewability of the eco-friendly natural materials, it is reasonable to use them as the natural adhesive components due to their safety when compared to the synthetic resins, as they sometimes have better properties such as higher resistance against moisture absorption and other advantages. One of these materials is an adhesive based on plant proteins where there is not only a higher bonding strength but a higher resistance against moisture absorption for external applications, which is mainly due to chemical modification. Many researchers have conducted experiments to improve the bonding strength and the resistance against protein water absorption through various physical [14], chemical [15] and biological [16] treatment methods. Moreover, chemical modifications using different methods such as acidic and alkaline treatment [17], graft copolymerization [18] and crosslinking [19] have been the most effective and most common methods.
The most recent research has focused on the prediction of the strength of composites using modeling instruments such as the response surface methodology (RSM) and the artificial neural network (ANN) when compared to the experimental approaches. The RSM is a statistical technique to design factorial experiments used to make mathematical models describing the effects of several quantitative factors on one or more responses [20], and to solve the multivariate equations simultaneously so to offer an optimal solution [21]. The RSM significantly reduces the number of experiments required to evaluate multiple parameters and their interactions so to generate more time and facility [21]. This method was used in the modeling of the physical and mechanical properties of different wood composites such as particleboard, fiberboard and glued laminated products. On the other hand, the ANN is a computational technique imitating human brain capacity. This modeling and simulation can handle very complicated complexes and nonlinear systems. When a data set is well-trained by a neural network, it can be used to predict the output of some new input data sets. This technique was not used to model the sandwich panel properties to predict and determine the optimal point of properties; however, various studies were performed on other wood-based panels and showed that the developed ANN models could evaluate the properties of the wood composite products with a good precision [22][23][24].
Previous studies have tried to use the RSM and the ANN methods separately to predict and optimize properties such as the strength of wood products. However, these methods were not studied simultaneously and their performance was not compared so to predict the current properties of the sandwich panels. Therefore, the present study attempts to evaluate, predict and optimize the bending strength of the sandwich panels by changing the formaldehyde to melamine-urea molar ratio (MR) used in the bonding line, in the sandwich panel thickness (PT) and in the modified protein to MUF resin weight ratio (WR) used in the bonding line by comparing the outputs of the RSM and the ANN methods.

Materials and Methods
First, the modified adhesive protein from soybean flour was prepared. According to the experimental set up being used, the MUF resin was made separately with different formaldehyde to melamine-urea MRs. Polyurethane foam manually made in the laboratory was also used as the core. The face sheets were also chosen from the MDF with a thickness of 2.7 mm.

Preparing the Protein Adhesive
After buying edible soybean flakes from the market (containing 1.1 g fat, 0.21 g sodium, 3.13 g sugar, 7.6 g carbohydrate and 53.3 g protein in 100 g) and grinding them, the resulting flour was sifted by a 100-mesh size sieve to be used as an ingredient to make the protein adhesive. A combination of 259.333 g water, 9.33 g NaOH 99% and 1.766 g ethylene glycol (phase transferer) was heated after being loaded in a three-necked flask equipped with a mechanical mixer, thermometer and condenser to 70 • C. After adding the flour (116.66 g) to the diluted environment of the flask while the mechanical mixer was working at a speed of 650 rpm (to achieve a suitable distribution of the flour in the diluted environment), the mixture's temperature increased to 88-90 • C after 15 min. After maintaining these conditions for 2 h, and when the reaction was complete, the flask's temperature decreased to 35 • C in the ice bath. The obtained solution with a pH of 12.9, after being neutralized by HCl 20%, was sifted by a 35-mesh size sieve to remove the lumps (with the portion of nearly 1:20) which likely formed due to the new hydrogen bonds during the alkali treatment. Finally, the resulting solution was kept in the refrigerator (+3 • C).

Making Melamine-Urea-Formaldehyde (MUF) Resin
Powder melamine with a purity of 99.8%, urea with a purity of 76% (Khorasan Petrochemical Co., Bojnourd, Iran) and formalin with a density of 1.08 gr/cm 3 , a pH of 2.5-4 and a concentration of 37.5% (Dr. Mojallali Industrial Chemical Complex Co., Tehran, Iran) were used. According to the experimental design, the MUF resin was made with three different F:MU MRs (including 1.68:1, 1.805:1 and 1.93:1) using the three-necked flask equipped with the condenser, thermometer and magnetic mixer. First, the urea (79.48 g ≈ 0.93 mol) and the formalin (134 g ≈ 1.68 mol, 144 g ≈ 1.805 mol or 154 g ≈ 1.93 mol formaldehyde) were added to the flask, and NaOH (40%) was added in drips, and then the solution's pH (5.3-5.5) increased to 8-8.4. As the reaction's temperature increased to 55-60 • C in this pH, and kept for 30 min, the temperature increased to 60-65 • C. Then, butanol (2.50 g) and the whole melamine (9.08 g ≈ 0.07 mol) were immediately added to the solution. After maintaining these conditions for 30 min, the mixture's pH reached 5-5.5 using the ammonium chloride solution of 20%, while the reaction environment temperature also increased to 80-85 • C. At this temperature (80-85 • C), and after 5-10 min, the solubility of the solution with the water was tested. After ensuring complete dissolution of the mixture with water at the WR of 1:1.5, the mixture was cooled down to 60 • C and the remaining urea (6.9 g) was loaded. The heater was then turned off and the solution cooled down. At the end of the synthesis process, the reacted compound (with pH ≈ 7) was mixed with ammonia (a formation of an aqueous media containing about 20% ammonia) so to improve the obtained resin stability. Hence, three types of MUF resin were prepared with different molar ratios. To prepare the final adhesive, the protein adhesive to MUF resin WRs were prepared as the third factor at ratios of 30:70, 50:50 and 70:30, respectively. To mix the protein adhesive with the MUF resin effectively, a mechanical mixer with a high speed was used for 5 min.

Making the Core Foam
To make the polyurethane (PU) foam, polyol (Table 1) and isocyanate ( Table 2) were used (produced by Mokarrar Industrial Group Co., Tehran, Iran) with the WR of 1:1 combined with plant fibers both to improve the core's bending strength and to act as 5% of the total foam produced. After being loaded into a mechanical mixer at a speed of 750 rpm to achieve a homogeneous system, the mixture was mixed for 5 min. After putting the mixture into a detachable mold with the dimensions 300 × 100 × 100 mm, enough time was given to the complex to create foam. After treating for 72 h at room temperature and removing from the mold, the foams were cut with a circular saw to the thicknesses of 10 mm, 20 mm and 30 mm (according to one of the independent variables being examined), respectively. The optimal point of applicability of the fibers in the system was determined by trial and error so to distribute the fibers properly in the foam solution environment. The recycled fibers were also obtained from the MDF being used after being immersed into water for 72 h and being defibrated by a hammer mill. After drying, the fibers were sifted by a 20-mesh size sieve and then a 70-mesh size sieve was chosen so to be added to the foam. The foam density was fixed at 0.07 g/cm 3 . Table 1. Characteristics of polyol (plyMbk-327b2).

Appearance Characterization Thick and Yellowish
Viscosity at 25 • C 1200-1600 mPa.s Specific gravity 1.1 g/cm 3 Consumption time 6 months after producing

Making the Sandwich Panel
After applying the adhesive on the measured MDF sheet with the density of 0.75 g/cm 3 and the thickness of 2.7 mm, the PU foam of this thickness was put between two face sheets according to the test plan being used. The level of the resin used was 180 g/m 2 (based on the dry matter) and the concentration of the adhesive solution was 70%. After being put under the press, the complex was exposed to a temperature of 160 • C and 2 kg/cm 2 pressure for 5 min. After removing from the press and trimming with the width of 90 mm, the bending strength test was performed to determine the modulus of rupture (MOR, based on the Equation (1)) of the specimens by a mechanical test device (Load Cell 2-Ton, Sanaf Co. Ltd. Tehran, Iran) with a loading speed of 2 mm/s. The density of the panels with the core thicknesses of 10 mm, 20 mm and 30 mm was 0.31 g/cm 3 , 0.27 g/cm 3 and 20 g/cm 2 , respectively. where σ f is the bending strength (MPa), F is the load (N), L is the span length (mm), h is the specimen's thickness measured along with the force applied (mm) and b is the specimen's width (mm).

Response Surface Methodology
A three-level and three-factor center composite design (CCD) of the RSM was created to model the bending strength of the sandwich panel by combining the conditions of 20 treatments with 2 iterations and 34 specimens in total. This design was created to achieve a second-degree polynomial model (Equation (2)) describing the increase in the bending strength of the sandwich panel (MOR) as a function of three independent variables, including the formaldehyde to melamine urea molar ratio (MR), the sandwich panel thickness (PT) and the protein adhesive to MUF resin weight ratio (WR), respectively. In fact, this design makes it possible to ensure the prediction correctness when testing for more extensive conditions in which the model's complexity is not known by establishing different levels for any variable being examined [25]. When combining the conditions of three independent variables with three levels, three types of points can be identified in the matrix cube produced by this method: 8 factorial points, 6 axial points and 6 center points assuming 2 iterations for any factorial point and star. All experiments were performed randomly, and the data were analyzed using the Design-Expert software (Version 13.0, Stat-Ease Inc., Minneapolis, MN, USA) where β 0 , β i , β ii and β ij are the model's regression coefficients, while X i X j and Y are independent and dependent variables, respectively. The parameters acting on the response, their levels and their combination of conditions are given in Table 3.

ANN Modeling
An artificial neural network (ANN) was used to predict the nonlinear relationship between the input parameters (X 1 , X 2 , X 3 ) and the dependent variable (Y). A back-propagation neural network together with the Levenberg-Marquardt algorithm were used to model the MOR of the sandwich panel as a function of the independent variables being examined. The accuracy was obtained by modifying the number of layers and the neurons in the different layers of the ANN architecture. The three-layered underlying architecture of the ANN model includes 1 input layer with 3 neurons (X 1 , X 2 , X 3 ) and 1 hidden layer composed of 6 neurons and 1 output layer (MOR) ( Figure 1). The hit and trial method was used to ensure that the number of neurons needed in the range, from 1 to 15 in the hidden layer, so to minimize any deviation between the predicted and the experimental results.  The experimental data set used to develop the RSM model was employed to make an ANN model in which 70% (24 points) was used to train the network, 15% (5 points) was used for validation and the remaining 15% (5 points) was used to test the network.
To determine a suitable limit of the percentages stated, trial and error was used based on different statistics such the determination coefficient (R 2 ), the root mean square error (RMSE), absolute average deviation (AAD) and the mean absolute percentage error (MAPE). The output signal was produced by passing the weighted sum of the input variables to any neuron through an active function that is usually nonlinear (sigmoid, hyperbolic tangent, Gaussian, linear, bipolar linear and threshold linear function) and is given by the hidden layer in the ANN architecture (Equation (3)).
where Y t is the network output, n is the number of hidden nodes, m is the number of input nodes, f is the transfer function, β ij {i = 1, 2, . . . , m; j = 0, 1, . . . , n} are the weights resulting from the hidden nodes, α j {j = 0, 1, . . . , 0} are the vectors of the weights from the hidden nodes toward the output nodes and α 0 and β 0j are the arc weights resulting from the bias terms that are always equal to 1.
To clarify the performance of the ANN models, different statistics were used based on the R 2 , the RMSE, the AAD and the MAPE (Equations (4)-(7), respectively): where n is the number of experiments, i is the experiment's number and Y ei and Y pi denote the experimental and predicted MOR values, respectively.

Characterization
From the MUF resin (with a MR equal to 1.93:1) and the MUF resin containing the modified protein adhesive (with WRs of 70:30 and 30:70), samples were selected for further characterization. The chemical structures of these samples were separately investigated using the Fourier transform infrared spectrometer (FTIR) (Nexus 470, Thermo Nicolet, Madison, WI, USA)). The spectra were obtained in the ATR mode at room temperature with 20 scans from 4000 to 400 cm −1 at the resolution of 8 cm −1 . Each sample was tightly fixed to the ATR accessory in order to ensure sufficient contact of the infrared light with the surface of the sample. A scanning electron microscope (SEM) (Hitachi, Japan) was used to observe the morphology of the samples as well.

FTIR Analysis
The FTIR analysis of the MUF resin (with a MR equal to 1.93:1) and the MUF resin containing the modified protein adhesive (with the WRs of 70:30 and 30:70) is given in Figure 2. Comparing the spectra of the control MUF and the modified MUF, there are also differences despite their similarities. It is observed that a strong absorption peak at 3310-3340 cm −1 belongs to the -OH and amide-groups (type 1) and NH 2 -group resulting from a stretching vibration band of the N-H functional group and the hydrogen bond [26] between the carbonyl groups of the peptide linkage in the protein and the wood surface. The larger intensity and the area under the curve in this band for the specimens with more protein indicates the exposure of the modified protein structures to the MUF resin. The presence of the symmetrical stretching vibration -C-H bands and the stretching -O-CH 3 bands is confirmed by a strong absorption peak in the range 2900-2990 cm −1 [27]. The shift and the change in the intensity of the stretching of C-H peaks in the composite confirm some chemical interactions in the matrix, so that as more modified protein is added, and even when the MR decreases, the B-amide absorption peak decreases gradually in this range, thereby demonstrating that the reactive groups are more exposed and have better access to the reactive units in the MUF resin so to create a chemical bond. Subsequently, the shift and the change in the stretching C-H peak indicates that the -CH 2 groups play a significant role in the connection of protein with resin, and the elimination of the peak in the curve (C) shows the high intensity of this interaction. Since the relative intensity of the bands is proportional to the level of the groups in the protein polymer [28], after the normalization of the absorption peak related to CH-in the range of 2940 cm −1 , the comparison of the relative number of reactive groups through the absorption intensity of the related IR peaks in 1600 and 2400 cm −1 indicated that their absorption intensity changed completely, thus showing that the relative number of the reactive groups to create a chemical bond in the adhesive mixture containing the maximum modified protein strongly increased even when the formaldehyde to melamineurea MR is less. A strong peak in the range of 2400 cm −1 corresponding to the stretching vibration of the bridged -CH 2 group is deleted, thereby strongly demonstrating the presence of the methylene bridge formation [29]. Indeed, the results show that the free amino protein groups were able to react with formaldehyde and penetrate into the coagulated UF connecting structure. At the same time, the absorption peaks of the bending vibration band of the N-H functional group were also observed in the ranges 1550 cm −1 and 1380 cm −1 . In the composition containing the modified protein, there are typical infrared absorption bands in the ranges of 1630-1680 cm −1 (C = O stretching), 1530-1559 cm −1 (for N-H bands) and 1260-1420 cm −1 (for C-N bending bands) for amides of the types I, II and III, respectively [30,31]. The stretching vibration mode of the modified proteins at 1630 cm −1 belongs to the C = O bands in amide I, while the bending vibration mode of the modified protein at 1530 cm −1 belongs to the N-H group for amide II and 1380 cm −1 and 1036 cm −1 belong to the C-N stretching of CH 2 -N and C-N stretching of the methylene linkage (NCH 2 N), respectively. N-H and C = O in peptide bands formed an initial protein backbone [32]. For the MUF/MP composition, amide bands I, II and III form compared to the MUF spectra in which these bands do not exist, and are higher compared to the lower MP values and the higher formaldehyde to melamine-urea MR. Moreover, the change in the relative intensity between amide I and II can be observed in curves B and C, and the relative intensity between 1380 cm −1 and the amide III band can be also observed for the specimen containing more MP and the lower formaldehyde to melamine-urea MR, indicating the reaction between the MUF resin and the MP adhesive. In addition, these bands may be accompanied by the reaction between protein -NH 2 and the C = O groups of the MUF resin to form the -CH 2 -NH-CH 2 -OH bridges [33].

SEM Analysis
The comparison of the MDF as the face sheet interior surface of the broken panel after the test (Figure 3) indicates that protein application could create a more complete smearing surface so to improve the connection and the contact surface, such that the protein content increases to 70%, as well as the coating which also increases ( Figure 3B in comparison to Figure 3C). The cross section of the connection surface also shows that, as the modified protein increases, the interphase compression increases between the face sheets and the core ( Figure 3B-c in comparison with Figure 3C-c). Moreover, the better stress transfer and the distribution on a larger area reduces the stress concentration; thus, the bonding interphase layer shows a higher strength during the force exertion. According to the porous surface of the face sheets, the existence of more protein results in a less porous surface followed by a higher contact surface. With less protein application, the viscosity cannot be changed and can even decrease due to the excessive breaking of the polymer chains under the effect of the alkaline treatments. However, as the consumption of the protein adhesives dramatically increases (as a result of the preference of the agglomeration effect due to the easier access of the water molecules to the coarser and bulkier functional hydrophilic groups, such as aldehyde and carboxyl groups, compared to the smaller hydroxyl groups due to the depolymerization), the viscosity also increases. Therefore, less resin penetrates the core layer with a very high porosity ( Figure 3A) and the interior surface of the face sheets with a high nonuniformity ( Figure 3D) so that the glue line thickness behaves in a way that a more uniform glue line is formed. Thus, the starve joints cannot develop, which causes stress concentration (i.e., the regions where microcracks start).  It is well-known that the smoothness of the surface of the MDF is very important. Nevertheless, the surface of the MDF can have microporosity. Adding the modified protein to the glue line improves the intimated contacts of the adhesive with the surface fibers. It was determined that applying modified protein can provide increments of properties of the MDF due to the formation of straightened fiber and its ordered distribution after hot press processing, so that fibers were composed closely and coated fully by the adhesive [34]. Modified protein is a hydrophilic compound and can absorb a higher amount of water. As water acts as plasticizer, the mobility of the soybean protein polypeptide chains is improved so that it can interact easily with other polymers. According to the findings based on SEM analyses, due to the hydrophilicity of proteins, after the hot press the shape of the fiber in the MDF became flat and straight and the adhesion between resin coated fibers was significantly improved and led to better entanglements.

RSM Modeling
According to the results (Table 3), changes in the strength ranged from 0.5 MPa to 4.664 MPa. The maximum bending strength belongs the specimen made by the resin consumed in the glue line with the MR of 1.68:1, the PT of 1 cm and the WR of 70:30. However, the minimum strength belongs to the boards with the MR of 1.68:1, the PT of 3 cm and the WR of 30:70. Therefore, it is worthy to note that the strength is maximum in the boards with a similar and minimum MR (1.68:1), as the board thickness has decreased while the WR is at maximum. Therefore, this finding indicates the potential of the application of the modified protein adhesive as a replacement for the MUF resin, even when less MR is used.
To determine the desirability for calculating the required response, after completing the experiments the quadratic model was chosen among different models including linear, interactive and cubic models. For this purpose, two important statistics, including the sequential model sum of squared and lack of fit tests, were used to describe the model performance (Table 4). It became clear that the sequential model fitting of the quadratic model was significant, while the lack of fit of the quadratic model was not significant. Thus, for the rest of the analytic process, the quadratic model was chosen for further analysis of the results. The ANOVA results and the relationship between the independent and dependent variables are given in Table 5. A minimum p-value (confirming the significance of any regression coefficient of the independent variables) always shows the most significant factor. Based on these criteria, two tested models indicate a p-value with a high significance (p < 0.0001), thereby emphasizing that they are suitable to affect the MOR. In addition, the determination coefficient (R 2 ) is 0.9960. Likewise, the adjusted coefficients of determination (adj. R 2 ) are 0.9673. Hence, the results obtained confirm that the model is significant (p < 0.05) and show that the experimental data fit the second-order polynomial equations well [35]. Furthermore, the very small coefficients of variations (c.v. < 10) indicate a high precision as well as the suitable validity of the MOR experimental values [36]. The model can predict with an adequate precision in the range of the experimental variables. Moreover, it measures a signal to noise ratio with a value larger than 4 and an adequate precision, which is desirable [37]. The coefficient of 78.11 indicates an adequate signal. However, the linear, interactive and quadratic regression coefficients are significantly different, with very small p-values for the response (p < 0.05). In the manufacturing process, all the independent variables significantly affect the MOR (p < 0.05). In particular, the MOR is significantly affected by the direct effect of the independent variable, "the panel thickness" (x 2 ), and the square of this variable (x 2 2 ) so that the increase in the response is inversely affected.

The Interaction between the Independent Variables in the RSM Approach
The interactive effect between the process parameters was described by the contour plots ( Figure 4). The results showed that as the panel thickness decreases, the bending strength strongly increases while the simultaneous effect of the MR on the change in the bending strength is lessened and the maximum increasing effect belongs to the MR of 1.805:1 (the axial point), and the effect is decreasing more than or less than at this middle limit ( Figure 4A). However, according to the regression coefficients β i and their sign, the direct effect of both factors (MR and PT) is similar, and the direct effect of PT is much more than that of MR (but it is inverse according to the negative β i coefficient). According to the regression equation in Table 4, it is observed that although the interactive effect β ij is significant, the direct effect of the PT is much more than the interactive effect of both factors (more than 10 times). The changes in the curve's slope in Figure 4 have proved it. When the PT is at the middle (2 cm) (as a center point), the interactive effect of the MR and the WR indicates that, as the WR increases during the MR changes, the bending strength increases considerably ( Figure 4B). The effect of the MR is slight according to the intangible changes in the colorful background in the whole x-axis versus the y-axis. According to the equation given in Table 4, with the regression coefficient β ij equal to 0.1563 and the small regression coefficients β ii equal to 0.0813 for the MR and 0.3226 for the WR, it becomes clear that the MR effect is much less than the PT effect (by about 4 times). This can be inferred from the sharp decrease in the changes along the x-axis (denoting the MR) when compared to the intensity of the changes in the y-axis (denoting the PT). Based on the regression equation, the interactive effect of the PT and the WR is significant on the MOR, though it is decreasing ( Figure 4C). It is also observed that in the curve, a large area of the plot belongs to the regions offering the minimum strength when applying more MUF resin (70%) for the panels with a thickness of 2 cm. As the WR becomes maximum (70:30) and the thickness becomes minimum (10 mm), the change in the color intensity of the plot's area indicates a maximum strength where the MR is at the middle (1.805:1). It is also evident from the regression coefficient β ij belonging to x 2 × 3 that it is much bigger than other β ij coefficients, thus indicating that this interaction has a maximum effect on the MOR. Due to the negative coefficient, it can be stated that the positive effect of the thickness decrease on the increase in the bending strength is much more than the positive effect of the increase in the protein adhesive percentage. This can be also deduced by comparing the regression coefficients β i belonging to the direct effect of x 2 and x 3 .

The Accuracy of ANN in the Output Prediction
The total MSE for the training and test data indicates the prediction accuracy of the ANN [38]. The selection of six hidden neurons resulted in the minimum MSE and the maximum correlation coefficient R. Figure 5 shows the plot's performance including the errors of the training, validation and test. The best validation performance occurs in epoch 3 where the validation error is minimum (0.096). No significant overfitting has occurred up to iteration 3. After epoch 3, the test and validation set errors show a similar trend until the validation stops at epoch 8.

The Interaction between the Independent Variables in the ANN Model Approach
Three possible states can occur: when the WR (50:50) was constant and the MR and the PT changed; when the PT (2 cm) was constant and the MR and the WR changed; and when the MR (1.805:1) was constant and the PT and the WR changed. As a result of this process, the results obtained for the MOR by the ANN models are given in Figure 7 graphically without doing experiments. The results obtained from Table 3 and Figure 6 indicate that the model offered by the ANN can be used effectively as a prediction instrument to determine the suitable construction conditions to achieve the highest bending strength of the sandwich panel. It may be noted that these construction conditions are very important industrially since the development of any industry strongly depends on the decrease in the production cost and the increase in the efficiency to survive competition.   As observed in Figure 7, the maximum bending strength belongs to the panels wit an average MR (1.805:1) but with a minimum thickness (1 cm) and a maximum WR. Th minimum bending strength is where the MR is minimum and maximum, while the P and the WR are maximum and minimum, respectively. The trend of the changes in th MOR of the sandwich panel in the model offered by the ANN are similar to the mode As observed in Figure 7, the maximum bending strength belongs to the panels with an average MR (1.805:1) but with a minimum thickness (1 cm) and a maximum WR. The minimum bending strength is where the MR is minimum and maximum, while the PT and the WR are maximum and minimum, respectively. The trend of the changes in the MOR of the sandwich panel in the model offered by the ANN are similar to the model offered by the RSM.
As the core thickness becomes maximum, deflection increased so that the specimens with a lower thickness behaved similarly to a material without any core. In other words, the top and bottom face sheets and the core showed properties similar to a laminate composite plate [39]. This increased the horizontal shear stress due to the increase in the deflection, and delamination occurs easily in the natural axis. Visual observation confirmed that the failure location of the specimens with a high thickness was generally at the middle of the core, and as the thickness decreased the failure mode was mainly a combination of the failure of the face sheets and the glue line. The increase in the ratio of the thickness of the core to the face sheets resulted in the development of more microcracks [40]. The beginning of any failure was the result of the speed of the propagation of the microcracks. As a result of the increase in the propagation of the microcracks, the stress transfer was lessened from one surface to another surface while its uniform distribution decreased and concentrated in the cracked regions. Moreover, this can decrease the panel's bending strength. Generally, due to resin penetration from the glue line to the polyurethane's porous texture, as the thickness increased the ratio of the resin penetration depth to the core thickness decreased. Therefore, in practice it is more likely that microcracks would begin and develop more quickly.
Generally, as the core thickness decreases, the panel's compressive strength increases due to the development of a border layer adjacent to the face sheet where the adhesive improves the local strength of the core cell walls. This strengthening effect can affect the behavior of the damage in the face sheets, which causes the specimens with a low thickness to be dominated by a high ratio of the core cells strengthened by the face sheet failure mechanisms [2]. As the thickness increases, the ratio of these cells decreases, and the face sheet failure mechanism is less predominated. As the thickness reaches the average value (20 mm), the face sheets failure and the core failure become affected. Since the specimens with a lower thickness show a behavior similar to the material without any core (i.e., the MDF used in the face sheets), the core does not significantly affect the applied force and the rigidity of the specimens is greater than those specimens with a larger thickness, and they behave like a laminate composite [39]. However, an increase in the thickness to 30 mm can have a small interactive effect on the thin face sheets and the core, and more load is borne in the core due to the damage caused by the high indentation.
Different failure modes can occur in the sandwich panel, and one of the most important failure modes is core collapse due to the local compressive loading. It is difficult to see this failure visually, but its importance must be considered so to significantly decrease the loadcarrying capacity. Since the core is flexible, the surface layer is damaged before the core in the three-point loading. When the top surface is under compressive stress, the bottom surface of the specimen is simultaneously under tensile stress. The connection between the surface layer and the core is weak due to the relative stress conditions. The brittle surface layer breaks together with the interphase layer and deformation increases. Therefore, the specimen's load-bearing capacity decreases so that the deformation mechanism of the sandwich panel is strongly dominated by the mechanical properties of the surface layers that bear the main compression/tension of the sandwich panel in the deflection. The effect of this deformation increases as the thickness to length ratio decreases, and the stress occurs more in those specimens with a high thickness to length ratio while the bending strength decreases at a constant deformation. In other words, as the thickness increases there is more deformation, and a higher stress exists in the face sheets.
Based on the transformed cross-section approach, it can be assumed that a sandwich panel can act as an I-beam in the bending loading wherein the core and the surface layers act as a web and flange, respectively. On the other hand, the rupture strength of an I-beam is calculated by calculating the maximum bending moment and the moment of inertia (Equations (8)-(10)) [41]: where I is the beam's moment of inertia, w is the beam's upper flange width, d is the beam's height, f is the beam's flange thickness, b is the beam's lower flange width, P is the concentrated loading force applied to the middle of the beam and L is the span length. However, according to Equation (9), the maximum bending moment of the beam is used in the concentrated three-point loading. For the I-beams, the maximum stress or the strength of the beam is as follows: where P is the load applied to the beam, D is the bending stiffness of the beam, L is the span length and I is the moment of inertia. An inverse relationship is observed between the MOR and the moment of inertia. Also, as the thickness increases, the moment of inertia strongly increases (in the first term of Equation (8), which is related to the moment of inertia of the surface layer, the thickness has a power 3; and in the second term, which is related to the moment of inertia of the core, the thickness has a power 2). In Equation (10), the MOR is affected significantly and inversely by the panel thickness so that the minimum MOR is where the thickness is maximum.
Based on the regression equation in Table 5, not only the linear coefficients but also the quadratic coefficients of the interactive effect of the independent variable "PT, (x 2 )" are much bigger than the other coefficients. This indicates that the effect of the thickness on the response surface (MOR) has been much larger than the in other independent variables. This may be due to the emergence of the core failure at smaller loading values in the thicker specimens.
The effect of the resin viscosity on the glue line strength is very significant. The adhesive fluidity decreases as the viscosity increases, and the resin cannot coat the surface well. However, a low viscosity adhesive can extensively penetrate into the wood components and foam, which can lead to the absence of the adhesive layer on the surface. During the alkaline treatment and oxidation, the viscosity of the protein adhesive increases so that its fluidity significantly decreases due to the destruction of the internal hydrogen bonds of the protein and its exposure to the chemical groups [42]. By combining the adhesive with the MUF resin with different MRs affecting the MW directly, the fluidity can be differently affected. Since as the formaldehyde molar ratio increases, there is polymerization and the viscosity increases, combining the resin with the protein adhesive which results in an increase in the viscosity. However, a resin with a smaller formaldehyde MR will result in a smaller viscosity combined with the protein adhesive due to the decrease in the polymerization process and the decrease in viscosity. Therefore, the change in the viscosity and the corresponding change in the adhesive fluidity can offer a suitable level of the composition in which the glue line strength is maximum. It was visually observed in the bending experiments that, in the specimens in which the formaldehyde MR is at the middle level (1.805) and the protein consumption is maximum, the failure was completely in the range of the neutral axis. However, in the specimens in which the formaldehyde MR was minimum or maximum, the failure occurred mainly in the glue line between the surface layer and the core.
Proteins normally have a compact spherical structure without any unfolding that mainly forms the dense layers and sometimes the brittle particles due to adsorption [43], which leads to a weak interfacial strength or bonding strength in the adhesive. Unfolding the structure by different treatments can not only release the jointed and hidden polar groups but can create additional cohesion strength of the adhesive due to the inherency of the intermolecular linkages, which can increase the efficiency and create a stronger contact to the wood surface by absorption [44].
Studies have shown that if the oxidized protein is not used in the MUF resin, the solid content increases while the viscosity significantly decreases, and the protein modification does not meaningfully affect the solid content. As more protein is used, the solid content decreases while the viscosity significantly increases [42]. It is known that the protein's isoelectric point is in a pH ranging from 7 to 8 [45]. This pH range must lead to a higher negative charge on the protein molecule and in the unfolding of the 3D protein structure by splitting both the intra-and intermolecular bonds. Therefore, the unfolded molecules (due to the action of NaOH as a denaturing agent) can result in a greater entanglement of the chain via the higher number of functional groups in the protein molecule that are accessible for more intermolecular interactions (crosslinking), especially during the application of the final heating in the press and can also the increase in the dispersion viscosity of the adhesive [46]. This can decrease the penetration due to the lower fluidity in the face sheets and there will be more adhesive in the glue line while more adhesive can penetrate the porous material of the core foam and smear a larger volume of the layers close to the glue line in practice. As a result, the bending stress can be transferred deeply inside the foam and decrease.
After the alkaline treatment, the reactive groups are exposed through the destruction of the Van der Waals forces and hydrogen bonds and hydrophobic interactions between the molecules in the native protein, so that the extra formaldehyde and urea in the resin can modify, unfold and then stabilize the protein from additional refolding [47]. These molecules thereby increase with growing dimensions due to the formation of new crosslinking between the methyl groups of the MUF resin and the reactive amine and amide bonds of the protein and the hydroxymethyl formation [42]. Concurrently, hydroxymethyl can react with the urea, methylol urea, amine or hydroxymethyl protein hydrolyzed under acidic conditions so to form methylene ether or methylene linkages [48]. More frequent methylene ether bonds can create methylene bridges, branching reactions and crosslinking in the resin to form a 3D network. The results indicate that the activation energy of the aminoplast adhesive modified by the hydrolyzed protein decreases so that, due to the presence of the hydrolyzed protein in the network structure, the methylene ether bridges can easily return, and the methylene ether bridges decrease in the cured adhesive containing the modified protein [48].
Based on the rheological behavior of the adhesives containing low values of protein, the distance between the random coils is more than the radius of gyration. Therefore, entanglement does not occur, and the suspension viscosity is comparable with the solvent. As the protein increases, the relative viscosity increases gradually to a threshold concentration beyond which the changes in the slope of the viscosity increase and are sudden and fast, and the distance between the random coils and the radius of the gyration becomes almost equal. As the solution's concentration increases beyond the critical concentration, the entanglement of the protein polymer is dominated and the relative viscosity of the complex increases [33]. It is known that even at a low protein content, the dispersed protein concentration is more than the critical concentration and an entanglement structure is formed with a pseudoplastic behavior [33], and as the protein content increases, the flow shows a dilatant plastic behavior [49]. This event may demonstrate the potential of the formation of a physical phenomenon to use the urea and formaldehyde in the system as crosslinking agents to form chemical bonds such as methylene and methyl ether bonds.

Comparison of the RSM and ANN Models
The comparison of the output of the ANN and RSM models with the measured bending strength values of the sandwich panel is given in Figure 8. The results of the graphical comparisons indicate a similarity between the experimental outputs and the ANN and the RSM outputs. As observed, the predicted outputs overlap with the measured outputs in many cases. Therefore, it can be said that the proposed model is probably trained and there is a very good precision in the bending strength prediction of the specimens. methylene and methyl ether bonds.

Comparison of the RSM and ANN Models
The comparison of the output of the ANN and RSM models with the measure bending strength values of the sandwich panel is given in Figure 8. The results of th graphical comparisons indicate a similarity between the experimental outputs and th ANN and the RSM outputs. As observed, the predicted outputs overlap with th measured outputs in many cases. Therefore, it can be said that the proposed model probably trained and there is a very good precision in the bending strength prediction the specimens. The goodness-of-fit of the ANN and the RSM models was evaluated by the regression coefficient R 2 . Figures 9 and 10 offer the fit between the predicted and the measured outputs of the bending strength for these models graphically. According to these Figures, R 2 is 0.9936 and 0.9841 for the test data set and 0.9969 and 0.9960 for all data sets based on the ANN and the RSM models, respectively. The R 2 results of the ANN model indicate that the selected model agrees with at least 99.36% and 99.69% of the measured test and all data sets of the bending strength. However, the RSM model can explain 98.41% and 99.60% of the deviations of the bending strength for the test data set and all data sets, respectively. Therefore, the bending strength prediction by the ANN is more precise than the RSM. However, both the ANN and the RSM could offer an excellent estimate of the bending strength with high reliability [50]. that the selected model agrees with at least 99.36% and 99.69% of the measured test and all data sets of the bending strength. However, the RSM model can explain 98.41% and 99.60% of the deviations of the bending strength for the test data set and all data sets, respectively. Therefore, the bending strength prediction by the ANN is more precise than the RSM. However, both the ANN and the RSM could offer an excellent estimate of the bending strength with high reliability [50].  The general predictability of a model is usually determined by R 2 , but the model performance may not be described by R 2 alone. In the multiple linear regression model, the adjusted R 2 also calculates the variations of a dependent variable and usually measures the goodness-of-fit with a higher precision compared to R 2 . Overall, high R 2 and adjusted R 2 values do not always mean that the regression model is an efficient model. In a good model, the AAD value must be as small as possible while the RMSE value must be close to zero. High RMSE and AAD values mean a higher probability of errors in the prediction. According to Table 6, the AAD value (1.28) for the RSM is three times more than that for the ANN (0.43). In sum, the RSME value (1.98) for the RSM is more than three times the value for the ANN (0.56). Moreover, the MAPE value (3.45) for the RSM is more than three times the value for the ANN (1.04). The prepared ANN had higher R 2 and adjusted R 2 values, while its AAD, RMSE and MAPE values were less than those for the RSM. Therefore, the ANN shows a prediction capacity which is significantly higher than the RSM. that the selected model agrees with at least 99.36% and 99.69% of the measured test and all data sets of the bending strength. However, the RSM model can explain 98.41% and 99.60% of the deviations of the bending strength for the test data set and all data sets, respectively. Therefore, the bending strength prediction by the ANN is more precise than the RSM. However, both the ANN and the RSM could offer an excellent estimate of the bending strength with high reliability [50]. The general predictability of a model is usually determined by R 2 , but the model performance may not be described by R 2 alone. In the multiple linear regression model, the adjusted R 2 also calculates the variations of a dependent variable and usually measures the goodness-of-fit with a higher precision compared to R 2 . Overall, high R 2 and adjusted R 2 values do not always mean that the regression model is an efficient model. In a good model, the AAD value must be as small as possible while the RMSE value must be close to zero. High RMSE and AAD values mean a higher probability of errors in the prediction. According to Table 6, the AAD value (1.28) for the RSM is three times more than that for the ANN (0.43). In sum, the RSME value (1.98) for the RSM is more than three times the value for the ANN (0.56). Moreover, the MAPE value (3.45) for the RSM is more than three times the value for the ANN (1.04). The prepared ANN had higher R 2 and adjusted R 2 values, while its AAD, RMSE and MAPE values were less than those for the RSM. Therefore, the ANN shows a prediction capacity which is significantly higher than the RSM. The general predictability of a model is usually determined by R 2 , but the model performance may not be described by R 2 alone. In the multiple linear regression model, the adjusted R 2 also calculates the variations of a dependent variable and usually measures the goodness-of-fit with a higher precision compared to R 2 . Overall, high R 2 and adjusted R 2 values do not always mean that the regression model is an efficient model. In a good model, the AAD value must be as small as possible while the RMSE value must be close to zero. High RMSE and AAD values mean a higher probability of errors in the prediction. According to Table 6, the AAD value (1.28) for the RSM is three times more than that for the ANN (0.43). In sum, the RSME value (1.98) for the RSM is more than three times the value for the ANN (0.56). Moreover, the MAPE value (3.45) for the RSM is more than three times the value for the ANN (1.04). The prepared ANN had higher R 2 and adjusted R 2 values, while its AAD, RMSE and MAPE values were less than those for the RSM. Therefore, the ANN shows a prediction capacity which is significantly higher than the RSM. There has been very limited information on the application of the ANN to predict the mechanical properties of the wood-based sandwich panel so far, and most data are on other wood panels. The R 2 results obtained from the ANN sometimes show higher values when compared to other wood-based composites. Fernandez et al. [51] have obtained the R 2 value of 0.66, while Bardak et al. [52], Bardak et al. [53], Demirkir et al. [54] and Nazerian et al. [22] have offered the R 2 values 0.97, 0.94, 0.97 and 99.5, respectively. Satisfactory results were observed in these, and many other, studies based on wood products.
In Figure 11, the residual errors distribution is plotted for both optimization techniques. The variations of the residuals are completely small, while there is no significant difference between the ANN and the RSM. However, based on the R 2 , it can be stated that the ANN-based residual value can be less than that which is based on the RSM. This means that the experimental data have a fit with a high accuracy, as both the RSM and ANN models are applied. Therefore, while artificial intelligence methods such as the ANN can be trained to estimate the nonlinear functions, and can be used to estimate the experimental data, they offer an opportunity to match any experimental set up with the model construction simultaneously, such that they are so flexible that they make it possible to add new experimental data to a reliable model construction. Although the RSM is only limited to the second-order polynomials [55], they can offer a regression equation to predict and show the effect of the experimental parameters and their interactive effect on the response being examined compared to the ANN.

Conclusions
In this study, the effect of the formaldehyde to melamine-urea molar ratio (MR), panel thickness (PT) and the replacement percentage of the modified soybean protein weight ratio (WR) was investigated on the bending strength of the polyurethane foamed-cored sandwich panels using the RSM and the ANN. The most important conclusions of this study are as follows.
The best MR is 1.805:1, and as this ratio increases or decreases, the bending strength of the sandwich panel decreases even when the other production variables decrease or increase. As the PT increases, the bending strength decreases. The increase in the WR results in a continuous increase in the bending strength.
The application of the ANN and the RSM models based on the experimental results indicated that they are powerful and efficient instruments to predict the bending strength.
The performance of the RSM and the ANN models was evaluated based on the statistical indices such as R 2 , RMSE, AAD and MAPE. The calculated R 2 values 0.9969 and 0.9960, respectively, based on the ANN and the RSM models show that the selected model agrees with at least the 99.69% and 99.60% of the measured bending strength value and can explain the bending strength deviations.
The AAD is 1.28 for the RSM, while it is only 0.43 for the ANN; the RSME is 1.98 for the RSM, while it is only 0.56 for the ANN; the MAPE is 3.45 for the RSM, while it is only