Screening copper-based single-atom alloy catalysts for electrochemical nitrogen reduction

Abstract

Highly efficient and selective electrocatalyst exploration is crucial for advancing the green electrochemical technology of direct nitrogen-to-ammonia conversion through the electrocatalytic nitrogen reduction reaction (eNRR). In this study, we have systematically screened Cu based single-atom alloy (SAA) catalysts for electrochemical nitrogen reduction using density functional theory calculations and interpretable machine learning. We reveal that V-, Nb-, Mo-, Ta-, and W@Cu(100) SAAs exhibit exceptional eNRR performance. These candidates demonstrate robust structural stability, superior selectivity, and remarkably low limiting potentials (U_L > −0.30 V) that surpass those of conventional single-atom catalysts. By employing intrinsic physicochemical properties of SAA catalysts as feature descriptors, we have developed an adaptive boosting machine learning model that achieves accurate prediction of limiting potentials. SHapley Additive exPlanations analysis highlights the critical influence of two key parameters: the electron accumulation in adsorbed *N₂ species and the number of valence electrons of transition metals. Subsequent application of the sure independence screening and sparsifying operator algorithm further identifies an optimal combination of four essential features, corroborating the significance of ML-identified descriptors. This multiscale investigation establishes a robust framework for rational design of SAA catalysts, combining first-principles calculations with interpretable machine learning to accelerate the development of sustainable ammonia synthesis technologies.

---

1. Introduction

NH₃ is a vital chemical product widely used in energy, the chemical industry, and agriculture. It occupies an essential position in the national economy.^1–3 Currently, NH₃ is mainly produced industrially by the Haber–Bosch (HB) method.^4–7 This conventional process needs to be carried out at high temperatures (400–500 °C) and high pressures (150–300 atm) under harsh conditions and relies on fossil fuels as the primary energy source.^8,9 The Haber–Bosch process also emits large amounts of CO₂, which can lead to serious air pollution.^10,11 In recent years, the electrochemical nitrogen reduction reaction (eNRR) has attracted much attention because it utilizes renewable energy to drive the reaction at ambient temperature and pressure, providing a new idea for a low-carbon ammonia economy.^12–16 However, conventional electrochemical catalysts usually exhibit low ammonia yields and low selectivity due to the complex six-electron transfer process involved, as well as the high dissociation energy of the N≡N bond (941 kJ mol⁻¹) and the hydrogen evolution reaction (HER), which severely limit their practical applications.^3,17–19

Single-atom alloy (SAA) catalysts, as an extension of single-atom catalysts (SACs), form a unique “isolated active center-metal substrate” synergistic interface by embedding isolated transition metal (TM) atoms into the lattice of the host metal substrate.^20–24 This structure combines the advantages of the high intrinsic activity of atomically dispersed sites and the high electrical conductivity of the metal substrate, while precise optimization of the electronic structure can be achieved through coordination environment modulation. Liu et al. experimentally synthesized PdAu SAA catalysts, which were used to catalyze the selective hydrogenation of 1-hexynes.²⁵ It was shown that the doping of Pd can significantly enhance the reactivity of highly selective Au hydrogenation catalysts at low temperatures. As a class of noble metals, Cu exhibits high chemical stability and strong hydrogenation activity and has been widely used as the host metal substrate in SAAs.^26–29 Cai et al. synthesized Ni-doped Cu-based SAA catalysts to convert nitrate to NH₃ directly. The results showed that the Ni atom could modulate the proton transfer step of the NO₃⁻ electroreduction and reduce the limiting potential of the reaction, resulting in a Faraday efficiency of nearly 100%.³⁰ Du and co-workers found that by doping the Pd atom on the surface of Cu(100), the HER can be effectively suppressed, and the energy barrier of the rate-determining step (RDS) of the eNRR can be significantly reduced, resulting in an NH₃ Faraday efficiency of up to 97.1%.³¹ Furthermore, copper-based SAA catalysts are expected to break the scaling relationship and achieve high eNRR activity and selectivity by optimizing the activation of reactants and the binding strength of intermediates.^32–34 For example, Dai et al. used density functional theory to screen the “structure–performance” relationship of the eNRR driven by Cu-based SAAs, and the results showed that there is a quintuple degenerate d-electron state in SAA catalysts, and the d-electrons can be redistributed in the functionalized orbitals.³⁵

However, although SAA catalysts show great potential for eNRR applications, their structure–activity relationships are still not fully explored. Experiments generally rely on costly traditional trial-and-error methods to screen efficient catalysts.^21,36,37 To realize the rational design of Cu-based SAA catalysts, it is necessary to understand the relationship between their atomic and electronic structural features and catalytic activities. Usually, the adsorption energy of key intermediates of the eNRR on SAAs is a decisive factor affecting catalytic activity.^38,39 The linear relationship between the adsorption energies of eNRR intermediate species on the catalyst surface was determined by the computational hydrogen electrode (CHE) model⁴⁰ proposed by Nørskov et al.,^41,42 and the d-band theory has also been widely used to characterize the trend of change between the activity of transition metal-based catalysts and the adsorption energy of eNRR intermediate species.^43–45 However, numerous studies have shown that the d-band model cannot accurately predict the catalytic activity of SAAs due to their special geometrical structure and electronic characteristics.^35,46,47 In recent years, machine learning (ML) and high-throughput density functional theory (DFT) calculations have provided a new paradigm for the rational design of SAA catalysts.^48–53 By extracting the intrinsic characteristics of dopant atoms (e.g., electronegativity, d-electron number, and atomic radius) and the adsorption energies of reaction intermediates as descriptors, high-precision prediction models can be constructed to reveal the key factors affecting the performance of SAA catalysts. For example, Liu et al. used a regression machine learning algorithm to screen transition metal-doped silver-based SAA catalysts. They showed that the local coordination environment of the active sites significantly affected their catalytic activity for NO reduction reactions.⁵⁴ Ren and co-workers proposed a “single-atom saturation” model. They used it to quantitatively evaluate the adsorption strength of different intermediates on SAAs, which revealed the structure–activity relationships of SAA catalysts in various reactions.⁵⁵ This data-driven research strategy promotes the screening of high-performance SAA catalysts and provides a theoretical framework for a deeper understanding of the key factors affecting the eNRR.

This study investigated the eNRR performance of the copper-based SAAs by DFT calculations, and all the catalysts were screened based on structural stability, N₂ activation, catalytic activity, and selectivity. Based on the reaction free energy calculated using DFT, the adaptive boosting (Adaboost) ML model was used to construct a correlation between the intrinsic properties of SAAs and the limiting potential (U_L) of the eNRR. In addition, the SHapley Additive exPlanations (SHAP) analysis was used to provide a global interpretation of the ML model. Finally, we used the framework of sure independence screening and sparsifying operator (SISSO) to obtain an optimal descriptor with the four most critical features whose predicted U_L agrees with the results of DFT.

2. Computational method

The 4 × 4 × 4 supercell of the Cu(100) crystal plane was used as the host metal substrate model, and a TM atom was used to replace a surface Cu atom to construct the SAA catalyst model. The Vienna Ab Initio Simulation Package (VASP)⁵⁶ was used to complete all the spin-polarization DFT calculations in this paper. The projector augmented wave (PAW)⁵⁷ method with a cutoff-energy of 450 eV and the Perdew–Burke–Ernzerhof (PBE)⁵⁸ exchange–correlation functional of generalized gradient approximation (GGA) were used to describe the ion–electron interaction and the exchange–correlation energy, respectively. The conjugate gradient method was applied to all geometry relaxation, and all coordinates of the bottom two layers of atoms were fixed. The convergence threshold for energy and force was set as 10⁻⁵ eV and 0.02 eV Å⁻¹, respectively. Gamma-centered 3 × 3 × 1 and 6 × 6 × 1 k-points based on Monkhorst–Pack k point grids⁵⁹ were applied to structural relaxation and electronic property calculations, respectively. Grimme's DFT-D3 strategy⁶⁰ was used to correct the van der Waals interaction. For all spin-polarization DFT calculations, the solvent environment was simulated using the VASPsol model,⁶¹ where the bulk dielectric constant of the solvent was set to 80. The local electronic interactions for all 3d TMs were better described using the DFT+U method, and the corresponding U values are listed in Table S1.⁶² The density of states (DOS) was described more precisely using the tetrahedron method with Blöchl corrections and dispersion of 2000 points. A vacuum layer of 15 Å was set up in the z-direction to avoid periodic interactions and eliminate steric effects. The integrated crystal orbital Hamilton population (ICOHP) analysis of the bonding and antibonding orbitals was completed using Lobster software.⁶³ More detailed computational parameters and settings can be found in the SI. In this study, only systems with the lowest energy were considered in the following calculations and discussions.

The stability of SAA catalysts was assessed using the binding energy (E_bin) of the TM atom to the Cu(100) substrate, which was defined as:

E_bin = E_TM@Cu(100) − E_{Cu(100)_defective} − E_{TM_single} (1)

where E_TM@Cu(100), E_{Cu(100)_defective}, and E_{TM_single} denote the energies of SAA catalysts, Cu(100) surfaces with individual atomic defects, and individual TM atom, respectively. A more negative E_bin value indicates stronger binding of the TM atom to the Cu(100) substrate, i.e., a more stable SAA catalyst.

The ability of the TM atom to aggregate on the surface of Cu(100) can be assessed by calculating the cohesion energy (E_c) of the TM atom and the formation energy (E_f) of SAA catalysts. The more negative value of E_f indicates that the TM atom is more inclined to have been distributed in an atomically monodispersed mode on Cu(100).

E_f = E_bin − E_c (2)

E_c = E_{TM_bulk}/n − E_{TM_single} (3)

where E_{TM_bulk} and n denote the energies of the corresponding unit cell of TM and the number of atoms in the unit cell, respectively.

The adsorption free energy (ΔG_ads) of the intermediates on SAAs during the eNRR process was calculated using the following equation:

ΔG_ads = G_{adsorbate@SAAs} − G_SAAs − G_adsorbate (4)

where G_SAAs and G_{adsorbate@SAAs} represent the energy of SAA catalysts before and after intermediate adsorption, respectively, while G_adsorbate represents the energy of the intermediate.

The Gibbs free energy of reaction (ΔG) was calculated for each elementary step of the eNRR according to the CHE model, which is defined as:

ΔG = ΔE + ΔE_ZPE − TΔS + ΔG_pH (5)

where ΔE is the total energy difference before and after the adsorption of the reaction intermediate in each elementary step as calculated by DFT. ΔE_ZPE and TΔS are the zero-point energy correction and entropy change at 298.15 K, which can be obtained using the VASPKIT software package.⁶⁴ ΔG_pH represents the corrected value of the free energy in different pH environments and is obtained using the equation ΔG_pH = k_BT × pH × ln 10, where k_B is the Boltzmann constant and pH is set to zero for this study.

To evaluate the activity of SAAs in catalyzing the conversion of N₂ to NH₃, the limiting potential (U_L) was calculated:

U_L = −ΔG_max/e (6)

where ΔG_max is the free energy change of the potential-determining step (PDS) of the reaction process calculated directly by DFT.

This study employs the computational hydrogen electrode (CHE) model to evaluate the thermodynamic trends of the eNRR on the TM@Cu(100) SAA catalysts. We explicitly recognize that the CHE model is a thermodynamic approximation that implicitly assumes a linear relationship between the kinetic potential barriers of each elementary step and the change in reaction free energy (i.e., adhering to the BEP relationship). However, it must be noted that completely disregarding specific kinetic barriers may compromise prediction accuracy in some instances. Furthermore, the kinetics of competing hydrogen evolution reactions may also exert complex effects on the final selectivity. Nevertheless, the CHE model is a core theoretical tool in electrocatalytic thermodynamics research. Converting electrode potential (U) into hydrogen chemical potential enables precise correction of reaction free energy. Obtaining reaction free energy (ΔG) via the CHE model allows assessment of reaction feasibility. Consequently, this method efficiently enables reliable preliminary screening of large numbers of catalyst materials, which has been widely validated in previous research in this field.^51,54,65

3. Results and discussion

3.1. Chemical stability of SAA catalysts

As shown in Fig. 1a, TM atoms doped with Cu(100) formed SAAs (denoted as TM@Cu(100)) and pure Cu(100) catalysts were used as models to probe eNRR (DFT-optimized structures as in Fig. S1 and S2). A four-step screening framework was constructed to obtain SAA catalysts with high stability and eNRR activity, as shown in Fig. 1b. First, the E_bin values of the TM@Cu(100) system were calculated using DFT. The catalyst model with E_bin < 0 eV was selected for subsequent studies, which indicated that the system was thermodynamically stable. Second, the adsorption free energy of the N₂ molecule on the TM@Cu(100) surface, ΔG_*N2 < 0 eV, indicates that N₂ forms a stable chemisorption configuration, which is a critical initial step for the eNRR. Third, the TM@Cu(100) catalyst must exhibit high eNRR activity, i.e., a low limiting potential U_L. Finally, the adsorption free energy of H on the TM surface, ΔG_*H > ΔG_*N2, implies that the catalyst has high eNRR selectivity.

Fig. 1 (a) Structural model of TM@Cu(100) SAAs. (b) DFT-based four-step screening strategy. (c) The binding energy between TM atoms and the Cu(100) surface.

As a key basis for maintaining the catalytic efficiency of the eNRR, the chemical stability of TM@Cu(100) SAA catalysts was first investigated. In this work, the stability of TM@Cu(100) was evaluated using the binding energies E_bin of TM atoms and the Cu(100) surface in SAA catalysts, which are plotted in Fig. 1c. It can be seen that the E_bin values for all systems are negative, ranging from −1.46 to −7.61 eV, which indicates that the TM atoms and the Cu(100) substrate exhibit a strong binding capacity, which gives the SAA catalysts excellent thermodynamic stability. The Bader charge analysis of Table S2 shows a significant charge transfer between the TM atoms and Cu(100) surface, and all TM atoms except Pd, Ir, Pt, and Au have lost a certain amount of electrons. Moreover, the trends of charge transfer and E_bin are the same (Fig. S3), further emphasizing the strong interaction between the TM atoms and Cu(100) substrate. Subsequently, the formation energy E_f was further calculated to assess the possibility of TM atoms aggregating into clusters (E_f > 0 eV) on the Cu(100) surface. The results in Fig. S3 and Table S2 show that most of the TM-doped Cu(100) systems have negative formation energies, elucidating that TM tends to be distributed as a single atom on the Cu(100) surface. It is worth stating that the TM@Cu(100) models with positive E_f were also used for the subsequent systematic exploration of eNRR activity.

3.2. Activation of the N₂ molecule on the surface of SAAs

As an initial step of the eNRR, the strength of adsorption of N₂ on the surface of the SAA catalysts can significantly affect the NH₃ yield.^39,66 As shown in Fig. 2a and b, N₂ can be chemisorbed on the TM@Cu(100) surface sites (top, bridge, and hollow sites) using both end-on and side-on configurations, which correspond to the four possible reaction pathways of the eNRR (distal, alternating, consecutive, and enzymatic mechanisms).¹⁷ According to the DFT-optimized adsorption configurations (Fig. S5–S9), the N₂ adsorption at bridge and hollow sites spontaneously transforms into the top-site adsorption configuration, suggesting that the N₂ molecule prefers to adsorb at the top site on the surface of the SAA catalysts. Furthermore, whether the SAA catalysts effectively capture and activate the N₂ molecule can be assessed by analyzing the changes in TM–N and N–N bond lengths. If N₂ is chemisorbed on the TM@Cu(100) surface in end-on or side-on configurations and activated, then the TM–N bond length is close to R_TM + R_N,⁶⁷ and the N–N bond length is significantly larger than 1.115 Å (N–N bond length in a free N₂ molecule). From Fig. 2c, it can be found that the N–N bond lengths in N₂ adsorbed on the surfaces of the TM@Cu(100) (TM = Sc, Ti, Zn, Y, Zr, Ag, Cd, Lu, Hf, Au, and Hg) (both end-on and side-on configurations) are very close to 1.115 Å, which indicates that the N₂ is not activated. Combined with Table S3, it can be seen that the N₂ presents a physisorption configuration on the catalyst surface. For other TM@Cu(100) systems, the N–N bond length in the adsorbed *N₂ was significantly larger than that in the free N₂ molecule, and there was a significant electron transfer between the TM atom and *N₂ (Fig. 2c), which indicated that the N₂ molecule was effectively activated. The DFT-calculated adsorption free energies of N₂ (Fig. S10) showed that ΔG_*N2 (end-on) < ΔG_*N2(side-on) for the SAA systems (TM = V, Mn, Fe, Co, Nb, Mo, Tc, Ru, Rh, Ta, W, Re, Os, Ir, and Pt), which elucidates the tendency of N₂ to adsorb on the catalyst surface in an end-on configuration. It is worth stating that although ΔG_*N2(end-on) and ΔG_*N2(side-on) for TM@Cu(100) systems (TM = Cr, Ni, and Pd) are essentially the same, N₂ spontaneously transforms from the side-on configuration to the end-on configuration at this point (Table S3). In addition, as shown in Fig. 2d, ΔG_*N2(end-on) and the adsorption free energy of protons on the surface of SAAs were further compared, and the results showed that N₂ occupied the active sites on all catalysts better than protons.

Fig. 2 (a) The possible adsorption sites. (b) The chemisorption and physisorption configurations of N₂ on TM@Cu(100) SAAs. (c) The N–N bond length and electron accumulation of adsorbed *N₂. (d) The adsorption free energy of N₂ and H.

In general, the adsorbed *N₂ can bind to TM atoms through two interactions, σ-donation (transfer of a lone pair of electrons from N₂ to the d orbitals of TM atoms) and π-back donation (transfer of d electrons back from TM atoms to the antibonding orbitals of *N₂).⁴ The charge density difference analysis (Fig. S11–S13) shows significant charge transfer between the adsorbed *N₂ and TM atoms, indicating a strong interaction between TM and *N₂. Taking Co@Cu(100) as an example, as shown in Fig. 3a, when N₂ is adsorbed on the surface of Co@Cu(100), there is a significant electron accumulation between Co and N atoms and the back donation effect of N on Co atoms leads to a decrease in the strength of the N≡N bond. By comparison with the pure Cu(100) system, the project density of states (PDOS) in Fig. 3b shows that the d orbitals of the Co atom and the 2p orbitals of the N atom are significantly hybridized near the Fermi energy level and the π* orbitals of the adsorbed *N₂ are partially occupied. This result indicates that there is no significant π-back donation effect in the pure Cu(100) system, which explains the relatively stronger adsorption of N₂ on the Co@Cu(100) surface (the PDOS of other SAAs with adsorbed *N₂ is shown in Fig. S14–S16). Furthermore, the downshift of the antibonding state of *N₂ adsorbed on the pure Cu(100) surface compared to Co@Cu(100) SAAs, as well as having a more negative ICOHP (Fig. 3c), further confirms the presence of a more stable N≡N bond.

Fig. 3 (a) Charge density difference of adsorbed *N₂ and Co@Cu(100) SAAs. The project density of states (b) and crystal orbital Hamilton populations of (c) Cu(100) and Co@Cu(100) with adsorbed *N₂.

3.3. Reaction pathways and activity trends of the eNRR

According to the Sabatier principle, the strength of N₂ adsorption on the catalyst surface significantly affects the reaction rate of the eNRR and the desorption process of product NH₃.⁶⁸ Therefore, despite N₂ exhibiting stronger adsorption on the TM@Cu(100) surface in the end-on configuration, we still considered both end-on and side-on adsorption configurations of N₂ in subsequent studies. In addition, pure Cu(100) and 29 types of TM@Cu(100) catalysts were used in the subsequent DFT calculations. As shown in Fig. 4, there are four possible reaction pathways for electrochemical reduction of the N₂ molecule to NH₃, which are distal, alternating, consecutive, and enzymatic mechanisms.

Fig. 4 Schematic illustration of the possible eNRR mechanism on TM@Cu(100) SAAs, including alternating, distal, consecutive and enzymatic pathways.

In the eNRR, *N₂ adsorbed on the active site undergoes six electron transfer (H⁺ + e⁻) steps. As shown in Fig. S17–S21, we systematically investigated the free energy changes of the eNRR (the adsorption configuration of each intermediate is shown in Fig. S22 and S33). From Table S4, it can be found that for the first hydrogenation step, SAA catalysts other than some TM@Cu(100) (TM = Cr, Ag, Cd, Au, and Hg) exhibit lower reaction free energy than the pure Cu(100) system. Among them, most TM@Cu(100) systems showed a decrease of ΔG(*N₂ → *N₂H) values more than 0.50 eV. However, for the Sc and Y@Cu(100) systems, N₂ could not be completely reduced to NH₃. The above results suggest that the 22 types of SAA catalysts may possess superior eNRR performance than the Cu(100) system. Moreover, the connection between the physicochemical properties of *N₂ adsorbed on the surface of SAAs in the end-on configuration and ΔG(*N₂), ΔG(*N₂ → *N₂H) and ΔG(*NH₂ → *NH₃) was further investigated. Fig. 5 shows a clear linear relationship between reaction free energy and N≡N bond length and electron accumulation (δ_N2) in adsorbed *N₂. With the increase of bond length and δ_N2, ΔG(*N₂) and ΔG(*N₂ → *N₂H) become more negative, while ΔG(*NH₂ → *NH₃) gradually increases. The results indicate that the stronger the adsorption of TM@Cu(100) catalysts on the N₂ molecule, the more electrons are transferred from the TM atom to adsorbed *N₂, and the higher the degree of activation of the N≡N bond. Accordingly, carrying out the first hydrogenation step is easier, but carrying out the last becomes more challenging.⁶⁹ As shown in Fig. 6a and b, N₂ has a more negative adsorption energy on the W@Cu(100) surface than the pure Cu(100) system, and the reaction free energy for the first hydrogenation step is significantly lower. However, in contrast, the reaction free energy of the last hydrogenation step on the Cu(100) surface is substantially lower than that of W@Cu(100) SAAs (−0.44 vs. 0.18 eV).

Fig. 5 The relationships between adsorption free energies and (a) N–N bond length and (b) electron accumulation of adsorbed *N₂.

Fig. 6 The free-energy diagram of (a) W@Cu(100) and (b) Cu(100). (c) The limiting potentials (U_L) of TM@Cu(100) SAAs and Cu(100). (d) Comparison of calculated ΔG_ads(*N₂) and ΔG_ads(*H). The adsorption structure of the H atom is shown in Fig. S36.

Subsequently, the free energy changes of the complete eNRR process based on four reaction mechanisms were further investigated. Among them, TM@Cu(100) catalysts could completely reduce N₂ to NH₃, except for Sc and Y doping. Among the many reports, for the distal and consecutive mechanisms, the N–N bond breaking occurred mainly in the third hydrogenation step (*NNH₂ or *N*NH₂ → *N + NH₃), whereas for the alternating and consecutive mechanisms, N–N bond breaking occurs at the fifth hydrogenation step (*NH₂NH₂ or *NH₂*NH₂ → *NH₂ + NH₃).^70,71 However, although the production of the first NH₃ molecule occurs at different steps and the reaction free energies of the eNRR show different variations, the PDSs of SAA catalysts excluding Ta/W@Cu(100) are all the first hydrogenation steps. That is to say, SAAs can significantly change the reaction free energy but have less effect on the PDSs. As shown in Fig. 6c, different catalysts have significantly different limiting potentials (U_L) ranging from −0.03 to −2.81 V. Among them, the U_L values of some TM@Cu(100) catalysts (TM = V, Nb, Mo, Tc, Ru, Rh, Hf, Ta, W, Re, Os, and Ir) are significantly lower than the U_L value (−0.98 V) of Ru(0001), which has the highest eNRR activity among bulk metal materials.⁷² In addition, the eNRR activities of V, Nb, Mo, Ta, and W@Cu(100) SAA catalysts are higher than those of most of the SACs reported in experiments^73,74 and theoretical calculations,^70,75 which have a limiting potential of −0.30 V. The above results indicate that SAAs can also significantly affect the limiting potential of the eNRR. Furthermore, although all SAAs can chemisorb N₂, the end-on adsorption configuration of *N₂ has a much lower free energy (Fig. S10). However, as can be seen in Fig. 6c, for catalysts where PDSs are the first hydrogenation step, the adsorbed *N₂ on half of the systems is present in a side-on conformation, whereas *N₂ on the other systems is adsorbed in an end-on conformation, which is also in good agreement with the Sabatier principle.

One of the most critical features of the highly selective TM@Cu(100) SAA catalysts is their ability to inhibit the HER significantly.⁷⁶ Therefore, to evaluate the selectivity of the catalysts for the eNRR, we further compared ΔG_ads(*N₂) and ΔG_ads(*H) for SAAs other than Sc/Y@Cu(100). If ΔG_ads(*N₂)–ΔG_ads(*H) > 0 eV, then the eNRR is more favorable for NH₃ selectivity.⁷³ As shown in Fig. 6d, all systems except Lu@Cu(100) tended to adsorb N₂ rather than hydrogen and thus exhibited excellent eNRR selectivity.

3.4. Machine learning model

The traditional d-band theory and CHE model cannot accurately describe the eNRR performance of TM-based SAA catalysts, and the construction of high-precision ML models can help the component optimization and rational design of SAA catalysts.^51,77,78 Here, we further analyze the relationship between the intrinsic properties of SAAs and the limiting potential (U_L) of the eNRR based on high-throughput DFT calculations. Of these, 20 features (Table S5), including the elemental properties of the TM (van der Waals radius, electron affinity, and electronegativity), the DFT-calculated electronic structures of the SAAs (d-band centers, filling of a d-band, and electron accumulation of *N₂), and structural features of SAAs with adsorbed *N₂ (N–N and TM–N bond length) were used to predict the U_L (labeling of ML model) of the eNRR.^79,80 However, by calculating the Pearson correlation coefficient (Fig. S34) between all the features, it was found that most of the features were strongly correlated with each other (|p| > 0.80).^54,81 Therefore, we reduced the number of features to eight to avoid introducing redundancy problems during the ML training process (Table S6). As shown in Fig. 7a, the p-values between all the features are significantly less than 0.80, indicating that the selected features are reasonable and help improve the ML model's accuracy.

Fig. 7 (a) The Pearson correlation coefficient between features. (b) Comparison between DFT and Adaboost-predicted U_L. (c) The ranking of important features of the Adaboost model. (d) The gain value and SHAP value of features. (e) The SHAP summary plot of the Adaboost model. (f) Comparison between DFT and SISSO-predicted U_L.

Subsequently, common ML algorithms such as random forest, adaptive boosting, and extreme gradient boosting models were used to predict the U_L of the eNRR.^82–84 As can be seen in Fig. 7b and Fig. S35, the Adaboost model showed the best performance with a training R² of 0.97 and a testing R² of 0.86. Feature importance analysis (Fig. 7c) shows that electron accumulation of *N₂ (δ_N2) and the valence electron number of TM (N) are the two most critical features, with importance scores of 55.5% and 16.9%, respectively. This confirms that the valence electron distribution of TM and its interaction with adsorbed *N₂ can significantly affect the U_L. We further calculated the SHAP values and used them to interpret the contribution of each feature to the prediction results.^85,86 As shown in Fig. 7d, the SHAP values of δ_N2 and N are 0.60 and 0.21, respectively, which are significantly larger than those of the other features, ranging from 0.02 to 0.06. The results show that δ_N2 and N are the essential features, while TM–N bond length (d_M–N) and electron affinity (EA) are less important, consistent with the gain values of the Adaboost model. Moreover, the global feature interpretation is provided in Fig. 7e for the SHAP summary of each feature, and it can be found that the eigenvalues of δ_N2 and N are positively and negatively correlated with the corresponding SHAP values, respectively. However, for the other six features, there is no significant correlation between their eigenvalues and SHAP values. Therefore, δ_N2 and N play a key role in rationally regulating the eNRR performance of SAA catalysts.

Introducing the SISSO method makes it possible to identify the few most critical features and obtain their mathematical expressions, facilitating the rapid screening of catalysts.^87–89 Using the six most important features identified by the Adaboost model and the limiting potentials as inputs to the SISSO algorithm, this study ended up with a 2-dimensional descriptor of complexity 4, which is denoted as:

As shown in Fig. 7f, the R² and RMSE of eqn (7) are 0.89 and 0.18, respectively, indicating that it can predict the U_L of the catalyst well. The four key features identified by the SISSO algorithm—the electron accumulation of *N₂ (δ_N2), the valence electron number of transition metal (N), TM–N bond length (d_M–N), and electron affinity (EA)—are not independent statistical metrics. Together, they form a physical picture describing the electronic structure of the TM site and its interaction with reaction intermediates. (1) δ_N2 directly quantifies the additional electrons acquired by the N₂ molecule upon adsorption at the TM active site. It serves as a direct measure of N≡N triple bond activation. Significant electron accumulation fills the anti-bonding π orbitals of *N₂, effectively weakening the N≡N triple bond and laying the groundwork for subsequent hydrogenation steps. (2) N represents the intrinsic electronic property of the TM atom, determining its capacity to donate d electrons and form feedback bonds with reactants like *N₂. This is crucial for regulating the intrinsic electron-donating strength of the TM site. (3) d_M–N reflects the strength of interaction between the TM active site and adsorbed nitrogen-containing intermediates. According to quantum mechanical theory, shorter bond lengths typically indicate stronger chemical bonds for a given system's ground state. Thus, d_M–N is closely correlated with the adsorption energy of critical reaction intermediates. (4) EA measures the energy released when an atom gains an electron, reflecting an element's tendency to accept electrons. In electrocatalytic environments, TM active sites with higher EA may be more favorable for receiving electrons from the electrode surface under applied potentials. These sites can then efficiently transfer electrons to adsorbed *N₂ or reaction intermediates, facilitating critical electron transfer steps. Overall, the Adaboost and SISSO ML models are helpful to understand the intrinsic connection between the physicochemical properties of SAA catalysts and their eNRR performances and provide a rational design of efficient eNRR SAA catalysts with an essential basis.

