Quantitative Structure–Electrochemistry Relationship (QSER) Studies on Metal–Amino–Porphyrins for the Rational Design of CO2 Reduction Catalysts

The quantitative structure–electrochemistry relationship (QSER) method was applied to a series of transition-metal-coordinated porphyrins to relate their structural properties to their electrochemical CO2 reduction activity. Since the reactions mainly occur within the core of the metalloporphyrin catalysts, the cluster model was used to calculate their structural and electronic properties using density functional theory with the M06L exchange–correlation functional. Three dependent variables were employed in this work: the Gibbs free energies of H*, C*OOH, and O*CHO. QSER, with the genetic algorithm combined with multiple linear regression (GA–MLR), was used to manipulate the mathematical models of all three Gibbs free energies. The obtained statistical values resulted in a good predictive ability (R2 value) greater than 0.945. Based on our QSER models, both the electronic properties (charges of the metal and porphyrin) and the structural properties (bond lengths between the metal center and the nitrogen atoms of the porphyrin) play a significant role in the three Gibbs free energies. This finding was further applied to estimate the CO2 reduction activities of the metal–monoamino–porphyrins, which will prove beneficial in further experimental developments.


Introduction
The combustion of fossil fuels has contributed to the accumulation of CO 2 in the atmosphere, which has led to significant increases in the greenhouse effect [1]. To achieve a carbon-neutral economy, increasing efforts have been devoted to converting CO 2 into fuel molecules or value-added chemical products [2,3]. Among numerous strategies, the electrochemical reduction of CO 2 is particularly appealing due to its mild reaction conditions, its tunable external parameters, such as electrode potentials, and available operating devices [2][3][4]. In addition, it uses water as a hydrogen source, makes full use of environmentally friendly renewable energy sources (such as solar, tidal, and wind energy), and can simultaneously contribute to carbon recycling and utilization and renewable electricity storage [5,6]. However, challenges remain in the electrochemical CO 2 reduction reaction (CO 2 RR) since CO 2 is a linear molecule with strong C=O bonds and requires a large overpotential for its activation [5]. Moreover, current electrochemical catalysts consistently suffer from low Faradaic efficiency (FE), high overpotential, and poor current potential for its activation [5]. Moreover, current electrochemical catalysts consistently suffer from low Faradaic efficiency (FE), high overpotential, and poor current density because of the competing hydrogen evolution reaction (HER) [7][8][9]. Therefore, the rational design of electrochemical catalysts with excellent catalytic performance and selectivity is imperative.
Metalloporphyrin complexes have been used as catalysts for CO2 reduction on account of their unique structure and electronic properties [10][11][12][13][14][15][16]; moreover, by changing the ligands on porphyrin complexes, one can modulate the different catalytic activities of such complexes [17][18][19][20][21]. Liu et al. used density functional theory (DFT) to design a series of transition-metal-porphyrins (TM-PPs) and found that the TM-PP monolayers display excellent catalytic stability and electrochemical CO2 reduction selectivity [22]. In another study, Davethu et al. tested several DFT methods to elucidate the ligand effect in the CO2RR performance of iron-porphyrins, finding that second-coordination-sphere perturbations influence CO2 positioning on the metalloporphyrins [23]. Abdinejad et al. observed that Fe-amino-porphyrins can selectively reduce CO2 to carbon monoxide (CO) at ambient pressure and temperature with competitive turnover numbers (TONs) [24]. Further, Wu et al. investigated the electrochemical CO2 reduction performance of a heterogenized zinc-porphyrin complex and demonstrated that such an electrocatalyst delivers a turnover frequency as high as 14.4 site −1 s −1 and an FE as high as 95% for the electroreduction of CO2 to CO at −1.7 V versus the standard hydrogen electrode [25].
During the CO2 reduction reaction, the first protonation step is CO2 + H + + e − → C*OOH or O*CHO [26]. These two main intermediates (C*OOH and O*CHO) will lead to the products of CO and HCOOH, respectively [27,28]. In the meantime, a hydrogen evolution reaction (HER) will also occur and it will form H* intermediates; this reaction is a competitive reaction with CO2RR [29]. The three important intermediates are shown in Figure 1. In theoretical studies, the Gibbs free energy of H* formation, G(H*), has been used to determine the favorable reaction between the hydrogen evolution reaction and the CO2RR [30], while the Gibbs free energies G(C*OOH) and G(O*CHO) can be used to predict the favored products [31][32][33][34][35][36]. Thus, based on theoretical catalytic studies, G(H*), G(C*OOH), and G(O*CHO) can be calculated and compared to determine the product selectivity. The quantitative structure-electrochemistry relationship (QSER) method can be used to study the relationship between the structures and electrocatalytic activity of catalysts and is widely applied in electrocatalysis research [37,38]. In other words, this approach can be used to correlate a catalyst's activity with its structural and electronic properties. The quantitative structure-electrochemistry relationship (QSER) method can be used to study the relationship between the structures and electrocatalytic activity of catalysts and is widely applied in electrocatalysis research [37,38]. In other words, this approach can be used to correlate a catalyst's activity with its structural and electronic properties. In this work, the concept of QSER with the genetic algorithm combined with multiple linear regression (GA-MLR) was applied to develop the mathematical models for the Gibbs free energies of the three intermediates (H*, C*OOH, and O*CHO), and DFT was used to investigate the structural and electronic properties of the metal-porphyrin catalysts. We first obtained a series of 10 TM-PPs (TM = Sc, Ti, V, Cr, Mn, Fe, Co, Ni, Cu, and Zn) and their G(H*), G(C*OOH), and G(O*CHO) values from the literature [22]. Subsequently, we used the collected data to train the QSER models and then used them to predict the product selectivity of the CO 2 RR for a series of newly designed transition-metal monoaminosubstituted tetraphenylporphyrins (TM-Amino-TPPs, TM = Sc, Ti, V, Cr, Mn, Fe, Co, Ni, Cu, and Zn); the workflow is depicted in Figure 2. s 2023, 28, x FOR PEER REVIEW 3 of 14 In this work, the concept of QSER with the genetic algorithm combined with multiple linear regression (GA-MLR) was applied to develop the mathematical models for the Gibbs free energies of the three intermediates (H*, C*OOH, and O*CHO), and DFT was used to investigate the structural and electronic properties of the metal-porphyrin catalysts. We first obtained a series of 10 TM-PPs (TM = Sc, Ti, V, Cr, Mn, Fe, Co, Ni, Cu, and Zn) and their G(H*), G(C*OOH), and G(O*CHO) values from the literature [22]. Subsequently, we used the collected data to train the QSER models and then used them to predict the product selectivity of the CO2RR for a series of newly designed transition-metal monoamino-substituted tetraphenylporphyrins (TM-Amino-TPPs, TM = Sc, Ti, V, Cr, Mn, Fe, Co, Ni, Cu, and Zn); the workflow is depicted in Figure 2.

Correlation between DFT Structural and Electronic Properties and Gibbs Free Energies
In this work, the structural and electronic properties were obtained, consisting of the natural bond orbital (NBO) charges of the TM (X1); average NBO charges of Na and Nb (X2), Ca (X3), and Cb (X4); dipole moment (X5); polarizability (X6); average TM-Na and TM-Nb bond length (X7); and average Na-Nb distance (X8), as shown in Figure 3 and Table  1. The correlation of the three Gibbs free energies and the eight DFT properties of the TM-PP complexes have also been analyzed using the Pearson correlation coefficient, as shown in

Correlation between DFT Structural and Electronic Properties and Gibbs Free Energies
In this work, the structural and electronic properties were obtained, consisting of the natural bond orbital (NBO) charges of the TM (X1); average NBO charges of N a and N b (X2), C a (X3), and C b (X4); dipole moment (X5); polarizability (X6); average TM-N a and TM-N b bond length (X7); and average N a -N b distance (X8), as shown in Figure 3 and Table 1. The correlation of the three Gibbs free energies and the eight DFT properties of the TM-PP complexes have also been analyzed using the Pearson correlation coefficient, as shown in

GA-MLR QSER-Model Relationship
The GA method was used to identify the molecular properties that were highly correlated with the three Gibbs free energies. The occurrences of the populations of the eight descriptors for the three models are depicted in Figure S1, indicating that G(H*) and G(C*OOH) are highly correlated with X1, X3, X7, and X8, while G(O*CHO) has a high correlation with X1, X2, X6, and X7. The obtained QSER equations to predict the Gibbs free energy values of H* (Y1), C*OOH (Y2), and (O*CHO) (Y3) for TM-PPs with MLR and GA-MLR methods are shown in Tables S1-S3. Table 2 shows the MLR models of Y1, Y2, and Y3 are denoted as Equations (1)-(3), respectively, while the selected GA-MLR models of Y1, Y2, and Y3 are represented in Equations (1A)-(3A), respectively. The QSER models

GA-MLR QSER-Model Relationship
The GA method was used to identify the molecular properties that were highly correlated with the three Gibbs free energies. The occurrences of the populations of the eight descriptors for the three models are depicted in Figure S1, indicating that G(H*) and G(C*OOH) are highly correlated with X1, X3, X7, and X8, while G(O*CHO) has a high correlation with X1, X2, X6, and X7. The obtained QSER equations to predict the Gibbs free energy values of H* (Y1), C*OOH (Y2), and (O*CHO) (Y3) for TM-PPs with MLR and GA-MLR methods are shown in Tables S1-S3. Table 2 shows the MLR models of Y1, Y2, and Y3 are denoted as Equations (1)-(3), respectively, while the selected GA-MLR models of Y1, Y2, and Y3 are represented in Equations (1A)-(3A), respectively. The QSER models

GA-MLR QSER-Model Relationship
The GA method was used to identify the molecular properties that were highly correlated with the three Gibbs free energies. The occurrences of the populations of the eight descriptors for the three models are depicted in Figure S1, indicating that G(H*) and G(C*OOH) are highly correlated with X1, X3, X7, and X8, while G(O*CHO) has a high correlation with X1, X2, X6, and X7. The obtained QSER equations to predict the Gibbs free energy values of H* (Y1), C*OOH (Y2), and (O*CHO) (Y3) for TM-PPs with MLR and GA-MLR methods are shown in Tables S1-S3. Table 2 shows the MLR models of Y1, Y2, and Y3 are denoted as Equations (1A,2A,3A), respectively, while the selected GA-MLR models of Y1, Y2, and Y3 are represented in Equations (1B,2B,3B), respectively. The QSER models with the MLR method resulted in a low predictive ability of R 2 (CV), which cannot be used for the further prediction of new design catalysts. Whereas the GA-MLR technique could improve the statistical values of the predictive ability of R 2 (CV), these obtained QSER with GA-MLR models are used for predicting the new design of metal-monoamino-porphyrin catalysts. Table 2. Equations of the lines of best fit for predicting the values of Y1 (G(H*)), Y2 (G(C*OOH)), and Y3 (G(O*CHO)) for the TM-PPs. QSER GA-MLR of G(H*) model. Equation (1B) (see Table 2) represents the QSER model obtained from the GA-MLR method for the prediction of G(H*).  Table 2

Newly Designed Compounds and Predicted Activities
As is commonly known, the ligands of porphyrin significantly affect its activity. For example, Savéant et al. introduced phenolic hydroxyl groups to the second coordination sphere of an Fe-porphyrin, which greatly improved the selectivity of the CO2RR [39]. In addition, Nichols et al. utilized second-sphere effects to design a series of iron-porphyrins to promote CO2 reduction selectivity [40]. In this study, a series of new TM-PP complexes were designed, introducing phenyl monoamino ligands into the porphyrin complexes; these were denoted as TM-Amino-TPPs, as shown in Figure 6. The structure was optimized using the M06L/6-31G (d, p) basis set, which is similar to that used with the training data set of the TM-PP catalysts. According to the QSER models with GA-MLR, three properties of the TM-Amino-TPPs were collected, namely, the NBO charges of TM (X1), The average NBO charges of Ca (X3), and the average bond length of TM-Na and TM-Nb (X7), as presented in Table S5. These properties were then used with the models in Table 2 to predict the values of G(H*), G(C*OOH), and G(O*CHO) from the Equations (1A)-(3A), respectively. The Gibb free energy predictions of the newly designed TM-Amino-TPPs are listed in Table 3.

Newly Designed Compounds and Predicted Activities
As is commonly known, the ligands of porphyrin significantly affect its activity. For example, Savéant et al. introduced phenolic hydroxyl groups to the second coordination sphere of an Fe-porphyrin, which greatly improved the selectivity of the CO 2 RR [39]. In addition, Nichols et al. utilized second-sphere effects to design a series of iron-porphyrins to promote CO 2 reduction selectivity [40]. In this study, a series of new TM-PP complexes were designed, introducing phenyl monoamino ligands into the porphyrin complexes; these were denoted as TM-Amino-TPPs, as shown in Figure 6. The structure was optimized using the M06L/6-31G (d, p) basis set, which is similar to that used with the training data set of the TM-PP catalysts. According to the QSER models with GA-MLR, three properties of the TM-Amino-TPPs were collected, namely, the NBO charges of TM (X1), The average NBO charges of C a (X3), and the average bond length of TM-Na and TM-Nb (X7), as presented in Table S5. These properties were then used with the models in Table 2 to predict the values of G(H*), G(C*OOH), and G(O*CHO) from the Equations (1B,2B,3B), respectively. The Gibb free energy predictions of the newly designed TM-Amino-TPPs are listed in Table 3.   Comparing the X1, X3, and X7 values of the TM-PPs and the TM-Amino-TPPs, one can note that with the introduction of amino ligands, the NBO charges of the TM and the average TM-Na and TM-Nb bond length for the porphyrins generally decrease. The TM-Amino-TPPs have less positive NBO charges on the TM and shorter TM-N bond lengths, while the amino substituent increases the average NBO charges of C a for these porphyrins, except for the case of Sc-Amino-TPP. The values of the descriptors X1, X3, and X7 of the TM-PPs versus those of the TM-Amino-TPPs are plotted in Figure S2, showing that the TM NBO charges of Cr-and Mn-Amino-TPP decrease the most and the NBO charges of C of these same two porphyrins increase the most with the introduction of amino ligands. This indicates that the incorporation of amino ligands could redistribute the electrons of the porphyrin molecules; the Cr and Mn metal centers have a greater positive charge, and they can adsorb the CO 2 − anion more easily. Gibbs free energy is one of the critically important values for determining the preferable reaction pathways. Thus, in case that the G(H*) values are lower than that of G(C*OOH) and G(O*CHO), the HER reaction is preferable, while the candidate catalysts for the CO 2 reaction pathway, their G(C*OOH) or G(O*CHO) values are lower than that of G(H*). Therefore, in the newly designed TM-Amino-TPPs based on the QSER predictions (Table S4), the Sc-, Ti-, V-, Cr-, Co-, Cu-, and Zn-Amino-TPPs show the lower Gibbs free energies of O*CHO compared to H* and C*OOH as shown in Figure 7. Therefore, this can imply that these catalysts could favorably convert CO 2 to HCOOH as O*CHO as an intermediate for HCOOH products. On the other hand, the Mn-, Fe-, and Ni-Amino-TPPs catalysts showed lower Gibbs free energies of C*OOH ( Figure 7); thus, these new catalysts could support the reaction pathway of CO product formation.
Furthermore, considering the crucial property of the NBO charges of the metal center of TM-Amino-TPPs ( Figure S2), they are significantly decreased compared to that of TM-PPs. Therefore, the charges on the TM-PPs catalysts are redistributed with the introduction of monoamino phenyl ligands, resulting in a more positive charge on the metal center, which would become a more active center. Meanwhile, it is worth noting that the QSER trendy prediction for the Fe-Amino-TPP catalyst corresponds well with the experimental work of Abdinejad et al., as the Fe-Amino-TPP catalyst exhibited a good CO 2 -to-CO conversion, with an FE of 49% [24]. Based on the three types of Gibbs free energy QSERmodel predictions, the newly designed Mn-, Fe-and Ni-Amino-TPPs are selected for further reaction mechanism investigations as their predicted G(C*OOH) is lower than the G(O*CHO) and G(H*). Thus, they could benefit from catalyzing the CO 2 to CO conversion selectivity. However, the product between the CO and HCOOH was intensively studied, and it was found that the product selectivity switched by changing the metal center [41]. Therefore, the QSER results only provide the preliminary screening for potential or candidate catalysts for further theoretical reaction mechanism investigation, which need to include the activation energy and the transition state to determine the rate-determining step of the CO 2 RR reaction.

Data for the Training Set
The Gibbs free energies of reaction intermediates determine the most appropriate reaction pathway and products for the CO 2 RR [7]. The lower the Gibbs's free energy, the more likely it is that the reaction will occur. In addition, the HER competes with the CO 2 RR, which needs to be taken into consideration in predicting the reaction selectivity [42]. In this work, we used the Gibbs free energies of the CO 2 RR and HER reaction intermediates C*OOH, O*CHO, and H* as our target activities to build the QSAR models. The data on the TM-PPs were obtained from Liu et al.'s work and are partitioned into three training sets-the Gibbs free energies of the H*, C*OOH, and O*CHO intermediates-denoted as G(H*), G(C*OOH), and G(O*CHO), or Y1, Y2, and Y3, respectively, as presented in Table 4 [22].

