Bioinformatics analysis of Lavandula angustifolia NAD(P)H-quinone oxidoreductases (LaNDHs) was conducted using data obtained from the UniProt database ( Supplementary Table S1 ). The molecular weights of these enzymes vary from 11.30 kDa to 84.17 kDa ( Table 1 ). The number of amino acids in the LaNDHs proteins ranges from 101 to 739 ( Table 1 ). Their isoelectric points (pI) span from 4.19 to 9.53 ( Table 1 ). The instability index of these enzymes varies between 22.67 and 55.85 ( Table 1 ).

Using the amino acid sequences of LaNDH-2, LaNDH-11, LaNDH-4L1 and LaNDH-4L2 (The reasons for our choice of LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2 can be found in the following content.), we predicted their secondary structures using the PSIPRED (Buchan et al., 2024; Jones, 1999) and NPS@ server (Combet et al., 2000) tools, respectively ( Figure 1 ; Tables 1 and 2 ). The predicted secondary structures of LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2 are predominantly composed of alpha helices, accounting for 59.02%, 62.36%, 72.28%, and 74.26% of the residues, respectively ( Figure 1 ; Table 1 ). Additionally, each protein contains multiple strands and coils ( Figure 1 ; Table 1 ). The number of residues in the helices for LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2 are 301, 227, 73, and 75, respectively ( Figure 1 ; Table 1 ).

The three-dimensional (3D) structures of LaNDHs were predicted using AlphaFold2 (Wayment-Steele et al., 2023; Jumper et al., 2021). AlphaFold2 is a deep learning-based tool known for providing highly accurate and reliable protein structure predictions, which outperform traditional homology modeling techniques. To assess the quality of the predicted structures ( Figures 2 , Supplementary Figure S1 ), we employed the Ramachandran plot to analyze the dihedral angles of the protein backbones. These ensured they fell within acceptable regions, which indicates a valid protein conformation ( Supplementary Figure S2 ; Table 3 ). A high-quality model is expected to have more than 90% of its residues in the most favored regions. In the most favored region, the residual rates of LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2 all exceeded 94%, indicating that these models represent the highest quality structures among these LaNDHs ( Table 3 ). Consequently, we proceeded with further analysis using LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2.

For LaNDH-2, 94.0% of residues were in the most favored region, 6.0% in the additionally allowed region, and none in the generously allowed or disallowed regions ( Table 3 ). For LaNDH-11, 95.0% of residues were in the most favored region, 4.7% in the additionally allowed region, 0.3% in the disallowed region, and none in the generously allowed region ( Table 3 ). For LaNDH-4L1, 94.7% of residues were in the most favored region, 5.3% in the additionally allowed region, and none in the generously allowed or disallowed regions ( Table 3 ). For LaNDH-4L2, 94.7% of residues were in the most favored region, 4.3% in the additionally allowed region, 1.0% in the generously allowed region, and none in the disallowed region ( Table 3 ).

ProSA analysis of the models revealed Z-scores of -6.22, -3.58, -2.47, and -2.59 for LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2, respectively ( Figures 3a , Supplementary Figure S3 ). The overall quality factors of these models were 97.21, 95.66, 96.63, and 96.63, respectively ( Figure 3b ), further confirming the high quality of the predicted structures.

Using the predicted models ( Figure 2 ), we employed the GalaxyWEB program (Ko et al., 2012; Heo et al., 2013, 2016; Seok et al., 2021) to identify the active sites of LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2 ( Figure 4 ). The results revealed that the active site residues of LaNDH-2 include I352, L356, K417, S440, I451, and L454 ( Figures 4a, e ). For LaNDH-4L1, the active site residues were identified as S40, I43, N44, T47, and F48 ( Figures 4b, e ). For LaNDH-4L2, the active site residues include I37, L38, S40, V41, M43, N44, and T47 ( Figures 4c, e ). In the case of LaNDH-11, the active site residues consist of R229, L265, L266, S269, I323, and A324 ( Figures 4d, e ). These residues are highly likely to be involved in the catalytic process, potentially interacting with the substrate side chain atoms to form essential bonds.

To examine the expression profiles of LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2 across various tissues (leaves, stems, flowers, and roots), we conducted real-time quantitative polymerase chain reaction (RT-qPCR). The results indicated that the highest expression levels of LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2 were found in the leaves compared to other tissues ( Figure 5 ). Specifically, the expression of LaNDH-2 was upregulated by 1663.5-fold in leaves, 10.6-fold in flowers, 5.7-fold in stems, and 1.1-fold in roots ( Figure 5 ). LaNDH-11 expression was increased by 560.3-fold in leaves, 4.6-fold in flowers, 2.9-fold in stems, and 1.1-fold in roots ( Figure 5 ). For LaNDH-4L1, expression was upregulated by 388.0-fold in leaves, 7.5-fold in flowers, 6.8-fold in stems, and 1.1-fold in roots ( Figure 5 ). LaNDH-4L2 expression increased by 812.9-fold in leaves, 20.2-fold in flowers, 4.3-fold in stems, and 1.2-fold in roots ( Figure 5 ). These results suggest that LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2 are predominantly expressed in leaf tissue, implying their primary involvement in chloroplast-based photosynthetic processes.

We performed RT-qPCR analysis to assess the expression levels of LaNDHs (LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2) in response to cadmium ion (Cd²⁺), heat, and salt treatments in leaves ( Figure 6 ). The results revealed that the expression of LaNDHs in leaves was positively correlated with increasing Cd²⁺ concentrations ( Figures 6a, d, g, j ). Similarly, LaNDHs expression in leaves increased as the temperature rose from 25°C to 40°C ( Figures 6b, e, h, k ). Additionally, LaNDHs expression in leaves was upregulated with higher salt concentrations ( Figures 6c, f, i, l ). These findings suggest that cadmium ion, heat and salt stress influence the photosynthetic rate in lavender, providing evidence for the association between LaNDHs genes (LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2) and the photosynthetic process.

To establish a comprehensive light-responsive gene expression profile, we evaluated the expression levels of LaNDHs genes (LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2) in leaves under various light conditions (white, red, and blue) using RT-qPCR. The results showed that the expression levels of these genes were highest under blue light compared to other light conditions ( Figure 7 ). Specifically, for LaNDH-2, the expression was highest under blue light (217,898.9-fold), followed by white light (1,663.5-fold), and red light (111.4-fold) ( Figure 7a ). For LaNDH-4L1, the highest expression was observed under blue light (3,251.1-fold), followed by white light (561.6-fold), and red light (176.9-fold) ( Figure 7b ). For LaNDH-4L2, the expression peaked under blue light (1,702.3-fold), followed by white light (388.1-fold), and red light (256.3-fold) ( Figure 7c ). For LaNDH-11, the highest expression was found under blue light (2,786.4-fold), followed by white light (812.4-fold), and red light (345.7-fold) ( Figure 7d ). These findings underscore the significant role of light in regulating the expression of LaNDHs genes (LaNDH-2, LaNDH-11, LaNDH-4L1, and LaNDH-4L2), further supporting the connection between these genes and photosynthesis.

The amino acid sequences of Lavandula angustifolia NAD(P)H-quinone oxidoreductases (LaNDHs) ( Supplementary Table S1 ) were analyzed using ProtParam to predict their chemical properties and physicochemical parameters (Duvaud et al., 2021; Gasteiger, 2003).

Secondary structures were predicted using PSIPRED 4.0 (Buchan et al., 2024; Jones, 1999) and the NPS@ v2.16.0 (Combet et al., 2000). Three-dimensional structural predictions for LaNDHs were carried out with the AlphaFold2 v2.1.1 (Wayment-Steele et al., 2023; Jumper et al., 2021). Active site residues were identified using the GalaxyWEB program (Ko et al., 2012; Heo et al., 2013, 2016; Seok et al., 2021). Multiple sequence alignment was performed using the LSQKAB program within the CCP4 suite (Collaborative Computational Project N, 1994), and the root mean square deviation (RMSD) for Cα atoms was calculated. Structural visualizations were generated using PyMOL 2.3.4 (https://www.pymol.org/2/).

To validate the tertiary structures, Ramachandran plots for LaNDHs were generated using the PDBsum database (de Beer et al., 2014; Laskowski et al., 2017; Laskowski, 2022, 2004). This tool evaluates the quality of protein structures by detecting geometric errors, thereby enhancing the accuracy of the models. The Ramachandran plot specifically analyzes the stereochemical properties by displaying the dihedral angles of amino acid residues, identifying the allowed conformational regions, and highlighting any disallowed orientations.

On the other hand, ProSA (Protein Structure Analysis) is a commonly used tool for analyzing and validating predicted protein models (Wiederstein and Sippl, 2007). The z-score provides an overall assessment of model quality and is plotted against the z-scores of all experimentally determined protein structures in the current PDB. This plot distinguishes between structural types (e.g., X-ray, NMR) using color coding, enabling the evaluation of whether the z-score for the input structure falls within the expected range for native proteins of similar size.

To quantify the expression levels of the target gene under different light conditions, real-time quantitative polymerase chain reaction (RT-qPCR) was conducted using PowerUp SYBR Green Master Mix (Applied Biosystems). Plant tissue samples (roots, stems, leaves, and flowers) were collected, immediately flash-frozen in liquid nitrogen, and stored at -80°C for later analysis. The light treatments included white, red, and blue light, with red light having a maximum wavelength of 660 nm and blue light having a maximum wavelength of 450 nm. The light intensity was set at 100 µmol/(m·s). For cadmium ion (Cd²⁺) stress, concentrations of 0, 20, and 40 µM Cd²⁺ were applied, while temperature stress was tested at 25°C, 30°C, and 40°C. Salt stress was induced using 0, 200, and 300 mM NaCl, respectively. Total RNA was extracted using the Universal Plant Total RNA Extraction Kit (Bioteke, Beijing, China) according to the manufacturer’s protocol. cDNA synthesis was performed with the PrimeScript 1st Strand cDNA Synthesis Kit (Takara, Kyoto, Japan). Primer sequences are listed in Supplementary Table S2 . The PCR reaction volume was 20 μL, with the following conditions: 90°C for 5 min, followed by 40 cycles of 95°C for 10 s and 60°C for 30 s, and a final step of 95°C for 15 s and 60°C for 60 s. RT-qPCR was performed using an Applied Biosystems QuantStudio 5 instrument. Data were analyzed using the 2^-ΔΔCT method (Hawkins and Guest, 2017; Green and Sambrook, 2018), and relative expression levels were presented as log₂ values in histograms. Beta-actin gene is expressed at relatively constant levels in different tissues and cells and is used to detect changes in gene expression levels. Beta-actin was used as the reference gene, with expression normalized to untreated controls. A positive control was included for the beta-actin gene. A ratio greater than zero indicated upregulation, while a ratio less than zero indicated downregulation.

All experiments were conducted at least in triplicate. The data were expressed as mean ± SD. Statistical analysis was conducted using Origin 8.5, Microsoft Excel 2013 and SPSS 19.0. In the all statistical evaluations, p < 0.05 was considered statistically significant, and p < 0.01 was considered high statistically significant.