It should be explicitly noted that the DFT calculations and machine learning descriptor construction in this study were performed based on a low coverage model. This assumption was made to ensure the feasibility of large-scale computational screening and is widely adopted in related studies within this field.^50,90,91 We fully recognize that in actual electrochemical environments, reaction intermediates (such as N₂, *NNH, *NH₂, etc.), along with numerous solvated water molecules and electrolyte ions, form higher coverage on electrode surfaces. This complex interfacial environment may influence catalytic performance through multiple pathways. (1) Direct repulsive or attractive forces between intermediates may alter their adsorption energies. For instance, repulsive interactions between polar intermediates like *NH₂ may weaken their adsorption, altering reaction thermodynamics. (2) High coverage may shift the optimal reaction pathway. For example, thermodynamically unfavorable pathways involving association mechanisms at low coverage may become feasible at high coverage due to steric hindrance or changes in the electronic structure. (3) The impact of coverage on competitive hydrogen evolution reactions may be particularly pronounced. Competition between HER intermediates and eNRR intermediates for surface sites directly determines Faraday efficiency. At high H coverage, eNRR active sites may become blocked, leading to reduced selectivity. Despite these potential complexities, our screening results under low coverage conditions retain a significant guiding value. First, they successfully identify candidate materials exhibiting optimal binding strengths for key eNRR intermediates in their intrinsic properties—a core foundation for high-performance catalysts. Second, the key physical descriptors we identified (e.g., δ_N2, N, and d_M–N) provide crucial electronic and geometric structural insights for understanding activity trends. Finally, the most promising SAA catalysts screened in this work (e.g., Mo and W@Cu(100)) provide a theoretical foundation for subsequent studies, including high-coverage, explicitly solvated molecular dynamics simulations.

4. Conclusions

The eNRR performance of TM@Cu(100) SAA catalysts was intensively investigated by combining DFT calculations and ML methods. The calculated reaction free energy results show that the SAAs can significantly reduce the free energy changes of the first hydrogenation steps (PDSs) compared to the pure Cu(100) system. In addition, five types of TM@Cu SAAs (TM = V, Nb, Mo, Ta, and W) are recognized as excellent eNRR catalysts with strong structural stability, good selectivity, and high eNRR activity with a low limit potential (U_L > −0.30 V). Subsequently, we applied the Adaboost ML algorithm to construct a well connection between the physicochemical properties of SAAs and the U_L of the eNRR. We used the SHAP analysis to reveal the significant effects of δ_N2 and N features on the activity of SAAs. Finally, based on the SISSO algorithm, the importance of the features determined by the ML model was further confirmed, and an optimal feature equation was obtained. This study not only evaluated the eNRR potential of copper-based SAA catalysts, but also revealed the relationship between the intrinsic features of SAAs and their catalytic activities. We hope that our work can contribute to the rational design of SAA catalysts.

Author contributions

H. Z. Liu: conception, investigation, performing calculations, and writing – original draft. Y.-G. Wang: supervision, conception, and writing – review and editing. All authors contributed to the scientific discussion and preparation of the manuscript.

Conflicts of interest

There are no conflicts to declare.

Data availability

The data supporting this article have been included as supplementary information (SI). Supplementary information: the atomic structure and electronic property analysis of catalysts with or without an absorbate. See DOI: https://doi.org/10.1039/d5cp03528g.

Acknowledgements

This work was financially supported by the National Key R&D Program of China (No. 2022YFA1503102), the NSFC Center for Single-Atom Catalysis (No. 22388102), the Shenzhen Science and Technology Program (No. JCYJ20220818100410023) and the Guangdong Provincial Key Laboratory of Catalysis (No. 2020B121201002). The computational resource is supported by the Center for Computational Science and Engineering at SUSTech and the CHEM high-performance supercomputer cluster (CHEM-HPC) located in the Department of Chemistry, SUSTech.

References