Optimization Details and Properties
The frequency of all of the TM-PP and TM-Amino-TPPs complexes were calculated and subsequently optimized using DFT with the M06L exchange-correlation functional and the 6-31G (d, p) basis set in the Gaussian 09 program [43]. Based on the optimized structures, the structural and electronic properties were obtained, consisting of the natural bond orbital (NBO) charges of the TM (X1); average NBO charges of N a and N b (X2), C a (X3), and C b (X4); dipole moment (X5); polarizability (X6); average TM-Na and TM-Nb bond length (X7); and average N a -N b distance (X8), as shown in Figure 3 and Table 1.

QSER Technique Used in This Work
The central idea of genetic algorithms is that the region to be searched is coded into one or multiple strings. The genetic algorithm operates with a set of such strings, termed a population [44]. In QSER studies, to build a model with high correlation and predictive ability, the choice of descriptors is notably significant [45]. Hence, we employed the GA method to effectively select the appropriate variables. As the Gibbs free energies of intermediates play an important role in predicting electrocatalytic performances, we used these values for the H*, C*OOH, and O*CHO intermediates as three dependent variables [38]. Additionally, the correlation of the three Gibbs free energies and the eight DFT properties of the TM-PP complexes have also been analyzed using the Pearson correlation coefficient. The GA-MLR technique was implemented in the Materials Studio program to derive the QSER models [46]. Initially, the training data were fully imported with the maximum equation length set at three as the limited training data; the population and the maximum number of generations were set to 3000 and 500, respectively, and the mutation probability was 0.1.

Statistical Terms for QSER Analysis
The performance of the regression models was measured using Fischer's ratio (Fvalue); R 2 ; the cross-validation R-squared, or R 2 (CV); and the residual sum of squares (RSS) [47]. The F-value reflects the ratio between the variance explained by the model and the variance due to the error in the regression; high F-values indicate that the model is statistically significant. R 2 is a statistical measure of fit that indicates how much of the variation in a dependent variable is explained by the independent variable(s) in a regression model, as given by Equation (4). We used the leave-one-out validation technique to obtain the values of R 2 (CV), a key measure of the predictive power of a model, as calculated in Equation (5); the closer the value is to 1.0, the better the predictive power. The RSS measures the level of variance in the error term, or residuals, of a regression model, calculated according to Equation (6). Ideally, the sum of the squared residuals should be smaller than the sum of squares from the regression model's inputs. We also calculated the root-mean-square error (RMSE) to evaluate the accuracy of our models (Equation (7)); the closer the value is to 0, the greater the prediction accuracy of the models.
where y i is the predicted value and y is the average of y i .

Conclusions
In summary, the exploration of potential CO 2 RR candidate catalysts can be accelerated through molecular structural descriptors obtained from the DFT of the catalyst cluster model, and then the QSER concept. A series of transition-metal-coordinated porphyrins and their calculated CO 2 reduction activity in terms of Gibbs free energies of H*, C*OOH, and O*CHO were used as training data for the QSER model. The genetic algorithm combined with multiple linear regression (GA-MLR) techniques was used to manipulate the QSER mathematical models of all three Gibbs free energies. Based on our QSER models, both the electronic properties (charges of the metal and porphyrin) and the structural properties (bond lengths between the metal center and the nitrogen atoms of the porphyrin) play a significant role in the three Gibbs free energies. The obtained results showed a good predictive ability both (R 2 value) R 2 (CV), which is more than 0.8. Subsequently, QSER with GA-MLR models were then used to predict the Gibbs free energies of H*, C*OOH, and O*CHO of CO 2 RR for a series of newly designed transition-metal monoamino-substituted tetraphenylporphyrins (TM-Amino-TPPs, TM = Sc, Ti, V, Cr, Mn, Fe, Co, Ni, Cu, and Zn). The predicted three kinds of energies found that the Sc-, Ti-, V-, Cr-, Co-, Cu-, and Zn-Amino-TPPs showed lower Gibbs free energies of O*CHO compared to H* and C*OOH; therefore, they could have the potential to convert CO 2 to HCOOH. While considering the Mn-, Fe-, and Ni-Amino-TPPs, which exhibit lower Gibbs free energies of C*OOH than others, it can be deduced that CO 2 -CO conversion selectivity occurs on these catalysts. Overall, the DFT descriptors and the QSER with GA-MLR methods can be used to estimate or screen potential catalysts for CO 2 electroreduction.