Next Article in Journal
Effects of Gene–Lifestyle Interaction on Obesity Among Students
Previous Article in Journal
Can Global DNA Methylation Be Influenced by Polymorphisms in Genes Involved in Epigenetic Mechanisms? A Review
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Screening Key Genes Related to Nitrogen Use Efficiency in Cucumber Through Weighted Gene Co-Expression Network Analysis

1
College of Life Science, Nankai University, Tianjin 300071, China
2
Cucumber Research Institute, Tianjin Academy of Agricultural Sciences, Tianjin 300192, China
3
State Key Laboratory of Vegetable Biobreeding, Tianjin Academy of Agricultural Sciences, Tianjin 300381, China
4
College of Agricultural Science, Nankai University, Tianjin 300071, China
*
Author to whom correspondence should be addressed.
These authors contributed equally to this work.
Genes 2024, 15(12), 1505; https://doi.org/10.3390/genes15121505
Submission received: 31 October 2024 / Revised: 21 November 2024 / Accepted: 22 November 2024 / Published: 24 November 2024
(This article belongs to the Special Issue Genetics and Breeding of Horticulture Crops—2nd Edition)

Abstract

:
Cucumber (Cucumis sativus L.) is a crucial vegetable crop, requiring significant nitrogen fertilizer inputs. However, excessive nitrogen application not only impairs growth but also poses severe environmental risks. Thus, enhancing nitrogen use efficiency (NUE) in cucumber is imperative. For the identification of genes associated with NUE in cucumber, roots of high NUE and low NUE lines were analyzed under high nitrogen conditions. Using transcriptome sequencing through WGCNA, a total of 15,180 genes were categorized into 35 co-expression modules, with 5 modules being highly correlated with NUE. Based on differential expression within the five modules and the results of GO and KEGG enrichment analyses, 25 genes were identified as potentially related to NUE. Among these, CsaV4_1G002492 (GLR22), CsaV4_2G003460 (GLR35), CsaV4_3G000307 (NRT1.1), and CsaV4_7G001709 (UPS2) were homologous to genes in Arabidopsis known to directly participate in NUE related process. These four genes were chosen as key genes for further analysis. qRT-PCR analysis revealed that CsaV4_3G000307 and CsaV4_7G001709 were more active during the early stages of the high nitrogen treatment in the high NUE line. Conversely, CsaV4_1G002492 and CsaV4_2G003460 were more active in the low NUE line. Using transcriptomic analysis, a frameshift INDEL mutation was observed in CsaV4_3G000307 in the low NUE line, which impacted the compactness of the protein structure, potentially altering its function. Analysis of protein interactions of these four key genes predicted some potential interaction networks. This research offers critical insights into the genetic factors influencing NUE in cucumber, presenting potential targets for genetic modification or breeding programs.

1. Introduction