1. D. E. Canfield, A. N. Glazer and P. G. Falkowski, Science, 2010, 330, 192–196.
2. S. L. Foster, S. I. P. Bakovic, R. D. Duda, S. Maheshwari, R. D. Milton, S. D. Minteer, M. J. Janik, J. N. Renner and L. F. Greenlee, Nat. Catal., 2018, 1, 490–500.
3. G. Qing, R. Ghazfar, S. T. Jackowski, F. Habibzadeh, M. M. Ashtiani, C.-P. Chen, M. R. Smith, III and T. W. Hamann, Chem. Rev., 2020, 120, 5437–5516.
4. M.-A. Légaré, G. Bélanger-Chabot, R. D. Dewhurst, E. Welz, I. Krummenacher, B. Engels and H. Braunschweig, Science, 2018, 359, 896–900.
5. W. Guo, K. Zhang, Z. Liang, R. Zou and Q. Xu, Chem. Soc. Rev., 2019, 48, 5658–5716.
6. A. J. Martín, T. Shinagawa and J. Pérez-Ramírez, Chem, 2019, 5, 263–283.
7. T. Wang and F. Abild-Pedersen, Proc. Natl. Acad. Sci. U. S. A., 2021, 118, e2106527118.
8. J. W. Erisman, M. A. Sutton, J. Galloway, Z. Klimont and W. Winiwarter, Nat. Geosci., 2008, 1, 636–639.
9. C. Smith, A. K. Hill and L. Torrente-Murciano, Energy Environ. Sci., 2020, 13, 331–344.
10. B. H. R. Suryanto, H.-L. Du, D. Wang, J. Chen, A. N. Simonov and D. R. MacFarlane, Nat. Catal., 2019, 2, 290–296.
11. J. G. Chen, R. M. Crooks, L. C. Seefeldt, K. L. Bren, R. M. Bullock, M. Y. Darensbourg, P. L. Holland, B. Hoffman, M. J. Janik, A. K. Jones, M. G. Kanatzidis, P. King, K. M. Lancaster, S. V. Lymar, P. Pfromm, W. F. Schneider and R. R. Schrock, Science, 2018, 360, eaar6611.
12. A. J. Martín and J. Pérez-Ramírez, Joule, 2019, 3, 2602–2621.
13. Z. Chen, C. Liu, L. Sun and T. Wang, ACS Catal., 2022, 12, 8936–8975.
14. J. Wu, S. Wang, R. Ji, D. Kai, J. Kong, S. Liu, W. Thitsartarn, B. H. Tan, M. H. Chua, J. Xu, X. J. Loh, Q. Yan and Q. Zhu, ACS Nano, 2024, 18, 20934–20956.
15. Y. Wu, M. Zhang, Z. Zhao, W. Song, J. Yan, J. Wen, A. Mahmood Ali, X. Ge and H. Zhang, Energy Fuels, 2025, 39, 1507–1540.
16. Z. Yin, Z. Yang, Y. Bian, H. Zhao, B. Chen, Y. Liu, Y. Wang, Y. Deng and H. Wang, J. Mater. Chem. A, 2025, 13, 9643–9650.
17. X. Guo, J. Gu, S. Lin, S. Zhang, Z. Chen and S. Huang, J. Am. Chem. Soc., 2020, 142, 5709–5721.
18. Y. Tang, Y. Wang, X. Cheng and H. Zhang, J. Phys. Chem. C, 2023, 127, 24574–24582.
19. M. Ismael and M. Wark, Appl. Mater. Today, 2024, 39, 102253.
20. R. T. Hannagan, G. Giannakakis, M. Flytzani-Stephanopoulos and E. C. H. Sykes, Chem. Rev., 2020, 120, 12044–12088.
21. T. Zhang, A. G. Walsh, J. Yu and P. Zhang, Chem. Soc. Rev., 2021, 50, 569–588.
22. T. Shen, S. Wang, T. Zhao, Y. Hu and D. Wang, Adv. Energy Mater., 2022, 12, 2201823.
23. R. J. Bunting, F. Wodaczek, T. Torabi and B. Cheng, J. Am. Chem. Soc., 2023, 145, 14894–14902.
24. Y. Zhang, G. Zhang, Q. Hou, S. Zhao, S. Wang, E. Xu, L. Wang, X. Zhang, F. Li, Y. Yang and M. Wei, ACS Catal., 2025, 15, 1206–1216.
25. J. Liu, M. B. Uhlman, M. M. Montemore, A. Trimpalis, G. Giannakakis, J. Shan, S. Cao, R. T. Hannagan, E. C. H. Sykes and M. Flytzani-Stephanopoulos, ACS Catal., 2019, 9, 8757–8765.
26. J. D. Lee, J. B. Miller, A. V. Shneidman, L. Sun, J. F. Weaver, J. Aizenberg, J. Biener, J. A. Boscoboinik, A. C. Foucher, A. I. Frenkel, J. E. S. van der Hoeven, B. Kozinsky, N. Marcella, M. M. Montemore, H. T. Ngan, C. R. O’Connor, C. J. Owen, D. J. Stacchiola, E. A. Stach, R. J. Madix, P. Sautet and C. M. Friend, Chem. Rev., 2022, 122, 8758–8808.
27. H. Yin, Y. Peng and J. Li, Environ. Sci. Technol., 2023, 57, 3134–3144.
28. M. R. Hossain and M. Trenary, J. Phys. Chem. C, 2024, 128, 19204–19209.
29. H. T. Ngan and P. Sautet, J. Am. Chem. Soc., 2024, 146, 2556–2567.
30. J. Cai, Y. Wei, A. Cao, J. Huang, Z. Jiang, S. Lu and S.-Q. Zang, Appl. Catal., B, 2022, 316, 121683.
31. C. Du, S. Lu, J.-A. Wang, X. Wang, M. Wang, H. M. Fruehwald, L. Wang, B. Zhang, T. Guo, J. P. Mills, W. Wei, Z. Chen, Y. Teng, J. Zhang, C.-J. Sun, H. Zhou, R. D. L. Smith, B. Kendall, G. Henkelman and Y. A. Wu, ACS Catal., 2023, 13, 10560–10569.
32. R. Réocreux and M. Stamatakis, Acc. Chem. Res., 2022, 55, 87–97.
33. P. Agharezaei, N. G. Tomohiro, H. Kobayashi, H. Schlenz, M. Yamauchi and K. K. Ghuman, J. Phys. Chem. C, 2024, 128, 3703–3717.
34. S. Li, C.-W. Wang, X. Zhao, J.-S. Dang and J. Li, J. Phys. Chem. Lett., 2024, 15, 5088–5095.
35. T. Dai, Z. Wang, X. Lang and Q. Jiang, J. Mater. Chem. A, 2022, 10, 16900–16907.
36. Y. Da, R. Jiang, Z. Tian, X. Han, W. Chen and W. Hu, SmartMart, 2023, 4, e1136.
37. Y. Yao, H. Geng, W. Duan, L. Zhang, F. Xu and H. Shang, Nano Res., 2025, 18, 94907313.
38. L. M. Tsiverioti, L. Kavalsky and V. Viswanathan, J. Phys. Chem. C, 2022, 126, 12994–13003.
39. Y. Sun, X. Li, Z. Wang, L. Jiang, B. Mei, W. Fan, J. Wang, J. Zhu and J.-M. Lee, J. Am. Chem. Soc., 2024, 146, 7752–7762.
40. J. K. Nørskov, J. Rossmeisl, A. Logadottir, L. Lindqvist, J. R. Kitchin, T. Bligaard and H. Jónsson, J. Phys. Chem. B, 2004, 108, 17886–17892.
41. C. Fu, L. Luo, L. Yang, S. Shen, G. Wei, J. Yin and J. Zhang, J. Phys. Chem. C, 2021, 125, 17051–17057.
42. D. Gao, D. Yi, C. Sun, Y. Yang and X. Wang, ACS Appl. Mater. Interfaces, 2022, 14, 52806–52814.
43. A. Biswas, S. Nandi, N. Kamboj, J. Pan, A. Bhowmik and R. S. Dey, ACS Nano, 2021, 15, 20364–20376.
44. C. Vogt and B. M. Weckhuysen, Nat. Rev. Chem., 2022, 6, 89–111.
45. X. Li, P. Shen, Y. Luo, Y. Li, Y. Guo, H. Zhang and K. Chu, Angew. Chem., Int. Ed., 2022, 61, e202205923.
46. V. Fung, G. Hu and B. Sumpter, J. Mater. Chem. A, 2020, 8, 6057–6066.
47. J. Long, S. Chen, Y. Zhang, C. Guo, X. Fu, D. Deng and J. Xiao, Angew. Chem., Int. Ed., 2020, 59, 9711–9718.
48. G. Zheng, Y. Li, X. Qian, G. Yao, Z. Tian, X. Zhang and L. Chen, ACS Appl. Mater. Interfaces, 2021, 13, 16336–16344.
49. B. Dou, Z. Zhu, E. Merkurjev, L. Ke, L. Chen, J. Jiang, Y. Zhu, J. Liu, B. Zhang and G.-W. Wei, Chem. Rev., 2023, 123, 8736–8780.
50. R. Ding, J. Chen, Y. Chen, J. Liu, Y. Bando and X. Wang, Chem. Soc. Rev., 2024, 53, 11390–11461.
51. J. Sun, R. Tu, Y. Xu, H. Yang, T. Yu, D. Zhai, X. Ci and W. Deng, Nat. Commun., 2024, 15, 6036.
52. P. Zhao, H. Wang, Z. Huang, B. Fan, Y. Liu, Z. Yin, Y. Wang, X. Zheng, Y. Deng and X. Fan, Appl. Surf. Sci., 2023, 631, 157480.
53. Z.-X. Yin, Y.-D. Li, Y.-H. Ye, Y. Liu, M.-F. Li, Z.-J. Yang, X.-R. Zheng, H.-Z. Wang, Y. Wang and Y.-D. Deng, Rare Met., 2024, 43, 5781–5791.
54. J. Liu, S. Wang, Y. Tian, H. Guo, X. Chen, W. Lei, Y. Yu and C. Wang, Angew. Chem., Int. Ed., 2025, 64, e202414314.
55. C. Ren, Y. Cui, Q. Li, C. Ling and J. Wang, J. Am. Chem. Soc., 2025, 147(16), 13610–13617.
56. G. Kresse and J. Hafner, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 558–561.
57. G. Kresse and D. Joubert, Phys. Rev. B: Condens. Matter Mater. Phys., 1999, 59, 1758–1775.
58. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865–3868.
59. S. Froyen, Phys. Rev. B: Condens. Matter Mater. Phys., 1989, 39, 3168–3172.
60. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104.
61. K. Mathew, R. Sundararaman, K. Letchworth-Weaver, T. A. Arias and R. G. Hennig, J. Chem. Phys., 2014, 140, 084106.
62. H. Xu, D. Cheng, D. Cao and X. C. Zeng, Nat. Catal., 2024, 7, 219.
63. V. L. Deringer, A. L. Tchougréeff and R. Dronskowski, J. Phys. Chem. A, 2011, 115, 5461–5466.
64. V. Wang, N. Xu, J.-C. Liu, G. Tang and W.-T. Geng, Comput. Phys. Commun., 2021, 267, 108033.
65. A. A. Peterson, F. Abild-Pedersen, F. Studt, J. Rossmeisl and J. K. Nørskov, Energy Environ. Sci., 2010, 3, 1311–1315.
66. Z. Xue, X. Zhang, J. Qin and R. Liu, Nano Energy, 2021, 80, 105527.
67. P. Pyykkö, J. Phys. Chem. A, 2015, 119, 2326–2337.
68. M. Che, Catal. Today, 2013, 218–219, 162–171.
69. L. Xu, Y. Huang, H. Lin, X. Wei and F. Ma, Nano Res., 2025, 18, 94907289.
70. L. Li, X. Wang, H. Guo, G. Yao, H. Yu, Z. Tian, B. Li and L. Chen, Small Methods, 2019, 3, 1900337.
71. L. Cai, N. Zhang, B. Qiu and Y. Chai, ACS Appl. Mater. Interfaces, 2020, 12, 20448–20455.
72. X. Guo, S. Lin, J. Gu, S. Zhang, Z. Chen and S. Huang, Adv. Funct. Mater., 2021, 31, 2008056.
73. H. Tao, C. Choi, L.-X. Ding, Z. Jiang, Z. Han, M. Jia, Q. Fan, Y. Gao, H. Wang, A. W. Robertson, S. Hong, Y. Jung, S. Liu and Z. Sun, Chem, 2019, 5, 204–214.
74. Y. Chen, R. Guo, X. Peng, X. Wang, X. Liu, J. Ren, J. He, L. Zhuo, J. Sun, Y. Liu, Y. Wu and J. Luo, ACS Nano, 2020, 14, 6938–6946.
75. C. Choi, S. Back, N.-Y. Kim, J. Lim, Y.-H. Kim and Y. Jung, ACS Catal., 2018, 8, 7517–7525.
76. A. R. Singh, B. A. Rohr, M. J. Statt, J. A. Schwalbe, M. Cargnello and J. K. Nørskov, ACS Catal., 2019, 9, 8316–8324.
77. H.-L. Jin, Q.-N. Li, Y.-Y. Tian, S.-A. Wang, X. Chen, J.-Y. Liu and C.-H. Wang, Rare Met., 2024, 43, 5813–5822.
78. X. Duan, Y. Li, J. Zhao, M. Zhang, X. Wang, L. Zhang, X. Ma, Y. Qu and P. Zhang, J. Am. Chem. Soc., 2025, 147, 651–661.
79. D. Deng, B. Song, C. Li and L.-M. Yang, J. Phys. Chem. C, 2022, 126, 20816–20830.
80. T.-T. Shi, G.-Y. Liu and Z.-X. Chen, J. Phys. Chem. C, 2023, 127, 9573–9583.
81. H. Mai, T. C. Le, D. Chen, D. A. Winkler and R. A. Caruso, Chem. Rev., 2022, 122, 13478–13515.
82. S. Lin, H. Xu, Y. Wang, X. C. Zeng and Z. Chen, J. Mater. Chem. A, 2020, 8, 5663–5670.
83. J. Gong, S. Chu, R. K. Mehta and A. J. H. McGaughey, npj Comput. Mater., 2022, 8, 140.
84. Z. Li, C. Zhao, H. Wang, Y. Ding, Y. Chen, P. Schwaller, K. Yang, C. Hua and Y. He, Proc. Natl. Acad. Sci. U.S.A., 2024, 121, e2320232121.
85. S. Lu, Q. Zhou, Y. Guo and J. Wang, Chem, 2022, 8, 769–783.
86. J. Sun, A. Chen, J. Guan, Y. Han, Y. Liu, X. Niu, M. He, L. Shi, J. Wang and X. Zhang, Energy Environ. Mater., 2024, 7, e12693.
87. R. Ouyang, S. Curtarolo, E. Ahmetcik, M. Scheffler and L. M. Ghiringhelli, Phys. Rev. Mater., 2018, 2, 083802.
88. J. Wang, H. Xie, Y. Wang and R. Ouyang, J. Am. Chem. Soc., 2023, 145, 11457–11465.
89. R. Miyazaki, K. S. Belthle, H. Tüysüz, L. Foppa and M. Scheffler, J. Am. Chem. Soc., 2024, 146, 5433–5444.
90. C. Ren, S. Lu, Y. Wu, Y. Ouyang, Y. Zhang, Q. Li, C. Ling and J. Wang, J. Am. Chem. Soc., 2022, 144, 12874–12883.
91. S. Picello, E. Inico, C. Saetta, G. Di Liberto and G. Pacchioni, ACS Catal., 2025, 15, 11232–11242.

---