Cucumber is one of the top ten vegetables globally, widely cultivated across various regions, and plays a significant role as an economic crop with substantial impacts on global agriculture [1]. According to the Food and Agriculture Organization (FAO) (http://agris.fao.org/), the total production area of cucumbers worldwide was 2.174 million hectares with a total yield of 94.718 million tons in 2022. Cucumber is a crop with high nutrient demands, particularly with a significant requirement for nitrogen fertilizers [2]. The application of nitrogen fertilizers significantly increases cucumber yield and is crucial for achieving high productivity [3,4]. However, excessive application of nitrogen fertilizer inhibits growth of cucumbers and increase economic costs [5,6,7]. With increasing nitrogen fertilizer application, the NUE of cucumbers tends to decline [8,9]. In Chinese agricultural production, the average utilization rate of nitrogen fertilizers is only 30–35%, meaning that approximately 65–70% of the nitrogen is not effectively used by the crops [10]. This results not only in substantial economic losses but also in severe environmental issues, such as water eutrophication, groundwater pollution, greenhouse gas emissions, and soil acidification [11,12,13,14]. Consequently, optimizing nitrogen fertilizer application and enhancing NUE are imperative for sustainable agricultural development.
However, the processes of nitrogen response and efficiency are governed by intricate networks involving metabolic, developmental, and environmental signals. Considerations of nitrogen absorption, transport, assimilation, and extensive signaling networks are necessary for enhancing NUE. Robust modern strategies are employed for identifying key genes. Chu et al. used genomic positional mapping and pinpointed a gene critical for high nitrate utilization, NRT1.1b. A SNP mutation in this gene has been demonstrated to influence nitrogen use efficiency in rice [15]. Zhou et al., through comparative transcriptomic and metabolomic analysis of maize and rice leaf tissues, pinpointed the gene OsDREB1C, which was involved in augmenting photosynthesis and NUE [16]. Gao et al. demonstrated that analysis of the expression profile of CsNRT2.1 showed that CsNRT2.1 is a high affinity nitrate transporter [17]. Ming et al. analyzed the expression patterns and subcellular localization of the GS1 protein under different nitrogen conditions, confirming its role in enhancing low nitrogen tolerance in cucumber [18]. Transcriptomic analysis through WGCNA is considered a powerful means of mining key genes [19]. WGCNA categorizes a multitude of genes into different modules based on gene expression changes, reducing analytical complexity and associating co-expression patterns within modules with phenotypic differences, thus identifying gene modules and key genes that play crucial roles during phenotypic variations [20,21,22,23]. This method has been extensively applied to unravel molecular mechanisms in various plants, including cucumbers [24]. Examples include identifying key genes responsible for resistance against Podosphaera xanthii in cucumber cultivars [25], pinpointing genes crucial for phenylpropanoid biosynthesis during cucumber storage [26], determining central genes for salt response [27], and identifying genes significantly associated with flavonoid biosynthesis in cucumbers [28]. Moreover, SNPs and INDELs are genetic variations that can significantly influence phenotypic traits [15,29,30]. The accurate identification of these mutations is crucial for genetic research and breeding programs. The presence of an SNP or INDEL may lead to phenotypic changes, playing a significant role in genetic breeding [31]. For instance, SNPs and INDELs have been used to assist in selecting genes resistant to powdery mildew [32], genes resistant to stem blight [33], genes associated with germination [34], and key genes for parthenocarpy in cucumbers [35].
Previously, we have established laboratory techniques for identifying phenotypes under varying nitrogen concentrations during the seedling stage, and yielded materials with extreme genotypes. Five candidate genes associated with cucumber tolerance to low nitrogen were identified through genome-wide association studies [36]. In this study, two extreme phenotypic cucumber lines under high nitrogen conditions were used. Transcriptomic sequencing of root tissues at seven time points after high nitrogen treatment was performed. Using the WGCNA method, a co-expression network of genes related to nitrogen use efficiency under high nitrogen conditions was constructed. Specific modules responding to NUE in cucumber were selected, and their functions were annotated and analyzed to mine related key genes and interaction networks. This analysis offers new insights into the regulation mechanisms of NUE in cucumber.

2. Materials and Methods

2.1. Plant Materials and High Nitrogen Treatment

Two extreme phenotypic cucumber lines, L2 with high NUE and L8 with low NUE, were used as plant material, which were obtained through previous screening of 112 cucumber lines of different cucumber types at Cucumber Research Institute of Tianjin Academy of Agricultural Sciences, China. L2 and L8 both belong to the South China cucumber type. Seeds were initially soaked in warm water at 55 °C for 15 min, followed by immersion in room temperature water for 3.5 h. Subsequently, seeds were placed in a germination culture box maintained at 28 °C for 2–3 days to allow sprout development. The sprouts were affixed to sponge plugs and then were transferred to a plastic box and cultivated in a half-strength Hoagland nutrient solution (nitrate concentration of 7 mmol/L) at a pH of 6.2. The nutrient solution was refreshed every three days and aeration was conducted twice a day. Once seedlings reached the one true leaf stage, they were transferred to 700 mL black hydroponic plastic buckets containing 500 mL of half-strength Hoagland nutrient solution. When seedlings reached the two true leaf stage, they were subjected to a higher nitrate treatment with a nitrogen concentration of 20 mmol/L. The Hoagland nutrient solution recipe is available in the Supplementary Table S1. Root sampling was conducted at 0 h, 2 h, 4 h, 8 h, 12 h, 24 h, and 30 h of cultivation, with at least three biological replicates per treatment.

2.2. RNA Extraction and Transcriptome Sequencing

Total RNA was extracted from root tissues using the standard extraction method provided by Novogene Co. Ltd. The quality and integrity of RNA were assessed with an Agilent 2100 bioanalyzer. Libraries for transcriptomic sequencing were prepared using the standard NEB library construction kit [37] and were sequenced on the Illumina platform to generate paired-end reads. Raw sequence reads were processed and aligned to the cucumber reference V4 genome (http://www.cucumberdb.com/#/home, accessed on 25 September 2024) utilizing HISAT2 [38].

2.3. Weighted Gene Co-Expression Network Analysis

Transcriptomic data from 14 samples were analyzed using the Novamagic cloud platform with a soft threshold set at 10. Genes were hierarchically clustered based on their Topological Overlap Matrix (TOM) dissimilarity [19], and modules were identified using dynamic tree cutting techniques. Each module contained a minimum of 30 genes. Modules that exhibited highly similar expression profiles (dissimilarity less than 0.25) were merged to minimize redundancy and enhance the robustness of the findings.

2.4. Differential Expression Genes Analysis

Differentially expressed genes (DEGs) between the high NUE and low NUE lines at each time point were identified using the edgeR [39]. The p-values were adjusted using the Benjamini and Hochberg method. Corrected p-value of 0.05 and absolute log2(foldchange) of 0.8 were set modules as the thresholds for significantly differential expression. A DEGs Venn diagram was created using an R 4.4.1 package.

2.5. Identification and Enrichment Analysis of the Core Modules

Identification and enrichment analysis of the core modules were conducted by calculating the Pearson correlation coefficients and the significance of the p-values for the module eigengenes in relation to L2 and L8. Modules were selected as core if they exhibited a correlation coefficient (r) greater than 0.6 and a p-value less than 0.05. Differentially expressed genes (DEGs) within these modules were then analyzed using Gene Ontology (GO, http://www.geneontology.org/, accessed on 11 October 2024) and the Kyoto Encyclopedia of Genes and Genomes (KEGG, http://www.genome.jp/kegg/, accessed on 11 October 2024) to elucidate the biological functions and pathways predominantly represented. Bar charts were utilized to display the top 30 significant GO enrichment analysis results. Bubble charts were employed to show the top 20 KEGG pathway analysis results.

2.6. Real-Time Quantitative PCR Analysis

Primers were designed using NCBI (https://www.ncbi.nlm.nih.gov/, accessed on 1 October 2024) and synthesized by Sangon Biotech Co., Ltd. (Shanghai, China). Detailed primer sequence information is available in Supplementary Table S2. RNA extraction was performed using the NG312 kit from Beijing LABLEAD Company, Beijing, China, and cDNA was synthesized with a First-strand cDNA Synthesis Mix kit from the same company. All procedures were carried out according to the instructions provided. The qRT-PCR was conducted using the 2x Realab Green PCR Fast mixture reagent from Beijing LABLEAD Company, Beijing, China in a 20 µL system as per the manufacturer’s instructions. The protocol included an initial denaturation at 95 °C for 30 s, followed by 40 cycles of denaturation at 95 °C for 10 s and annealing at 60 °C for 30 s. A melting curve analysis was performed from 60 to 90 °C in 10 s intervals with a ∆T of 0.5 °C. An actin gene served as the internal control for data standardization, with relative expression calculated using the 2−∆∆Ct method. Variance analysis was conducted using R with a t-test to determine significance.

2.7. Transcriptomic SNPs and INDELs Analysis

SNPs and INDELs were analyzed using the Novamagic cloud platform. Sequence logos were generated to visualize these variations using the WebLogo tool, providing a graphical representation of sequence conservation and variation at identified loci [40,41].

2.8. Protein Structure and Interaction Prediction

Protein structure prediction was performed using AlphaFold3, with the default settings for speed and copies parameters [42]. The detailed protein sequence is in Supplementary File S1. STRING was used to predict protein interaction, employing default settings (https://cn.string-db.org/, accessed on 26 October 2024).

3. Results

3.1. Identification of Gene Co-Expression Modules Through WGCNA

An optimal soft-threshold power was selected for the construction of the gene co-expression network by assessing the scale-free topology. Various soft-threshold powers, ranging from 1 to 20, were evaluated (Figure 1). A power of 10 was determined to be optimal, corresponding to a fit index that exceeded 0.8 and ensuring stabilization and the attainment of a scale-free topology within the network.
Initial clustering based on the TOM revealed an extensive array of color blocks, reflecting a complex pattern of gene co-expressions that complicated further analysis (Supplementary Figure S1). The clustering led to the identification of 35 gene co-expression modules, which varied significantly in size. The smallest module, labeled grey, contained 48 genes, while the largest, labeled turquoise, included over 4000 genes, demonstrating the diversity of gene co-expression across the sampled conditions.

3.2. Correlation Between the Modules and NUE

To elucidate the roles of the 35 identified modules in the two lines, associations were explored by calculating Pearson correlation coefficients between module eigengenes and the L2 (high NUE)/ L8 (low NUE) lines. As shown in Figure 2, a detailed heatmap illustrates the strength and significance of these correlations, highlighting each module’s potential contribution to NUE. Five out of the 35 modules identified were found to be significantly correlated with the extreme NUE lines. The purple (correlation coefficient 0.99), skyblue (correlation coefficient 0.95), cyan (correlation coefficient 0.72), darkorange (correlation coefficient 0.71), and lightcyan (correlation coefficient 0.68) modules, contain 241, 86, 200, 101, and 195 genes, respectively.
Modules such as darkorange and skyblue, which show positive correlations, were considered to be active in the high NUE line, whereas modules like purple, cyan, and lightcyan, showing negative correlations, were considered to be active in the low NUE line, all of which were targeted for further comprehensive functional studies.

3.3. Comprehensive Analysis of the Core Modules

3.3.1. Analysis of the Purple Module

The purple module, encompassing 241 genes, displayed distinct expression patterns across the seven different time points in response to high nitrogen. Figure 3 illustrates that eigengene expressions differ significantly between the high NUE and low NUE cucumber lines. A peak in expression at 4 h in the high NUE line corresponded to the lowest expression in the low NUE line, highlighting significant differences in high nitrogen response.
The maximum number of DEGs within this module was observed at 4 h, totaling 145 genes (Supplementary Figure S2). The time points at 2 and 8 h each recorded 140 DEGs, while the lowest count of 129 DEGs occurred at 24 h. Seventy-four genes were differentially expressed at all time points, underlining their possible vital role in nitrogen metabolism.
GO (Supplementary Figure S3 and Table S3) and KEGG enrichment analyses (Supplementary Figure S4 and Table S4), elucidated the biological functions and pathways of these DEGs. The enrichment analysis highlighted several crucial biological processes and molecular functions related to NUE. Key NUE-related functions and pathways included glycine, serine, and threonine metabolism, aminoacyl-tRNA biosynthesis, and arginine biosynthesis and metabolism. Fifteen genes were identified, with nine exhibiting differential expression at four or more time points. These nine genes were selected as candidate genes for further analysis.

3.3.2. Analysis of the Lightcyan Module

The lightcyan module, encompassing 195 genes, exhibited diverse expression patterns across the seven time points (Figure 4). Significant differences in expression between the high NUE and low NUE lines were observed, which was particularly notable at 4 h where the low NUE line displayed remarkable peak eigengene expressions in contrast to the lower expressions in the high NUE line.
Within this module, the highest number of DEGs was observed at 4 h, totaling 147 genes (Supplementary Figure S5). This peak was followed by 98 DEGs at 24 h, with the lowest count of 29 DEGs at 2 h. While DEG counts at time points 0, 8, 12, and 30 h were recorded at 35, 63, 82, and 55, respectively. Notably, 74 genes were differentially expressed at four or more time points, with 7 genes differentially expressed across all time points.
GO enrichment identified 11 genes with NUE-related functions (Supplementary Figure S6 and Table S5) and KEGG pathways enrichment identified 10 genes (Supplementary Figure S7 and Table S6). After excluding duplicates, the enrichment analysis highlighted 19 genes with several crucial biological processes and molecular functions related to NUE, including arginine and proline metabolism, as well as glycine, serine, and threonine metabolism. Among these, 10 genes exhibited differential expression at four or more time points, emphasizing their potential roles in nitrogen metabolism. These 10 genes were selected as candidate genes for further analysis.

3.3.3. Analysis of the Darkorange Module

The darkorange module encompassed 101 genes. Notably higher eigengene expressions were observed at 0 and 2 h in the high NUE line, marking these periods as critical for gene activity within this module. In contrast, the low NUE line exhibited negative values at these times (Figure 5).
The number of DEGs peaked at 0 h with 64 genes, followed by 58 at 2 h (Supplementary Figure S8). The lowest differential expression was noted at 4 h with 23 genes, with consistent counts of 30 genes observed at both 8 and 12 h, and varying counts of 29 and 39 genes at 24 and 30 h, respectively. A total of 35 genes were identified as differentially expressed at four or more time points, where 12 exhibited differential expression across all seven time points.
GO enrichment identified 10 genes with NUE-related functions (Supplementary Figure S9 and Table S7) and KEGG pathways enrichment identified 5 genes (Supplementary Figure S10 and Table S8). These analyses highlighted crucial pathways such as nitrogen metabolism and arginine biosynthesis, which were essential for understanding nitrogen utilization in plants. After excluding duplicates, 14 genes critical for nitrogen utilization were identified, where 5 genes showed differential expression across four or more time points, underscoring their significant roles in nitrogen metabolism. These five genes were selected as candidate genes from this module for further analysis.

3.3.4. The Analysis of the Cyan Module

The cyan module comprised 200 genes. Marked differences in eigengene expressions between the high NUE and low NUE lines were particularly prominent at 4 h, where the low NUE line showed a peak in positive values, which is in stark contrast to the negative values observed in the high NUE line (Figure 6).
The module exhibited a peak in differential expression at 4 h with 62 genes, and a secondary peak at 24 h with 47 genes (Supplementary Figure S11). The lowest differential expression was noted at 2 h with only 16 DEGs, and variable DEG counts were observed at other time points—29 at 0 h, 22 at 8 h, 35 at 12 h, and 36 at 30 h. Of these, 17 genes were identified as differentially expressed at four or more time points, with one gene expressed across all seven time points.
GO (Supplementary Figure S12 and Table S9) and KEGG pathways (Supplementary Figure S13 and Table S10) revealed the involvement of these genes in crucial metabolic pathways, including amino acid biosynthesis and arginine and proline metabolism. These pathways played vital roles in the plant metabolic framework, indirectly influencing nitrogen assimilation and utilization. However, there were no genes differentially expressed across four or more time points and therefore no genes were selected as candidate genes.

3.3.5. The Analysis of the Skyblue Module

The skyblue module comprised 86 genes. The expression heatmap revealed no distinct time-specific variations within this module (Figure 7). A peak in differential expression was observed at 12 h involving 31 genes, closely followed by 29 genes at 4 h (Supplementary Figure S14). The lowest number of differentially expressed genes, totaling 20, was recorded at 8 h. Counts of differentially expressed genes at 0 and 30 h included 28 genes each, with 25 and 26 genes at 2 and 24 h, respectively. Importantly, 17 genes in this module were identified as differentially expressed across four or more time points, while 6 genes differentially expressed at all time points.
GO enrichment identified two genes with NUE-related functions (Supplementary Figure S15 and Table S11) and KEGG pathways enrichment identified one gene (Supplementary Figure S16 and Table S12). Crucial functions and pathways included nitrogen compound transport and arginine and proline metabolism. Among these genes, only one gene showed differential expression across multiple time points. This one gene was selected as a candidate gene from this module for further analysis.

3.4. Information and Functional Annotation of 25 Selected Genes

Table 1 provides detailed information about the candidate genes from the five modules, including gene descriptions, species affiliations, and their respective expression modules. These genes have been studied across multiple species, such as in Arabidopsis thaliana, Oryza sativa, Cucumis melo, Cucumis sativus, and so on. We further screened 25 genes by reviewing related studies to determine whether the genes or their gene family members had been previously studied for functions related to NUE. The results identified genes directly related to NUE such as CsaV4_1G002492 (GLR22) and CsaV4_2G003460 (GLR35), which are glutamate receptors potentially involved in nitrogen signal transduction; CsaV4_3G000307 (NRT1.1/NPF6.3), which is involved in nitrogen absorption and transport; and CsaV4_7G001709 (UPS2), which is involved in nitrogen recycling. These four genes were selected as key genes for further analysis.

3.5. Expression Profile Analysis of the Key Genes

As shown in Figure 8a, after high nitrogen treatment, CsaV4_3G000307 was significantly upregulated at 2 and 4 h in the high NUE line, whereas no significant up-regulation was observed in the low NUE line, suggesting that CsaV4_3G000307 may play a more active role in the initial high nitrogen response in the high NUE line. A remarkable down regulation was observed in both lines at time point 8 h. The expression level recovered again at time points 12 h and 24 h in both the high NUE and low NUE lines. At the time point 30 h, the expression level declined remarkably in the low NUE line but maintained normal expression in the high NUE line.
As illustrated in Figure 8c, CsaV4_1G002492 showed no differences in expression before high nitrogen treatment. A significant down-regulation in the high NUE line was observed at all time points except at 8 h, and a significant up-regulation in the low NUE line was observed at all time points. Figure 8e revealed that CsaV4_2G003460 exhibited a higher expression in the low NUE line at 0, 2, 8, and 30 h than that in high NUE line. According to Figure 8g, CsaV4_7G001709 showed up-regulation at 2 h and 4 h, but down-regulation at other time points in the high NUE line, while in low NUE line, it was up-regulated at 4 h, but down-regulation at other time points.
This analysis demonstrates the distinct expression profiles of CsaV4_3G000307, CsaV4_7G001709, CsaV4_2G003460, and CsaV4_1G002492 in the high NUE and low NUE lines. CsaV4_3G000307 and CsaV4_7G001709 were suggested to be related to enhancing NUE, whereas CsaV4_2G003460 and CsaV4_1G002492 appeared to be reducing NUE. Linear regression analysis confirmed a high consistency between qPCR results and transcriptomic data (Figure 8b,d,f,h), verifying the reliability of the transcriptomic data.

3.6. Protein Interaction Prediction for CsaV4_3G000307

Four interacting proteins were identified through protein interaction predictions for CsaV4_3G000307 (Csa_3G027720) (Figure 9a). Notably, NIR and the protein of Csa_2G372190 were recognized within the nitrogen metabolism pathway, as indicated in the KEGG database (Figure 9b). The protein of Csa_3G730930 was linked to the ammonium transmembrane transporter complex. The protein of Csa_2G372190 was characterized as a non-specific serine/threonine protein kinase within the protein kinase superfamily. The protein of Csa_3G150160 was characterized as glutamine synthetase, which plays a crucial role in nitrogen assimilation.

3.7. Protein Interaction Prediction for CsaV4_1G002492

Five interacting proteins were identified through protein interaction predictions for CsaV4_1G002492 (Csa_1G418780) (Figure 10a). Csa_3G644800 contained the Cupin 1 domain. The proteins of Csa_2G296030 and Csa_7G447890 are characterized as having undefined functions. The protein of Csa_6G118360 possessed a C2H2-type zinc finger domain, which is typically involved in the transcriptional regulation of plants. The protein of Csa_7G453500 contained a plastocyanin-like domain, which is related to photosynthesis. Figure 10b highlighted the enriched functions of these proteins, including No apical meristem (NAM) protein, Lipoxygenase, Phytoene dehydrogenase activity, and the Cupin domain.

3.8. Protein Interaction Prediction of CsaV4_2G003460

Five interacting proteins were identified through protein interaction predictions for CsaV4_2G003460 (Csa_2G418930) (Figure 11a). Among them, the proteins of Csa_5G606590, Csa_1G039020, and Csa_1G012100 contained the t-SNARE coiled-coil homology domain, typically crucial in vesicular trafficking processes like endocytosis and exocytosis, essential for cellular transport mechanisms. The proteins of Csa_4G062380 and Csa_2G145880 were identified as serine hydroxymethyltransferases, directly involved in amino acid metabolism. As illustrated in Figure 11b, the enriched gene functions primarily included the glycine biosynthetic process from serine, L-serine catabolic process, folic acid metabolic process, tetrahydrofolate interconversion, and vesicle-mediated transport.

3.9. Protein Interaction Prediction of CsaV4_7G001709

Four interacting proteins were identified through protein interaction predictions for CsaV4_7G001709 (Csa_7G343330) (Figure 12a). Among these, the protein of the gene Csa_3G828970 has undefined characteristics. The proteins of Csa_3G894470, Csa_6G381850, and Csa_4G573860 all contained the Aa_trans (amino acid transporter) domain, which typically plays a crucial role in amino acid transport processes vital for nitrogen assimilation and metabolism. As shown in Figure 12b, the functions enriched for these proteins included amino acid transmembrane transport, nitrogen compound transport, and organic substance transport, which are all closely related to nitrogen metabolism.

3.10. Transcriptomic SNP and INDEL Analysis of the Key Genes

The structural integrity of the four genes was analyzed at the transcriptional level for the presence of SNPs and INDELs. No SNPs or INDELs were found in CsaV4_2G003460, CsaV4_1G002492 and CsaV4_7G001709 in either the high NUE or low NUE lines. However, an INDEL was identified in CsaV4_3G000307, causing a frame-shift mutation that occurs at the 522nd leucine exclusively in the low NUE line (Figure 13a and Figure 14a). Structural analysis for both the high NUE and low NUE lines revealed that the frame-shift mutation impacted the compactness of the protein, potentially altering its functionality (Figure 13b and Figure 14b).

4. Discussion

Nitrogen is essential for synthesizing chlorophyll, amino acids, and nucleic acids, making it indispensable for plant growth and development [61,62,63]. However, a significant amount of nitrate nitrogen is lost from the soil due to plant uptake and leaching by irrigation or rain, reducing nitrogen use efficiency [64]. By improving plant NUE, we can significantly reduce cultivation costs and mitigate nitrate pollution.
In this study, root transcriptomic responses to high nitrogen were evaluated across multiple time points in two contrasting cucumber genotypes, high NUE and low NUE, selected based on extensive field and laboratory research. A multitude of DEGs were identified at seven time points within both lines, affirming the validity of the material selection and laying a solid foundation for subsequent analyses. To isolate genes associated with NUE, WGCNA was employed to categorize 15,180 genes into 35 modules, ranging from the smallest gray module containing 48 genes to the largest turquoise module encompassing over 4000 genes. This modular division reflected the complexity of gene expression in the high NUE and low NUE lines. Notably, five modules showed expression patterns significantly exclusive to either the high NUE or low NUE lines. The purple module exhibited the highest correlation coefficient at 0.99, followed by skyblue at 0.95, darkorange at 0.71, cyan at 0.72, and lightcyan at 0.68. Among these, the purple module was particularly notable with 241 genes, 83.4% of which were differentially expressed between the two lines, whereas the cyan module, comprising 200 genes, showed only 53.5% differential expression, reducing its significance. Despite the skyblue module being highly correlated with the high NUE line, it contained the fewest genes, 86, with 63.95% showing differential expression, which suggests that a few key genes play a critical role in enhancing this module’s relevance to NUE. Indeed, gene CsaV4_7G001709 from the skyblue module was annotated as UPS2, which is a protein rich in nitrogen and directly related to NUE [57]. Additionally, genes CsaV4_2G003460 and CsaV4_1G002492 from the lightblue module are annotated as glutamate receptor 3.5 and 2.2 and are directly involved in nitrogen metabolism [43,44]. CsaV4_3G000307 and CsaV4_7G001709 may play positive regulatory roles in the NUE process. In contrast, CsaV4_2G003460 and CsaV4_1G002492 typically showed higher expression levels in the low NUE line, indicating their potential negative regulatory effects on NUE. Additional analysis explored whether these expression differences stemmed from structural changes within the genes, such as SNPs or INDELs, since these genetic variations are acknowledged as potential factors in phenotypic segregation [47].
Transcriptome analysis revealed a frame-shifting INDEL mutation exclusively in CsaV4_3G000307 within the low NUE line; no missense mutation SNPs were detected in any of other the three genes across either line. Both missense mutation SNPs and frame-shifting INDELs can modify amino acid sequences, potentially affecting protein functionality [48], thus structural predictions were conducted for CsaV4_3G000307’s sequences before and after the mutation. Post-mutation, the protein structure appeared looser, which might influence protein binding and functionality. However, whether this INDEL actually causes functional changes requires further verification, such as whether the gene also exhibits the same INDEL in the genome, whether the presence of this INDEL in the transcriptome genuinely affects protein coherence, and whether the occurrence of this INDEL is due to alternative splicing or sequencing errors; these aspects all necessitate subsequent verification.
The homolog gene of CsaV4_3G000307 in Arabidopsis thaliana, AtNPF6.3/NRT1.1, plays a pivotal role in detecting external nitrogen levels and transporting extracellular nitrogen into cells [49]. Positioned on the plasma membrane, AtNPF6.3 enables nitrate transport from the extracellular environment into the cytoplasm by facilitating phosphorylation, which activates nitrate ion channels [50]. This function is conserved across species, as demonstrated by the nitrate transport activity of AtNPF6.3 homologs in rice and maize [15,51]. Consequently, CsNPF6.3/NRT1.1 is likely to be integral to the nitrate transport process in cucumber as well. Furthermore, AtNPF6.3 is implicated in various other functions, including root development regulation via auxin transport, nitrate signal transduction, nutrient signal integration, and plant–microbe interaction regulation, highlighting its essential role in plant stress resistance and developmental processes [52,53,54]. NRT1.1 has been implicated in several regulatory networks; for example, STOP1 activation enhances NRT1.1-mediated nitrate absorption [55], while NRT1.1 collaborates with NLP7 in nitrogen sensing and transport [56]. It also partners with NRT2.1 to facilitate nitrate transport [65], and is negatively regulated by abscisic acid signaling, which phosphorylates NRT1.1 via SnRK2s, inhibiting nitrate uptake [66]. Interaction predictions for CsaV4_3G000307 identified four potential interacting genes or proteins: Csa_2G372190, Csa_3G730930, and Csa_3G150160 are identified as NIRs, where Csa_3G150160 serves as glutamine synthetase (GS). Once nitrate is transported to relevant sites, enzymes such as nitrate reductase (NR), NIR, GS, and glutamate synthase (GOGAT) transform it into amino acids that fulfill diverse biological functions [67]. Csa_3G730930 acts as an AMT, which is essential for the nitrogen needs of plants, particularly in flood-prone or acidic soils where ammonium is prevalent [68]. The protein of Csa_2G372190, which is identified as a non-specific serine/threonine protein kinase, is part of the protein kinase superfamily [69,70,71]. The protein kinase translated by Csa_2G372190 may stimulate the expression of certain pathways. If there are also differences in the expression of this gene between high NUE and low NUE lines, it could be related to the differential expression of CsaV4_3G000307 in the two lines. Additionally, Csa_3G730930 and CsaV4_3G000307 may have a synergistic effect on nitrogen absorption. The variable expression of CsaV4_3G000307 in the high NUE and low NUE lines might lead to changes in nitrogen absorption and transport efficiency, thereby impacting the role of proteins such as GS and NIR in the nitrogen assimilation process.
The proteins of CsaV4_1G002492 and CsaV4_2G003460 belong to GLRs. In Arabidopsis, overexpression of AtGLR2 leads to symptoms of Ca2+ deficiency [45], and Ca2+ is closely related to the absorption and transport of nitrogen within cells [56]. AtGLR1.1 is a regulator of carbon and nitrogen metabolism in Arabidopsis, controlling seed germination by influencing abscisic acid (ABA) to regulate carbon and nitrogen metabolism [46]. The enrichment of CsaV4_1G002492’s predicted interaction proteins, such as NAM protein, lipoxygenase, and phytoene dehydrogenase are directly or indirectly related to nitrogen metabolism [72,73,74]. Among them, the protein of Csa_6G118360 possessed a C2H2-type zinc finger domain, which typically involved in the transcriptional regulation of plants and the regulation of Csa_6G118360 may affect the differential expression of the CsaV4_1G002492 gene between the two extreme lines. The enriched interaction gene functions of CsaV4_2G003460 include the glycine biosynthetic process from serine, L-serine catabolic process, folic acid metabolic process, tetrahydrofolate interconversion, and vesicle-mediated transport, which are also directly or indirectly related to nitrogen metabolism [75,76,77]. The proteins of CsaV4_1G002492 and CsaV4_2G003460 may have a synergistic effect in the process of nitrogen metabolism.
The protein of CsaV4_7G001709 is UPS2. In soybean, the UPS1 protein plays a crucial role in managing the allocation of ureides from the source to the sink tissues. Enhancing the expression of UPS1 significantly increases the nitrogen supply to the sink tissues of soybean, which in turn contributes to higher grain yields [58,59,75]. In wheat, TaUPS1 and TaUPS2.1 are involved in nitrogen distribution and NUE. Three of the predicted interaction proteins (the proteins of Csa_3G894470, Csa_6G381850 and Csa_4G573860) contained the Aa_trans (amino acid transporter) domain and the functions enriched for these proteins include amino acid transmembrane transport, nitrogen compound transport, and organic substance transport, participating in nitrogen metabolism [78,79,80]. Overall, the enrichment of these functions indicates that the associated genes play a direct role in regulating plant nitrogen nutrition and metabolism, enhancing the four key genes potential significance in NUE.
Future research will employ transgenic technologies such as gene silencing, knockout, and over-expression to elucidate the roles of the key genes that regulate NUE in cucumbers. The investigation will extend to mining upstream and downstream genes and performing network analysis, enhancing the understanding of the genetic networks of NUE in cucumber.

5. Conclusions

The study analyzed transcriptomic data through WGCNA combined with other methods, identifying four key genes that are potentially instrumental in enhancing NUE. This investigation offers potential loci for molecular breeding and lays the groundwork for ensuing genetic and molecular inquiries dedicated to enhancing NUE in cucumbers.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/genes15121505/s1, Table S1: The Hoagland nutrient solution; Table S2: The primer sequence information; Table S3: GO Result of Purple Module; Table S4: KEGG Result of Purple Module; Table S5: GO Result of Lightcyan Module; Table S6: KEGG Result of Lightcyan Module; Table S7: GO Result of Darkorange Module; Table S8: KEGG Result of Darkorange Module; Table S9: GO Result of Cyan Module; Table S10: KEGG Result of Cyan Module; Table S11: GO Result of Skyblue Module; Table S12: KEGG Result of Skyblue Module; Figure S1: Gene co-expression modules identified through WGCNA; Figure S2: Venn Diagram of DEGs in Purple Module Across Time Points; Figure S3: GO enrichment analysis of DEGs in Purple Module; Figure S4: KEGG Pathway enrichment analysis of DEGs in Purple Module; Figure S5: Venn Diagram of DEGs in Lightcyan Module; Figure S6: GO enrichment analysis of DEGs in Lightcyan Module; Figure S7: KEGG enrichment analysis of DEGs in Lightcyan Module; Figure S8: Venn Diagram of DEGs in Darkorange Module; Figure S9: GO enrichment analysis of DEGs in Darkorange Module; Figure S10: KEGG enrichment analysis of DEGs in Darkorange Module; Figure S11: Venn Diagram of DEGs in Cyan Module; Figure S12: GO enrichment analysis of DEGs in Cyan Module; Figure S13: KEGG enrichment analysis of DEGs in Cyan Module; Figure S14: Venn Diagram of DEGs in Skyblue Module; Figure S15: GO enrichment analysis of DEGs in Skyblue Module; Figure S16: KEGG enrichment analysis of DEGs in Skyblue Module; File S1: The protein sequence of the key genes.

Author Contributions

Conceptualization, S.D.; methodology, A.W. and N.W.; software, C.L.; formal analysis, L.M.; resources, Y.H., N.L. and Z.C.; writing—original draft, L.M. and A.W.; supervision, S.D.; project administration, S.D.; funding acquisition, A.W. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the National Natural Science Foundation of China (Grant No. 32372708) and Tianjin Seed Industry Innovation Major Special Project (23ZXZYSN00010).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The raw Illumina sequence reads have been deposited into the National Center for Biotechnology Information under sequence read archive (SRA) PRJNA1175851.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Mallick, P.K. Evaluating Potential Importance of Cucumber (Cucumis sativus L.-Cucurbitaceae): A Brief Review. Int. J. Appl. Sci. Biotechnol. 2022, 10, 12–15. [Google Scholar] [CrossRef]
  2. Bai, L.; Deng, H.; Zhang, X.; Yu, X.; Li, Y. Gibberellin Is Involved in Inhibition of Cucumber Growth and Nitrogen Uptake at Suboptimal Root-Zone Temperatures. PLoS ONE 2016, 11, e0156188. [Google Scholar] [CrossRef] [PubMed]
  3. Hu, X.; Zhang, J.; Liu, W.; Wang, Q.; Wang, T.; Li, X.; Lu, X.; Gao, L.; Zhang, W. CsNPF7.2 Has a Potential to Regulate Cucumber Seedling Growth in Early Nitrogen Deficiency Stress. Plant Mol. Biol. Rep. 2020, 38, 461–477. [Google Scholar] [CrossRef]
  4. McAllister, C.H.; Beatty, P.H.; Good, A.G. Engineering nitrogen use efficient crop plants: The current status. Plant Biotechnol. J. 2012, 10, 1011–1025. [Google Scholar] [CrossRef] [PubMed]
  5. Zhang, W.; Wu, L.; Ding, Y.; Yao, X.; Wu, X.; Weng, F.; Li, G.; Liu, Z.; Tang, S.; Ding, C.; et al. Nitrogen fertilizer application affects lodging resistance by altering secondary cell wall synthesis in japonica rice (Oryza sativa). J. Plant Res. 2017, 130, 859–871. [Google Scholar] [CrossRef]
  6. Wu, L.; Zhang, W.; Ding, Y.; Zhang, J.; Cambula, E.D.; Weng, F.; Liu, Z.; Ding, C.; Tang, S.; Chen, L.; et al. Shading contributes to the reduction of stem mechanical strength by decreasing cell wall synthesis in japonica rice (Oryza sativa L.). Front. Plant Sci. 2017, 8, 881. [Google Scholar] [CrossRef]
  7. Haque, M.A.; Haque, M.M. M. Growth, yield and nitrogen use efficiency of new rice variety under variable nitrogen rates. Am. J. Plant Sci. 2016, 7, 612–622. [Google Scholar] [CrossRef]
  8. Hua, B.; Liang, F.; Zhang, W.; Qiao, D.; Wang, P.; Teng, H.; Zhang, Z.; Liu, J.; Miao, M. The Potential Role of bZIP55/65 in Nitrogen Uptake and Utilization in Cucumber Is Revealed via bZIP Gene Family Characterization. Plants 2023, 12, 3228. [Google Scholar] [CrossRef]
  9. Vitousek, P.M.; Naylor, R.; Crews, T.; David, M.B.; Drinkwater, L.; Holland, E.; Johnes, P.; Katzenberger, J.; Martinelli, L.; Matson, P.; et al. Nutrient imbalances in agricultural development. Science 2009, 324, 1519–1520. [Google Scholar] [CrossRef]
  10. Li, H.; Dai, M.W.; Dai, S.L.; Dong, X.J. Current Status and Environment Impact of Direct Straw Return in China’s Cropland—A Review. Ecotoxicol. Environ. Saf. 2018, 159, 293–300. [Google Scholar] [CrossRef]
  11. Ma, X.; Tan, Z.; Cheng, Y.; Wang, T.; Cao, M.; Xuan, Z.; Du, H. Water-Nutrient Coupling Strategies That Improve the Carbon, Nitrogen Metabolism, and Yield of Cucumber under Sandy Cultivated Land. Land 2024, 13, 958. [Google Scholar] [CrossRef]
  12. Cui, Z.; Chen, C.; Chen, Q.; Huang, J. Difference in the Contribution of Driving Factors to Nitrogen Loss with Surface Runoff between the Hill and Plain Agricultural Watersheds. J. Geophys. Res. Biogeosci. 2024, 129, e2023JG007931. [Google Scholar] [CrossRef]
  13. Zhang, S.; Hou, X.; Wu, C.; Zhang, C. Impacts of Climate and Planting Structure Changes on Watershed Runoff and Nitrogen and Phosphorus Loss. Sci. Total Environ. 2020, 706, 134489. [Google Scholar] [CrossRef] [PubMed]
  14. Ding, N.; Tao, F.; Chen, Y. Effects of Climate Change, Crop Planting Structure, and Agricultural Management on Runoff, Sediment, Nitrogen, and Phosphorus Losses in the Hai-River Basin since the 1980s. J. Clean. Prod. 2022, 359, 132066. [Google Scholar] [CrossRef]
  15. Hu, B.; Wang, W.; Ou, S.; Tang, J.; Li, H.; Che, R.; Zhang, Z.; Chai, X.; Wang, H.; Wang, Y.; et al. Variation in NRT1.1B contributes to nitrate-use divergence between rice subspecies. Nat. Genet. 2015, 47, 834–838. [Google Scholar] [CrossRef]
  16. Wei, S.; Li, X.; Lu, Z.; Zhang, H.; Ye, X.; Zhou, Y.; Li, J.; Yan, Y.; Pei, H.; Duan, F.; et al. A transcriptional regulator that boosts grain yields and shortens the growth duration of rice. Science 2022, 377, eabi8455. [Google Scholar] [CrossRef]
  17. Li, Y.; Li, J.; Yan, Y.; Liu, W.; Zhang, W.; Gao, L.; Tian, Y. Knock-Down of CsNRT2.1, a Cucumber Nitrate Transporter, Reduces Nitrate Uptake, Root length, and Lateral Root Number at Low External Nitrate Concentration. Front. Plant Sci. 2018, 9, 722. [Google Scholar] [CrossRef]
  18. Xin, M.; Qin, Z.-w.; Yang, J.; Zhou, X.-y.; Wang, L. Functional analysis of the nitrogen metabolism-related gene CsGS1 in cucumber. J. Integr. Agric. 2021, 20, 1515–1524. [Google Scholar] [CrossRef]
  19. Zhang, B.; Horvath, S. A general framework for weighted gene co-expression network analysis. Stat. Appl. Genet. Mol. Biol. 2005, 4, 17. [Google Scholar] [CrossRef]
  20. Horvath, S. Weighted Network Analysis: Applications in Genomics and Systems Biology; Springer: New York, NY, USA, 2011. [Google Scholar]
  21. Ivliev, A.E.; ’t Hoen, P.A.C.; Sergeeva, M.G. Coexpression network analysis identifies transcriptional modules related to proastrocytic differentiation and sprouty signaling in glioma. Cancer Res. 2010, 70, 10060–10070. [Google Scholar] [CrossRef]
  22. Hu, Y.; Wu, G.; Rusch, M.; Lukes, L.; Buetow, K.H.; Zhang, J.; Hunter, K.W. Integrated cross-species transcriptional network analysis of metastatic susceptibility. Proc. Natl. Acad. Sci. USA 2012, 109, 3184–3189. [Google Scholar] [CrossRef] [PubMed]
  23. Giotti, B.; Joshi, A.; Freeman, T.C. Meta-analysis reveals conserved cell cycle transcriptional network across multiple human cell types. BMC Genom. 2017, 18, 30. [Google Scholar] [CrossRef] [PubMed]
  24. Zinati, Z.; Nazari, L. Deciphering the Molecular Basis of Abiotic Stress Response in Cucumber (Cucumis Sativus L.) Using RNA-Seq Meta-Analysis, Systems Biology, and Machine Learning Approaches. Sci. Rep. 2023, 13, 12942. [Google Scholar] [CrossRef] [PubMed]
  25. Meng, X.; Yu, Y.; Song, T.; Yu, N.; Cui, N.; Ma, Z.; Chen, L.; Fan, H. Transcriptome Sequence Analysis of the Defense Responses of Resistant and Susceptible Cucumber Strains to Podosphaera Xanthii. Front. Plant Sci. 2022, 13, 872218. [Google Scholar] [CrossRef] [PubMed]
  26. Wang, J.; Tian, P.; Sun, J.; Li, B.; Jia, J.; Yuan, J.; Li, X.; Gu, S.; Pang, X. CsMYC2 Is Involved in the Regulation of Phenylpropanoid Biosynthesis Induced by Trypsin in Cucumber (Cucumis Sativus) during Storage. Plant Physiol. Biochem. 2023, 196, 65–74. [Google Scholar] [CrossRef]
  27. Karimi, M.; Pakdel, M.H.; Lashaki, K.B.; Soorni, A. Identification of hub salt-responsive genes in Cucumis sativus using a long non-coding RNA and mRNA interaction network. Hort. Environ. Biotechnol. 2022, 63, 539–556. [Google Scholar] [CrossRef]
  28. Wang, J.; Jia, J.; Sun, J.; Pang, X.; Li, B.; Yuan, J.; Chen, E.; Li, X. Trypsin Preservation: CsUGT91C1 Regulates Trilobatin Biosynthesis in Cucumis Sativus during Storage. Plant Growth Regul. 2023, 100, 633–646. [Google Scholar] [CrossRef]
  29. López-Kleine, L.; Leal, L.; López, C. Biostatistical approaches for the reconstruction of gene co-expression networks based on transcriptomic data. Brief. Funct. Genom. 2013, 12, 457–467. [Google Scholar] [CrossRef]
  30. Zhang, T.; Gu, M.; Liu, Y.; Lv, Y.; Zhou, L.; Lu, H.; Liang, S.; Bao, H.; Zhao, H. Development of novel InDel markers and genetic diversity in Chenopodium quinoa through whole-genome re-sequencing. BMC Genom. 2017, 18, 685. [Google Scholar] [CrossRef]
  31. Roy, S.C.; Lachagari, V.B.R. Assessment of SNP and InDel Variations Among Rice Lines of Tulaipanji × Ranjit. Rice Sci. 2017, 24, 336–348. [Google Scholar] [CrossRef]
  32. Xu, Q.; Shi, Y.; Yu, T.; Xu, X.; Yan, Y.; Qi, X.; Chen, X. Whole-Genome Resequencing of a Cucumber Chromosome Segment Substitution Line and Its Recurrent Parent to Identify Candidate Genes Governing Powdery Mildew Resistance. PLoS ONE 2016, 11, e0164469. [Google Scholar] [CrossRef] [PubMed]
  33. Han, J.; Dong, S.; Shi, Y.; Miao, H.; Liu, X.; Beckles, D.M.; Gu, X.; Zhang, S. Fine mapping and candidate gene analysis of gummy stem blight resistance in cucumber stem. Theor. Appl. Genet. 2022, 135, 3117–3125. [Google Scholar] [CrossRef] [PubMed]
  34. Cao, M.; Li, S.; Deng, Q.; Wang, H.; Yang, R. Identification of a major-effect QTL associated with pre-harvest sprouting in cucumber (Cucumis sativus L.) using the QTL-seq method. BMC Genom. 2021, 22, 249. [Google Scholar] [CrossRef] [PubMed]
  35. Wu, Z.; Zhang, T.; Li, L.; Xu, J.; Qin, X.; Zhang, T.; Cui, L.; Lou, Q.; Li, J.; Chen, J. Identification of a stable major-effect QTL (Parth 2.1) controlling parthenocarpy in cucumber and associated candidate gene analysis via whole genome re-sequencing. BMC Plant Biol. 2016, 16, 182. [Google Scholar] [CrossRef] [PubMed]
  36. Li, B.; Wei, A.; Tong, X.; Han, Y.; Liu, N.; Chen, Z.; Yang, H.; Wu, H.; Lv, M.; Wang, N.N.; et al. A Genome-Wide Association Study to Identify Novel Candidate Genes Related to Low-Nitrogen Tolerance in Cucumber (Cucumis sativus L.). Genes 2023, 14, 662. [Google Scholar] [CrossRef]
  37. Parkhomchuk, D.; Borodina, T.; Amstislavskiy, V.; Banaru, M.; Hallen, L.; Krobitsch, S.; Lehrach, H.; Soldatov, A. Transcriptome analysis by strand-specific sequencing of complementary DNA. Nucleic Acids Res. 2009, 37, e123. [Google Scholar] [CrossRef]
  38. Mortazavi, A.; Williams, B.A.; McCue, K.; Schaeffer, L.; Wold, B. Mapping and quantifying mammalian transcriptomes by RNA-Seq. Nat. Methods 2008, 5, 621–628. [Google Scholar] [CrossRef]
  39. Robinson, M.D.; McCarthy, D.J.; Smyth, G.K. edgeR: A Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 2010, 26, 139–140. [Google Scholar] [CrossRef]
  40. Crooks, G.E.; Hon, G.; Chandonia, J.M.; Brenner, S.E. WebLogo: A sequence logo generator. Genome Res. 2004, 14, 1188–1190. [Google Scholar] [CrossRef]
  41. Schneider, T.D.; Stephens, R.M. Sequence Logos: A New Way to Display Consensus Sequences. Nucleic Acids Res. 1990, 18, 6097–6100. [Google Scholar] [CrossRef]
  42. Jumper, J.; Evans, R.; Pritzel, A.; Green, T.; Figurnov, M.; Ronneberger, O.; Tunyasuvunakool, K.; Bates, R.; Žídek, A.; Potapenko, A.; et al. Highly accurate protein structure prediction with AlphaFold. Nature 2021, 596, 583–589. [Google Scholar] [CrossRef] [PubMed]
  43. Liu, K.; Sutter, B.M.; Tu, B.P. Autophagy sustains glutamate and aspartate synthesis in Saccharomyces cerevisiae during nitrogen starvation. Nat. Commun. 2021, 12, 57. [Google Scholar] [CrossRef] [PubMed]
  44. Guo, X.; Nan, Y.; He, H.; Ma, B.-L.; McLaughlin, N.B.; Wu, X.; Chen, B.; Gao, Y. Post-flowering nitrogen uptake leads to the genotypic variation in seed nitrogen accumulation of oilseed rape. Plant Soil 2021, 461, 281–294. [Google Scholar] [CrossRef]
  45. Kim, S.A.; Kwak, J.M.; Jae, S.K.; Wang, M.H.; Nam, H.G. Overexpression of the AtGluR2 Gene Encoding anArabidopsis Homolog of Mammalian Glutamate Receptors Impairs CalciumUtilization and Sensitivity to Ionic Stress in TransgenicPlants. Plant Cell Physiol. 2001, 42, 74–84. [Google Scholar] [CrossRef] [PubMed]
  46. Kang, J.; Turano, F.J. The putative glutamate receptor 1.1 (AtGLR1.1) functions as a regulator of carbon and nitrogen metabolism in Arabidopsis thaliana. Proc. Natl. Acad. Sci. USA 2003, 100, 6872–6877. [Google Scholar] [CrossRef]
  47. Seo, J.-H.; Dhungana, S.K.; Kang, B.-K.; Baek, I.-Y.; Sung, J.-S.; Ko, J.-Y.; Jung, C.-S.; Kim, K.-S.; Jun, T.-H. Development and Validation of SNP and InDel Markers for Pod-Shattering Tolerance in Soybean. Int. J. Mol. Sci. 2022, 23, 2382. [Google Scholar] [CrossRef]
  48. Barton, H.J.; Zeng, K. New Methods for Inferring the Distribution of Fitness Effects for INDELs and SNPs. Mol. Biol. Evol. 2018, 35, 1536–1546. [Google Scholar] [CrossRef]
  49. Shaofen, J.; Jingsong, L.; Qiong, L.; Qiang, L.; Chunyun, G.; Zhenhua, Z. NRT1.1 Regulates Nitrate Allocation and Cadmium Tolerance in Arabidopsis. Front. Plant Sci. 2019, 10, 384. [Google Scholar] [CrossRef]
  50. Sun, J.; Zheng, N. Molecular Mechanism Underlying the Plant NRT1.1 Dual-Affinity Nitrate Transporter. Front. Physiol. 2015, 6, 386. [Google Scholar] [CrossRef]
  51. Wen, Z.; Tyerman, S.D.; Dechorgnat, J.; Ovchinnikova, E.; Dhugga, K.S.; Kaiser, B.N. Maize NPF6 proteins are homologs of Arabidopsis CHL1 that are selective for both nitrate and chloride. Plant Cell 2017, 29, 2581–2596. [Google Scholar] [CrossRef]
  52. Wang, W.; Hu, B.; Li, A.; Chu, C. NRT1.1s in plants: Functions beyond nitrate transport. J. Exp. Bot. 2020, 71, 4373–4379. [Google Scholar] [CrossRef] [PubMed]
  53. Fang, X.Z.; Fang, S.Q.; Ye, Z.Q.; Liu, D.; Zhao, K.L.; Jin, C.W. Dual-Affinity Nitrate Transport/Signalling and Its Roles in Plant Abiotic Stress Resistance. Front. Plant Sci. 2021, 12, 715694. [Google Scholar] [CrossRef] [PubMed]
  54. Rolly, N.K.; Yun, B.-W. Regulation of Nitrate (NO3) Transporters and Glutamate Synthase-Encoding Genes under Drought Stress in Arabidopsis: The Regulatory Role of AtbZIP62 Transcription Factor. Plants 2021, 10, 2149. [Google Scholar] [CrossRef] [PubMed]
  55. Ye, J.Y.; Tian, W.H.; Zhou, M.; Zhu, Q.Y.; Du, W.X.; Zhu, Y.X.; Liu, X.X.; Lin, X.Y.; Zheng, S.J.; Jin, C.W. STOP1 Activates NRT1.1-Mediated Nitrate Uptake to Create a Favorable Rhizospheric pH for Plant Adaptation to Acidity. Plant Cell 2021, 33, 3658–3674. [Google Scholar] [CrossRef] [PubMed]
  56. Jia, Y.; Qin, D.; Zheng, Y.; Wang, Y. Finding Balance in Adversity: Nitrate Signaling as the Key to Plant Growth, Resilience, and Stress Response. Int. J. Mol. Sci. 2023, 24, 14406. [Google Scholar] [CrossRef]
  57. Takagi, H.; Watanabe, S.; Tanaka, S.; Matsuura, T.; Mori, I.C.; Hirayama, T.; Shimada, H.; Sakamoto, A. Disruption of ureide degradation affects plant growth and development during and after transition from vegetative to reproductive stages. BMC Plant Biol. 2018, 18, 287. [Google Scholar] [CrossRef]
  58. Lu, M.Z.; Carter, A.M.; Tegeder, M. Altering ureide transport in nodulated soybean results in whole-plant adjustments of metabolism, assimilate partitioning, and sink strength. J. Plant Physiol. 2022, 269, 153613. [Google Scholar] [CrossRef]
  59. Thu, S.W.; Lu, M.-Z.; Carter, A.M.; Collier, R.; Gandin, A.; Sitton, C.C.; Tegeder, M. Role of ureides in source-to-sink transport of photoassimilates in non-fixing soybean. J. Exp. Bot. 2020, 71, 4495–4511. [Google Scholar] [CrossRef]
  60. Meng, X.; Zhang, Z.; Wang, H.; Nai, F.; Wei, Y.; Li, Y.; Wang, X.; Ma, X.; Tegeder, M. Multi-scale analysis provides insights into the roles of ureide permeases in wheat nitrogen use efficiency. J. Exp. Bot. 2023, 74, 5564–5590. [Google Scholar] [CrossRef]
  61. Wen, B.; Xiao, W.; Mu, Q.; Li, D.; Chen, X.; Wu, H.; Li, L.; Peng, F. How does nitrate regulate plant senescence? Plant Physiol. Biochem. 2020, 157, 60–69. [Google Scholar] [CrossRef]
  62. Luo, L.; Zhang, Y.; Xu, G. How does nitrogen shape plant architecture? J. Exp. Bot. 2020, 71, 4415–4427. [Google Scholar] [CrossRef] [PubMed]
  63. Mur, L.A.J.; Simpson, C.; Kumari, A.; Gupta, A.K.; Gupta, K.J. Moving nitrogen to the centre of plant defence against pathogens. Ann. Bot. 2017, 119, 703–709. [Google Scholar] [CrossRef] [PubMed]
  64. Wang, X.; Zou, C.; Gao, X.; Guan, X.; Zhang, Y.; Shi, X.; Chen, X. Nitrate leaching from open-field and greenhouse vegetable systems in China: A meta-analysis. Environ. Sci. Pollut. Res. 2018, 25, 31007–31016. [Google Scholar] [CrossRef] [PubMed]
  65. Lay-Pruitt, K.S.; Takahashi, H. Integrating N Signals and Root Growth: The Role of Nitrate Transceptor NRT1.1 in Auxin-Mediated Lateral Root Development. J. Exp. Bot. 2020, 71, 4365–4368. [Google Scholar] [CrossRef] [PubMed]
  66. Su, H.; Wang, T.; Ju, C.; Deng, J.; Zhang, T.; Li, M.; Tian, H.; Wang, C. Abscisic Acid Signaling Negatively Regulates Nitrate Uptake via Phosphorylation of NRT1.1 by SnRK2s in Arabidopsis. J. Integr. Plant Biol. 2021, 63, 597–610. [Google Scholar] [CrossRef]
  67. Fortunato, S.; Nigro, D.; Lasorella, C.; Marcotuli, I.; Gadaleta, A.; de Pinto, M.C. The Role of Glutamine Synthetase (GS) and Glutamate Synthase (GOGAT) in the Improvement of Nitrogen Use Efficiency in Cereals. Biomolecules 2023, 13, 1771. [Google Scholar] [CrossRef]
  68. Li, H.; Hu, B.; Chu, C. Nitrogen Use Efficiency in Crops: Lessons from Arabidopsis and Rice. J. Exp. Bot. 2017, 68, 2477–2488. [Google Scholar] [CrossRef]
  69. Nagarajan, S.N.; Lenoir, C.; Grangeasse, C. Recent Advances in Bacterial Signaling by Serine/Threonine Protein Kinases. Trends Microbiol. 2022, 30, 553–566. [Google Scholar] [CrossRef]
  70. Beekmann, K.; de Haan, L.H.J.; Actis-Goretta, L.; van Bladeren, P.J.; Rietjens, I.M.C.M. Effect of Glucuronidation on the Potential of Kaempferol to Inhibit Serine/Threonine Protein Kinases. J. Agric. Food Chem. 2016, 64, 1256–1263. [Google Scholar] [CrossRef]
  71. Roskoski, R. Cyclin-dependent Protein Serine/Threonine Kinase Inhibitors as Anticancer Drugs. Pharmacol. Res. 2019, 139, 471–488. [Google Scholar] [CrossRef]
  72. Alhabbar, Z.; Yang, R.; Juhasz, A.; Xin, H.; She, M.; Anwar, M.; Sultana, N.; Diepeveen, D.; Ma, W.; Islam, S.; et al. NAM gene allelic composition and its relation to grain-filling duration and nitrogen utilisation efficiency of Australian wheat. PLoS ONE 2018, 13, e0205448. [Google Scholar] [CrossRef] [PubMed]
  73. Zougari, A.; Guy, S.; Planchon, C. Genotypic lipoxygenase variation in soybean seeds and response to nitrogen nutrition. Plant Breed. 1995, 114, 313–316. [Google Scholar] [CrossRef]
  74. Shang, C.; Wang, W.; Zhu, S.; Wang, Z.; Qin, L.; Alam, M.A.; Xie, J.; Yuan, Z. The responses of two genes encoding phytoene synthase (Psy) and phytoene desaturase (Pds) to nitrogen limitation and salinity up-shock with special emphasis on carotenogenesis in Dunaliella parva. Algal Res. 2018, 32, 1–10. [Google Scholar] [CrossRef]
  75. Zimmermann, S.E.; Benstein, R.M.; Flores-Tornero, M.; Blau, S.; Anoman, A.D.; Rosa-Téllez, S.; Gerlich, S.C.; Salem, M.A.; Alseekh, S.; Kopriva, S.; et al. The phosphorylated pathway of serine biosynthesis links plant growth with nitrogen metabolism. Plant Physiol. 2021, 186, 1487–1506. [Google Scholar] [CrossRef] [PubMed]
  76. Grubbs, R.C.; Trossbach, J.; Houghton, B.C.; Hitchcock, F.A. Effects of Folic Acid on Respiratory and Nitrogen Metabolism. J. Appl. Physiol. 1949, 2, 327–342. [Google Scholar] [CrossRef]
  77. Tong, S.; Zhao, L.; Zhu, D.; Chen, W.; Chen, L.; Li, D. From formic acid to single-cell protein: Genome-scale revealing the metabolic network of Paracoccus communis MA5. Bioresour. Bioprocess. 2022, 9, 55. [Google Scholar] [CrossRef]
  78. Feng, L.; Yang, T.; Zhang, Z.; Li, F.; Chen, Q.; Sun, J.; Shi, C.; Deng, W.; Tao, M.; Tai, Y.; et al. Identification and characterization of cationic amino acid transporters (CATs) in tea plant (Camellia sinensis). Plant Growth Regul. 2018, 84, 57–69. [Google Scholar] [CrossRef]
  79. The, S.V.; Snyder, R.; Tegeder, M. Targeting Nitrogen Metabolism and Transport Processes to Improve Plant Nitrogen Use Efficiency. Front. Plant Sci. 2021, 11, 628366. [Google Scholar] [CrossRef]
  80. Franklin, O.; Cambui, C.A.; Gruffman, L.; Palmroth, S.; Oren, R.; Näsholm, T. The carbon bonus of organic nitrogen enhances nitrogen use efficiency of plants. Plant Cell Environ. 2017, 40, 25–35. [Google Scholar] [CrossRef]
Figure 1. Soft-threshold power selection for gene co-expression network construction. Left panel: red line represents scale-free topology fit index across varying soft-threshold powers. Right panel: mean connectivity changes with increasing soft-threshold power.
Figure 1. Soft-threshold power selection for gene co-expression network construction. Left panel: red line represents scale-free topology fit index across varying soft-threshold powers. Right panel: mean connectivity changes with increasing soft-threshold power.
Genes 15 01505 g001
Figure 2. Correlation heatmap between gene co-expression modules and L2/L8 lines. Red indicates positive correlations and blue denotes negative correlations. Cells contain Pearson correlation coefficients and p-values in parentheses.
Figure 2. Correlation heatmap between gene co-expression modules and L2/L8 lines. Red indicates positive correlations and blue denotes negative correlations. Cells contain Pearson correlation coefficients and p-values in parentheses.
Genes 15 01505 g002
Figure 3. Heatmap of the gene expression of the eigengenes in the purple module. Color ranges from green (low) to red (high).
Figure 3. Heatmap of the gene expression of the eigengenes in the purple module. Color ranges from green (low) to red (high).
Genes 15 01505 g003
Figure 4. Gene expression of the eigengenes in the lightcyan module. Color ranges from green (low) to red (high).
Figure 4. Gene expression of the eigengenes in the lightcyan module. Color ranges from green (low) to red (high).
Genes 15 01505 g004
Figure 5. Gene expression of the eigengenes in the darkorange module. Color ranges from green (low) to red (high).
Figure 5. Gene expression of the eigengenes in the darkorange module. Color ranges from green (low) to red (high).
Genes 15 01505 g005
Figure 6. Gene expression of the eigengenes in the cyan module. Color ranges from green (low) to red (high).
Figure 6. Gene expression of the eigengenes in the cyan module. Color ranges from green (low) to red (high).
Genes 15 01505 g006
Figure 7. Gene expression of the eigengenes in the skyblue module. Color ranges from green (low) to red (high).
Figure 7. Gene expression of the eigengenes in the skyblue module. Color ranges from green (low) to red (high).
Genes 15 01505 g007
Figure 8. Expression profiles and the linear regression analysis of the four key genes. The expression profiles of (a) CsaV4_3G000307; (c) CsaV4_1G002492; (e) CsaV4_2G003460; and (g) CsaV4_7G001709. The linear regression analysis of (b) CsaV4_3G000307; (d) CsaV4_1G002492; (f) CsaV4_2G003460; and (h) CsaV4_7G001709. ‘r’ indicates the Pearson correlation coefficient and ‘p’ denotes significance. One ⭐ indicates a significance difference with p < 0.05, two ⭐ indicate p < 0.01, three ⭐ indicate p < 0.001, four ⭐ indicate p < 0.0001.
Figure 8. Expression profiles and the linear regression analysis of the four key genes. The expression profiles of (a) CsaV4_3G000307; (c) CsaV4_1G002492; (e) CsaV4_2G003460; and (g) CsaV4_7G001709. The linear regression analysis of (b) CsaV4_3G000307; (d) CsaV4_1G002492; (f) CsaV4_2G003460; and (h) CsaV4_7G001709. ‘r’ indicates the Pearson correlation coefficient and ‘p’ denotes significance. One ⭐ indicates a significance difference with p < 0.05, two ⭐ indicate p < 0.01, three ⭐ indicate p < 0.001, four ⭐ indicate p < 0.0001.
Genes 15 01505 g008aGenes 15 01505 g008b
Figure 9. (a) Protein interaction network for the gene CsaV4_3G000307. Nodes represent proteins, each encapsulating all proteins from a single gene locus, including splice isoforms and post-translational modifications. Edges denote specific and meaningful protein–protein interactions, highlighting the precision of the depicted associations. (b) Heatmap of Gene Enrichment Analysis. This illustrates False Discovery Rate (FDR) values across various categories, with a color gradient from lighter green (lower FDR values) to darker blue (higher values) indicating the significance of gene enrichment. Circles on the heatmap represent gene counts in each category, varying sizes correspond to gene numbers: smaller circles for two genes, larger circles for three genes.
Figure 9. (a) Protein interaction network for the gene CsaV4_3G000307. Nodes represent proteins, each encapsulating all proteins from a single gene locus, including splice isoforms and post-translational modifications. Edges denote specific and meaningful protein–protein interactions, highlighting the precision of the depicted associations. (b) Heatmap of Gene Enrichment Analysis. This illustrates False Discovery Rate (FDR) values across various categories, with a color gradient from lighter green (lower FDR values) to darker blue (higher values) indicating the significance of gene enrichment. Circles on the heatmap represent gene counts in each category, varying sizes correspond to gene numbers: smaller circles for two genes, larger circles for three genes.
Genes 15 01505 g009
Figure 10. (a) Protein structural and interaction network for the gene CsaV4_1G002492. Nodes are synonymous with proteins. (b) Heatmap of Gene Enrichment Analysis.
Figure 10. (a) Protein structural and interaction network for the gene CsaV4_1G002492. Nodes are synonymous with proteins. (b) Heatmap of Gene Enrichment Analysis.
Genes 15 01505 g010
Figure 11. Protein interaction network for the gene CsaV4_2G003460. (a) Nodes are synonymous with proteins; (b) Heatmap of Gene Enrichment Analysis.
Figure 11. Protein interaction network for the gene CsaV4_2G003460. (a) Nodes are synonymous with proteins; (b) Heatmap of Gene Enrichment Analysis.
Genes 15 01505 g011
Figure 12. (a) Protein interaction network for the gene CsaV4_7G001709. Nodes are synonymous with proteins. (b) Heatmap of Gene Enrichment Analysis.
Figure 12. (a) Protein interaction network for the gene CsaV4_7G001709. Nodes are synonymous with proteins. (b) Heatmap of Gene Enrichment Analysis.
Genes 15 01505 g012
Figure 13. Protein sequence and structural analysis in the high NUE line. (a) The protein sequence logo, where the conservation of amino acids at specific positions within the sequence is represented by the height of letters. The red box marks the start position of the frame-shift mutation. (b) The three-dimensional model of the protein structure. The accompanying green matrix represents the predicted interaction strength or structural stability between amino acid residues within the protein, with darker green indicating higher prediction accuracy or stability, and lighter green indicating lower accuracy or stability.
Figure 13. Protein sequence and structural analysis in the high NUE line. (a) The protein sequence logo, where the conservation of amino acids at specific positions within the sequence is represented by the height of letters. The red box marks the start position of the frame-shift mutation. (b) The three-dimensional model of the protein structure. The accompanying green matrix represents the predicted interaction strength or structural stability between amino acid residues within the protein, with darker green indicating higher prediction accuracy or stability, and lighter green indicating lower accuracy or stability.
Genes 15 01505 g013
Figure 14. Protein sequence and structural analysis in the low NUE line. (a) The protein sequence logo. The red box marks the start position of the frame-shift mutation. (b) The three-dimensional model of the protein structure.
Figure 14. Protein sequence and structural analysis in the low NUE line. (a) The protein sequence logo. The red box marks the start position of the frame-shift mutation. (b) The three-dimensional model of the protein structure.
Genes 15 01505 g014
Table 1. Gene information and functional annotation of the 25 selected genes.
Table 1. Gene information and functional annotation of the 25 selected genes.
Gene IDGene DescriptionOrganism SpeciesModuleRelated
References
CsaV4_1G000167BALDH—Benzaldehyde dehydrogenase,
mitochondrial
Antirrhinum majuslightcyan
CsaV4_1G002491TAL2—Transaldolase 2Streptomyces coelicolordarkorange
CsaV4_1G002492GLR22—Glutamate receptor 2.2Arabidopsis thalianalightcyan[43,44,45,46]
CsaV4_1G002522NAS4—Probable nicotianamine synthase 4Arabidopsis thalianalightcyan
CsaV4_1G003819TENAC—Bifunctional TH2 protein,
mitochondrial
Arabidopsis thalianapurple
CsaV4_2G000093HXK2—Hexokinase-2Arabidopsis thalianapurple
CsaV4_2G000502NAS3—Nicotianamine synthase 3Arabidopsis thalianalightcyan
CsaV4_2G000612MDL3— (R)-mandelonitrile lyase 3Prunus serotinadarkorange
CsaV4_2G003460GLR35—Glutamate receptor 3.5Arabidopsis thalianalightcyan[43,44,45,46]
CsaV4_2G003519AB14B—ABC transporter B family member 14Arabidopsis thalianalightcyan
CsaV4_3G000307NPF6.3—Protein NRT1/PTR FAMILY 6.3Arabidopsis thalianadarkorange[47,48,49,50,51,52,53,54,55,56]
CsaV4_3G001073AB9C—ABC transporter C family member 9Arabidopsis thalianalightcyan
CsaV4_3G003070SCP35—Serine carboxypeptidase-like 35Arabidopsis thalianapurple
CsaV4_3G004703SOX—Probable sarcosine oxidaseArabidopsis thalianapurple
CsaV4_3G004765ISOA3—Isoamylase 3, chloroplasticArabidopsis thalianapurple
CsaV4_4G000627ISPF—2-C-methyl-D-erythritol
2,4-cyclodiphosphate synthase, chloroplastic
Arabidopsis thalianapurple
CsaV4_4G000995RIBA1—Bifunctional riboflavin biosynthesis protein, chloroplasticArabidopsis thalianalightcyan
CsaV4_4G002807ZIP3—Zinc transporter 3Arabidopsis thalianalightcyan
CsaV4_5G000080CHIA—Acidic endochitinaseCucumis sativuspurple
CsaV4_6G000223PER52—Peroxidase 52Arabidopsis thalianadarkorange
CsaV4_6G000691RBL13—RHOMBOID-like protein 13Arabidopsis thalianapurple
CsaV4_6G002589ACCO1—1-aminocyclopropane-1-carboxylate oxidase 1Cucumis melolightcyan
CsaV4_6G002776PIP21—Probable aquaporin PIP2-1Oryza sativa subsp. japonicadarkorange
CsaV4_7G000732XYN4—Endo-1,4-beta-xylanase 4Arabidopsis thalianapurple
CsaV4_7G001709UPS2—Ureide permease 2Arabidopsis thalianaskyblue[57,58,59,60]
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Ma, L.; Wei, A.; Liu, C.; Liu, N.; Han, Y.; Chen, Z.; Wang, N.; Du, S. Screening Key Genes Related to Nitrogen Use Efficiency in Cucumber Through Weighted Gene Co-Expression Network Analysis. Genes 2024, 15, 1505. https://doi.org/10.3390/genes15121505

AMA Style

Ma L, Wei A, Liu C, Liu N, Han Y, Chen Z, Wang N, Du S. Screening Key Genes Related to Nitrogen Use Efficiency in Cucumber Through Weighted Gene Co-Expression Network Analysis. Genes. 2024; 15(12):1505. https://doi.org/10.3390/genes15121505

Chicago/Turabian Style

Ma, Linhao, Aimin Wei, Ce Liu, Nan Liu, Yike Han, Zhengwu Chen, Ningning Wang, and Shengli Du. 2024. "Screening Key Genes Related to Nitrogen Use Efficiency in Cucumber Through Weighted Gene Co-Expression Network Analysis" Genes 15, no. 12: 1505. https://doi.org/10.3390/genes15121505

APA Style

Ma, L., Wei, A., Liu, C., Liu, N., Han, Y., Chen, Z., Wang, N., & Du, S. (2024). Screening Key Genes Related to Nitrogen Use Efficiency in Cucumber Through Weighted Gene Co-Expression Network Analysis. Genes, 15(12), 1505. https://doi.org/10.3390/genes15121505

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop