Comparative Proteomic Analysis Unveils Critical Pathways Underlying the Role of Nitrogen Fertilizer Treatment in American Elderberry.

American elderberry (Sambucus nigra subsp. canadensis) is a rapidly growing specialty crop in Missouri and eastern North America. Nitrogen (N) is a major nutrient involved in plant growth and development. However, proteome changes for different genotypes of elder in response to varying levels of N-treatment remain undefined. To reveal plant responses to N, comparative proteomic analyses were performed to determine consistent changes in three genotypes of elderberry leaves (Adams II, Bob Gordon and Wyldewood) grown under different N-fertilizer treatments. 165 proteins separated by two dimensional gel electrophoresis showed significant differences in abundance (p < 0.05 and greater than 2-fold). Principal component analysis of the abundance profiles of these proteins revealed Bob Gordon as a distinct genotype. The 165 proteins were identified by mass spectrometry and showed similar functional distributions in these genotypes underlying the N-treatment. Among the proteins identified, 23 are mainly involved in photosynthesis, protein metabolism and redox homeostasis. Their abundance profiles were not altered upon exposure to N or genotype. These results provide novel insights into plant responses to fertilizer treatment at the proteome level and could lead to a better understanding of molecular mechanisms of elderberry growth.

1. Introduction

American elderberry (Sambucus nigra subsp. canadensis) is a specialty crop grown in many parts of North America. Overall elderberry cultivation is expanding in the U.S. with the potential for further growth [1]. The cultivation of elder plants requires the use of fertilizers; most commonly nitrogen (N)–based. Therefore, improving and understanding the role of N in elder production is needed [2]. Fertilizers typically contain N in the form of nitrate and ammonium compounds, which rapidly dissolve in water, resulting in loss to the soil, and low nitrogen use efficiency (NUE) [3]. A better understanding of the molecular interactions and mechanisms that underlie elder response to different N-treatments may aid in more efficient utilization of N-fertilizer as well as furthering the understanding of elderberry growth.

Several studies have examined gene expression changes in plants in response to N treatment. Bi and co-workers found a set of N-responsive genes in three genotypes of maize under N-limiting conditions by high throughput RNA sequencing [4]. Transcriptome analyses have identified differentially expressed genes under N-deprivation conditions and yielded a complex pathway in the model plant Arabidopsis thaliana [5]. However, the exact function of these genes remains unknown. Proteomic approaches can aid in the understanding of changes and functions of the plant proteome to N [6,7,8,9]. Proteomics can complement many functional genomics approaches [10] and provide useful information about post translational modifications (PTM) with annotation of the corresponding gene [11]. While several studies have found a proportional relationship between transcript levels and the abundance of proteins [12,13,14], it is still difficult to predict the relationship between cellular protein concentrations and the abundance of corresponding mRNAs [15]. Part of the difficulty is likely due to PTMs and the effects of micro RNAs on expression of proteins [16]. A comparative proteomics approach could provide a better assessment of differentially abundant proteins and metabolic processes.

Here, we present a comparative proteomics investigation of the impact of N-fertilizer treatment on American elderberry (Sambucus nigra subsp. canadensis). Three genotypes are used in the investigation; including an heirloom genotype (Adams II) and the more recently developed genotypes (Bob Gordon [17] and Wyldewood [18]). The selection of suitable N-fertilizer usage is important to ensure agriculture productivity and sustainability [2]. Elder leaves were chosen for the study due to their sensitivity to varying environmental conditions [19]. Two-dimensional electrophoresis (2-DE), followed by in-gel digestion, and mass spectrometry (MS) allowed us to identify 165 proteins that showed altered abundance upon exposure to different levels of N-fertilizer treatment. Protein functional studies helped to further our understanding of protein profile changes in response to treatment. These results provide novel insight into cellular changes which could lead to new N-treatment management strategies to improve NUE and plant production in the future.

2. Materials and Methods

Acrylamide, phenol (Tris-equilibrated), CHAPS, 3-(4-heptyl) phenyl-3-hydroxypropyl) dimethylammoniopropanesulfonate (C7BzO), dithiothreitol (DTT), iodoacetamide (IAA), EDTA and 2-hydroxyethyl disulfide (2-HED) were obtained from Sigma-Aldrich (St. Louis, MO, USA). Ammonium persulfate, bis-acrylamide, urea (electrophoresis grade), thiourea (electrophoresis grade), ammonium acetate, glycerol, Coomassie Brilliant Blue G-250, agarose (low EEO) and 2-mercaptoethanol were obtained from Thermo Fisher Scientific (Houston, TX, USA). N,N,N′,N′-tetramethylethylenediamine (TEMED) was obtained from BioRad Laboratories, Inc. (Hercules, CA, USA). EZQ protein quantification kits were obtained from Molecular Probes (Invitrogen, Carlsbad, CA, USA). IPG strips and IPG buffer 3–10 were obtained from GE Healthcare (Piscataway, NJ, USA). Modified porcine trypsin was acquired from Promega (Madison, WI, USA).

Study materials (leaves) were harvested from an active research trial in southwest Missouri that incorporated three elderberry genotypes and four N-fertilization treatments (applied 9 May 2014 and the two previous springs), as previously described [20]. Fresh, fully mature elder leaves from the middle portions of branches were randomly collected from 36 of the 48 plots in August 2014, promptly refrigerated, and transported under refrigeration to the laboratory for analysis. About 1000 cm3 of leaves were collected from each plot. The leaves were collected after fruit had been harvested. The middle portions of the compound elder leaves were then washed, ground under liquid nitrogen, and 100 mg of each replicate was aliquoted into 2 mL microcentrifuge tubes and stored at −80 °C prior to protein extraction.

The proteins were extracted according to a method described by Hurkman and Tanaka with some modifications [21]. Briefly, proteins were isolated by phenol and precipitated overnight by five volumes of pre-chilled 0.1 M ammonium acetate in 100% methanol. Protein concentration was determined by EZQ quantitation kits (Life Technologies, Grand Island, NY, USA) using the manufacturer’s instructions. Proteins were separated by 2-DE as previously described [22]. All 2-DE stained gels were scanned using a UMAX PowerLook 2100 XL scanner (UMAX Systems GmbH, Willich, Germany). At least triplicate gels were performed for each treatment, and a total of 36 gel images were analyzed by ImageMaster 2D Platinum software version 7.0 (GE Healthcare). To compensate for gel-to-gel variation, the volume of each spot was normalized as a relative volume as previously described [23].

Three replicates were performed for each genotype at each N-fertilizer treatment. One-way analysis of variance (ANOVA), part of ImageMaster 2D Platinum (Version 7.0; GE Healthcare), was initially performed for each group to identify statistical significance by using the relative volume (vol %). A significance level of 95% was chosen. To evaluate the quality of spots, several parameters were set manually. The 2-fold difference calculated as the ratio of the average spot volume of each treatment to that of the control was used; spots were observed in at least 80% of spot maps and repeated on triplicates. Once the spots were considered to be statistically significant, the average values of these spots were analyzed by Tukey’s honest significant difference test, and then filtered by two-way ANOVA by genotype, two-way ANOVA by N-treatment, and/or two-way ANOVA interaction (p < 0.05) using the SAS software package (Version 9.4, SAS Institute, Cary, NC, USA). Only spots of interest were selected for further MALDI-TOF/TOF or LC-MS/MS analysis. Additionally, principal component analysis (PCA) was performed using the SAS software.

Six pooled gels were generated by pooling equal amounts (200 µg) of three biological replicates onto one gel following the same steps as described above. Each protein spot of interest was excised from the pooled gel, trypsin digested and identified by MALDI-TOF/TOF (4700 Proteomics Analyzer, Applied Biosystems, Foster City, CA, USA) in the positive ion reflector mode as previously described [22]. For peptide and protein identification of MALDI-TOF/TOF data, the resulting peptide peak lists were submitted to the MASCOT database search engine against the National Center for Biotechnology Information non-redundant (NCBInr) database. A homology search was performed due to the limited sequences of Sambucus nigra genome (164 sequences as of May 2018). The following parameters were selected: Viridiplantae (Green Plant) as taxonomy, trypsin as digesting enzyme, one missed cleavage, fixed modification of carbamidomethyl (Cys), variable modification of oxidation (Met), precursor ion mass error tolerance of 100 ppm, MS/MS fragment ion mass error tolerance of 0.1 Da, peptide charge of 1+, monoisotopic and MALDI-TOF/TOF as instrument. Confident protein identifications were defined as: (1) the highest protein score on the database searching report, (2) a minimum of two matched peptides, (3) less than 15% deviation between theoretical and experimental Mr and pI values (gel-based). An in-house BlastP search at NCBI was performed to verify all matches and update annotations and identification of all hypothetical or unknown proteins (sequence similarity > 80%).

In cases where the MASCOT search did not reveal confident identification, we performed LC-MS/MS experiments as previously described [24]. For peptide and protein identification by LC-MS/MS, raw files were analyzed and quantified using the SEQUEST algorithm implemented with a Scorcerer2 integrated Data Appliance (SageN Research, Milpitas, CA, USA) against a local copy of the Viridiplantae database downloaded in FASTA format via file transfer protocol from NCBInr (released in November, 2013; 2,355,794 proteins). Homology searches included full trypsin specificity (KR/P), two missed cleavage sites, peptide mass tolerance of 50 ppm, fragment ion mass tolerance of 1 Da, carbamidomethylation of Cys, as a static modification and Met oxidation as variable modification. Search results were verified by Scaffold viewer V4.0.5 (Proteome Software Inc., Portland, OR, USA). For each identification, a minimum of two matched peptides; peptide threshold (p < 0.05) and protein threshold (p < 0.001); cross-correlation factor (Xcorr) 2.0, 3.0 and 3.5 for the charge states and +2, +3 and +4 respectively, and minimum Delta CN (Delta correlation) of 0.1. Results were then filtered using the Protein- and PeptideProphet [25,26] implemented in the Scaffold software to achieve a peptide and protein global false discovery rate of <5 and 0.1%, respectively. An extra BlastP search was also performed to update the best matches and align all hypothetical or unknown proteins (sequence similarity > 80%).

Hierarchical clustering was constructed using the software PermutMatrix [27]. The average abundance of identified proteins, which were classified in the same orthologous group and displayed similar trends of abundance profiles across all genotypes, were taken and the fold change calculated between each group and control, followed by log2 transformation for heat map representation. The dissimilarities were calculated based on Euclidean distances and hierarchical clustering was carried out according to Ward’s method. The trees were generated using the multiple-fragment heuristic algorithm (MF) as a seriation rule.

The function of the identified proteins was sought by transferring original sequences to the Arabidopsis genome and orthologous genes using the Mercator web-based pipeline (http://mapman.gabipd.org/app/mercator), which divides protein into 35 hierarchical, nonredundant functional classes using MapMan bin codes [28]. The cellular locations of identified proteins were determined by searching the best matched Arabidopsis orthologous proteins in the SUBA database (Version 3) [29].

3. Results

3.1. Comparative Proteome Analyses and Differentially Abundant Proteins in Elderberry Leaves

To investigate the effect of N-fertilizer treatment upon protein profiles of elder leaves, three different applications were applied to the plants (56, 112 and 169 kg N/ha). Control plants had no N-fertilizer applied. The three N-fertilizer applications can be thought of as a low, medium, and high concentration. 2-D gel electrophoretic separation of the proteins extracted from leaves of three genotypes (Adams II, Bob Gordon and Wyldewood) of American elderberry was performed. A combination of 12 treatments of genotype × N-concentration as a factorial design was developed, along with three replicates derived from separate leaf extractions giving a total of 36 images. The correlation coefficients among the biological triplicates were shown to vary from 0.9863 to 0.9993 (Figure S1). Replicate gel images were then digitally integrated to build a “master gel” for each treatment condition. Three different groups were established to identify differentially abundant proteins; Group 1: 0 (or control) vs 56 kg N/ha, Group 2: 0 vs 112 kg N/ha, and Group 3: 0 vs 169 kg N/ha. On average, 300 protein spots were initially detected on each 2-DE gel with some distinctive changes within genotypes or N-fertilizer conditions. A one-way ANOVA analysis was applied to all matched spots in each group (p < 0.05). Once the protein spots were determined to be differentially abundant, the spots were selected using the following stringent criteria; at least 2-fold change, correlation of all spot maps, and appearance on all replicates. Of the 3564 spots initially detected among the master gels, only 165 “high class” protein spots fit these criteria (Table 1).

Table 1

Frequency of differentially abundant protein spots detected from 2-DE gel analysis among three elderberry genotypes (Adams II, Bob Gordon and Wyldewood) and three N fertilizer treatments a.

N Concentration Comparison	Initial Detected Protein Spots b	Differentially Abundant Protein Spots c	High Class Spots d	Effect e				MS/MS Analyzed Spots f
				G	N	G and N	G × N
0 vs 56	306	86	61	30	4	11	16	28
0 vs 112	295	69	53	21	2	11	19	37
0 vs 169	281	68	51	19	4	8	20	36
Total	882	223	165	70	10	30	55	101

a Protein spots that exhibited differential abundance were identified by MALDI TOF/TOF or LC-MS/MS. b Protein spots were detected by ImageMaster. c The number of protein spots was determined by a one-way ANOVA analysis. d Manually excluded spots that had less than a 2-fold change or absence within any of the triplicate gels. e Four effects were separated based on p values (G: genotype as main effect; N: N as main effect; G and N: genotype and N as main effect, but no interaction; G × N: interaction effect between genotype and N.) f Included those spots of interest that were excised, digested and analyzed by MS.

For the 165 protein spots which met the above criteria, a post hoc Tukey test was conducted to determine statistical differences in abundance changes between each treatment and control. A two-way ANOVA analysis was performed on the 165 protein spots to assess any significant differences due to either genotype or N-treatment on protein abundance (Table S1). The results of the two-way ANOVA analysis revealed that genotype had a significant effect on the expression of 70 protein spots, whereas N-treatment had no effect. For 10 of the protein spots N-treatment resulted in significant changes in abundance, but genotype had no effect on these spots. Both genotype and N-treatment had a significant effect on the expression level of 30 protein spots. Genotype and N-treatment had an interaction effect on the abundance level of 55 protein spots. Based on the analysis, we were able to identify significant protein abundance level changes and the conditions resulting in abundance change.

The representative 2-DE gel maps are shown in Figure S2. Pick locations of spots that have significant volume changes based on the variance of group means are labeled. The proteins that have same location (Mr and pI) throughout three groups are considered to be the same proteins which reduced the number of spots for MS analysis to 101.

3.2. Protein Identification

15 proteins were successfully identified by MALDI-TOF/TOF, whereas the remaining 86 proteins were determined by LC-MS/MS (Tables S2 and S3). Notably, 39 proteins were interpreted either as unknown and hypothetical proteins, or proteins without specific function in the database. This indicates that the comprehensive function of many N-responsive proteins from Sambucus nigra subsp. canadensis have yet to be elucidated. To gain further information about these unknown proteins, we performed a BLASTP search of their homologues using original sequences as queries. The corresponding homologues with the highest homology are shown in Tables S4–S6. All 39 proteins shared more than 80% positive identity with homologues at the amino acid level, showing that they might have similar function. The relative expression level for the 101 identified proteins are also listed in Tables S4–S6. In particular, 63 proteins displayed different experimental Mr and/or pI and, confirmed by the identification of at least two unique peptides, were identified as protein isoforms. For example, isoforms were observed for the 20 kDa chaperonin family protein (CPN21; spots 49, 56, 61 and 90), the ATP synthase CF1 α subunit (ATP-SCF1S; spots 13, 14 and 15), and catalase (CAT; spots 18 and 19).

3.3. Hierarchical Clustering of Differentially Abundant Proteins in the Three Treatment Groups

Overall, there are 52 up-regulated spots and 73 down-regulated spots for Adams II. 42 up-regulated spots and 71 down-regulated spots are detected for Bob Gordon. 55 up-regulated spots and 70 down-regulated spots are detected for Wyldewood. The number of differentially abundant proteins for down-regulated was greater than up-regulated across the three genotypes at each N-fertilizer rate (except Adams II at 169 kg N/ha and Wyldewood at 112 kg N/ha). For Bob Gordon, there were about twice as many down-regulated differentially abundant proteins, indicating that Bob Gordon is dissimilar to Adams II or Wyldewood.

Protein abundance analysis for a given fertilizer treatment and genotype can highlight proteins which may respond to N during plant growth and development. 13 proteins showed the same trend of regulation (up or down regulation) for a given N-treatment across the three genotypes. In Group 1, ATP-dependent zinc metalloprotease FTSH 1 (FTSH1; spot 11), CA2 (spot 62) and actin 11 (ACT11; spot 86) were down-regulated. In Group 2, eukaryotic initiation factor 4A-8 (EIF4A8; spot 22) and citrate synthase family protein (CS; spot 23) were up-regulated, while heat shock protein 70 (HSP70; spot 5), 30S ribosomal protein S1 (RPS1; spot 21), probable plastid-lipid-associated protein 6 (PAP6; spot 53) and AMT (spot 87) were down-regulated. In Group 3, 20S proteasome subunit PBA1 (PBA1; spot 98) was up-regulated. In contrast, V-type proton ATPase catalytic subunit A (VHA-A; spot 12), sedoheptulose-1,7-bisphosphatase (SHB; spot 97) and RPI (spot 99) were down-regulated (Table S1).

Also, a hierarchical clustering of 165 differentially abundant proteins was performed with the aim of evaluating similarity of protein regulation in leaf organs across the different genotypes in each N-treatment group as shown in Figure 1. In Figure 1a, the low N-treatment, 61 proteins were distributed in four main clusters. Cluster A included nine proteins abundances that were mainly down-regulated in the three genotypes and all Adams II protein abundances in this cluster are down-regulated. Cluster B comprised nine protein groups with the largest abundance values of the genotypes. In Cluster D of Figure 1a, most of the protein abundances are down-regulated and all the Bob Gordon proteins abundances in this cluster are down-regulated. In Figure 1b, the medium N-treatment, 53 proteins fell into three clusters. Cluster B included seven proteins and represents the largest expression level. Cluster C proteins were mainly down-regulated in the three genotypes, but the downregulation in Bob Gordon was the strongest. In Figure 1c, the high N-treatment, two main clusters were recognized based on 51 protein expression profiles. Cluster A was mainly down-regulated in three genotypes, while cluster B was up-regulated, where the Adams II and Bob Gordon expressions were larger than that of Wyldewood.

Figure 1

(a–c) Hierarchical clustering of 165 differentially abundant proteins (p < 0.05) for each N fertilizer concentration (number following A (Adams II), B (Bob Gordon), W (Wyldewood)). The abundance of each protein was transformed using log2. Every colored box represents a protein that is up-regulated (red) or down-regulated (green) in a certain treatment. Gray color means + ∞. Hierarchical clusters (A, B, C, and D) were established as described in the text.

3.4. Principle Component Analysis (PCA) of Each Group

Global correlation of each group is visualized by plotting the two largest sources of variance, PCA1 and PCA2, which orthogonally divide the spots maps among expression data (Figure S3). In general, the first two components of Group 1 (3 genotypes, low N-treatment), Group 2 (3 genotypes, medium N-treatment), and Group 3 (3 genotypes, high N-treatment) contained 93.64%, 95.04% and 94.32% of the total data variance, respectively. The largest variance (PCA 1) represents technical and biological variation. The second largest variance (PCA 2) represents applied treatment variation. PCA reveals that the largest variation is not due to genotype or N-fertilizer treatment. PCA of the 165 differentially abundant proteins in the three groups demonstrated that, with and without N fertilizer treatment, Bob Gordon has a protein expression profile distinct from that of Adams II and Wyldewood. This result suggests that Bob Gordon is a distinct phenotype. The correlation of Adams II and Wyldewood demonstrated a great similarity among different N fertilizer treatments. However, PCA could not distinguish between control and each N-treatment of Adams II, Bob Gordon or Wyldewood. The same trends are shown across the three groups in the PCA analysis, demonstrating a satisfactory reproducibility of the experimental processes.

3.5. Functional Classification and Cellular Localization of Differentially Abundant Proteins

The 101 identified proteins were classified into functional groups including the largest group photosynthesis (30.7%), followed by protein metabolism (15.8%) and glycolysis (10.9%), redox homeostasis (8.9%), amino acid metabolism (8.9%), tricarboxylic acid (TCA) cycle (5.9%), oxidative pentose phosphate (OPP) (4.0%) pathway and miscellaneous (14.9%) (Figure 2). Proteins involved in photosynthesis were classified into two categories, including Calvin cycle and light reaction. Proteins involved in protein metabolism included those involved in protein assembly and cofactor ligation, folding, synthesis, degradation and PTM. The 101 identified proteins were assigned to six major cellular locations (Tables S4–S6). Most proteins were found in the chloroplast and cytosol. The remaining proteins were discovered in the mitochondrion, peroxisome, Golgi apparatus and/or nucleus. Both functional classification and cellular localization results show that a large number of chloroplast proteins are involved in photosynthesis, which were likely greatly affected by N-treatment.

Figure 2

The functional category distribution of the 101 identified proteins. Proteomes identified in the leaf were matched to their Arabidopsis orthologues and classified into 16 major functional groups using the Mercator software (gabipd.org/biotools/mercator/). Proteins involved in photosynthesis and protein metabolism were further sub-classified as indicated by arrows. OPP: oxidative pentose phosphate; TCA: tricarboxylic acid; PTM: post translational modification.

Using the same approach, similar distributions were observed for up- or down-regulated differentially abundant proteins for each genotype and N-treatments (Figure S4). The functional classification of N characteristic-proteins exhibited expression profiles associated with key roles in photosynthesis, protein metabolism, glycolysis and redox homeostasis. The functional classes are also included for each identified protein in Table S4.

3.6. Expression Profiles of Unique Proteins in Response to N Across Three Genotypes

A total of 23 unique proteins representing 58 isoforms appeared in all three groups. This profile was not altered upon exposure to N or genotype. Table 2 summarized the identified proteins, spot number, and functional class. Not surprisingly, photosynthesis and protein metabolism comprise the majority of functional classes for the leaves.

Table 2

Expression profiles of unique proteins in response to N across three genotypes.

Protein Name	Spot Number	Functional Class
RubisCO activase α form precursor (RCA)	20 and 35	Photosynthesis (34.8%)
Glyceraldehyde-3-phosphate dehydrogenase B (GAPB)	26, 28 and 75
Phosphoglycerate kinase (PGK)	27, 30 and 31
Oxygen-evolving enhancer protein 1 (OEE1)	41, 44, 57 and 58
Chlorophyll a–b binding protein of LHCII type 1-like (LHCII type1)	46
Chloroplast pigment-binding protein CP26 (CP26)	48 and 52
Ferredoxin-NADP reductase, leaf isozyme (LFNR)	64, 69, 70, 91 and 92
Glyceraldehyde-3-phosphate dehydrogenase A (GAPA)	65, 66, 73, 76 and 101
Chaperonin 60 subunit β 1 (CPN 60)	10	Protein Metabolism (17.4%),
Elongation factor-GTP-binding domain-containing protein (EF)	32
Photosystem II stability/assembly factor HCF136 (HCF136)	37 and 40
20 kDa Chaperonin family protein (CPN21)	49, 56, 61 and 90
Catalase (CAT)	18 and 19	Redox Homeostasis (13.0%)
Monodehydroascorbate reductase family protein (MDAR)	29
Ascorbate peroxidase (APX)	51, 54, 55 and 59
Phosphoglucomutase (PGM)	6, 7, 8 and 9	Glycolysis (8.7%)
Glyceraldehyde-3-phosphate dehydrogenase C subunit1 (GAPC1)	74, 77, 78, 81 and 89
Transketolase family protein (TK)	2, 3 and 4	OPP pathway (4.3%)
Glutamine synthetase cytosolic isozyme 1-1-like (GS)	33 and 34	N-metabolism (4.3%)
Putative lactoylglutathione lyase (LGL)	36	Biodegradation of xenobiotics (4.3%)
4-Nitrophenylphosphatase (4-NPPase)	43	Minor CHO metabolism (4.3%)
Guanine nucleotide-binding protein subunit β-like protein (GNBP)	72	Development (4.3%)
Monomeric lectin SNAlm precursor (SNAlm)	42	Not assigned (4.3%)

Among these, 18 unique proteins appeared as 23 isoforms which agrees well with the two-way ANOVA results (genotype × N interaction, p < 0.05). In some instances, the same proteins (or isoforms) were assigned to different effects across three groups, suggesting the concentration may change. Thus these consistently shown proteins are represented in important pathways and functional classes. We performed hierarchical clustering to compare trends of proteins enrichment in response to all N treatments among the different genotypes. These protein groups fell into four main clusters (Figure 3). Four photosynthesis related proteins (LHCII type 1, PGK, RCA and HCF136) were distributed in cluster A along with other proteins such as EF, TK and CPN21. These proteins were particularly up-regulated in Wyldewood at 56 and 169 kg N/ha. Cluster B comprised most glyceraldehyde-3-phosphate dehydrogenase proteins, CAT, 4-NPPase, one isoform of GS and GNBP, which were mainly down-regulated across genotypes and N treatments. Cluster C contained LFNR and GAPA proteins which were mainly up-regulated at 169 kg N/ha regardless of genotype analyzed. Finally, in cluster D, APX, GAPB, SNAlm, one isoform of GS, one isoform of CPN21 and LGL were grouped. Also, expression profiles from all genotypes and N treatments were clustered and are represented at the top of the columns in Figure 3. Bob Gordon under different N-treatments were clustered together, indicating similarity of expression of Bob Gordon among N fertilizers. In summary, the global abundance of proteins related to OPP pathway, detoxification, carbohydrate recognition, protein metabolism (elongation, folding and assembly) and redox homeostasis were upregulated at 169 kg N/ha compared to control, while proteins associated with development, glycolysis, photosynthesis, minor CHO metabolism and ammonium assimilation were down-regulated.

Figure 3

Hierarchical clustering of orthologous N characteristic-protein groups consistently enriched in three N groups compared to control (p < 0.05) in different genotypes. The color gradient range of the heat map indicates the proportional up-regulation (red) and down-regulation (green) of protein abundance (fold difference, log2-transformed). Hierarchical clusters (A, B, C, and D) were established as described in the text.

Figure 4 summarizes the major pathways identified and the protein abundance changes with treatment. Our results showed that 13 proteins representing 21 isoforms had similar up- or down-regulated change patterns in response to all N-treatment groups. For example, the protein GLO (spot 36) was up-regulated in both Adams II and Bob Gordon. In contrast, CPN 60 (spot 10) and LHCII type 1 (spot 46) were down-regulated in both Adams II and Bob Gordon. For each genotype, three isoforms of PGK (spots 27, 30 and 31) and EF (spot 32) were down-regulated in Adams II. Two isoforms of CAT (spots 18 and 19), MDAR (spot 29), two isoforms of HCF136 (spots 37 and 40), 4-NPPase (spot 43) and GNBP (spot 72) were down-regulated in Bob Gordon. Three isoforms of GAPB (spots 26, 28 and 75), two isoforms of GS (spots 33 and 34) and two isoforms of CP26 (spots 48 and 52) were down-regulated in Wyldewood (Table S1).

Figure 4

Schematic representation of key pathways that regulate N fertility. These differentially abundant proteins (interaction p < 0.05) across different genotypes and N treatments were also classified in functional categories using MapMan bin codes (http://mapman.gabipd.org). The cellular locations of each protein were determined by the SUBA database (Version 3). Arrows indicate up- and down-regulation of proteins within each set of isoforms compared to control group at a given N condition. Red squares in the functional category indicate an overall regulation of the category at 169 kg N/ha compared to control (0 kg N/ha) in different genotypes, while green squares means down regulation.

On the other hand, nine proteins appeared as 36 isoforms, and consistently exhibited opposite expression patterns within each set of isoforms across all genotype and N groups. These proteins included three isoforms of TK (spots 2, 3 and 4), four isoforms of PGM (spots 6, 7, 8 and 9), two isoforms of RCA (spots 20 and 35), four isoforms of OEE1 (spots 41, 44, 57 and 58), four isoforms of APX (spots 51, 54, 55 and 59), five isoforms of LFNR (spots 64, 69, 70, 91 and 92), five isoforms of GAPA (spots 65, 66, 73, 76 and 101), five isoforms of GAPC1 (spots 74, 77, 78, 81 and 89) and four isoforms of CPN21 (spots 49, 56, 61 and 90). In addition, one protein, SNAlm (spot 42), which was identified in European elder (Sambucus nigra subsp. nigra) [30], showed opposite expression patterns in different N groups regardless of genotype analyzed. Three isoforms of GAPB (spots 26, 28 and 75), two isoforms of GS (spots 33 and 34) and two isoforms of CP26 (spots 48 and 52) were also consistently observed to have opposite abundance levels for Adams II and Bob Gordon. Two isoforms of CAT (spots 18 and 19) and two isoforms of HCF136 (spots 37 and 40) were found for Adams II and Wyldewood. Three isoforms of PGK (spots 27, 30 and 31) were represented for Bob Gordon and Wyldewood (Table S1).

4. Discussion

It is essential to develop a more comprehensive understanding of how different plants respond to different levels of N-fertilizer in the field. In this study, we exploited a comparative proteomics approach to identify specificity and dynamic changes in leaf proteomes of American elderberry affected by N-treatment (Table 1). Protein functional classification with their Arabidopsis orthologues was adopted to enable us to evaluate the expression level of the same proteins in different genotypes. Understanding how American elderberry responds to different levels of fertilizer not only permits us to better manage its utilization, but also to understand the growth and development of elder, possibly improving the yields of different genotypes with optimized choice of N-fertilizer usage.

2-DE coupled with MALDI-TOF/TOF and LC-MS/MS and advanced statistical methods were used in this study as an efficient proteomics approach, allowing us to observe spatial difference in protein expression. 165 protein spots (approximately 5% of all detectable protein spots) showed significant changes in response to N-treatment compared to control. PCA analysis revealed that the protein expression profile of Bob Gordon differs from that of the other two genotypes, indicating that Bob Gordon is a distinct phenotype and sensitive to N-treatment. This is consistent with the observation that Bob Gordon has the highest yield among the three genotypes under same growing conditions [31]. Our study suggests that Bob Gordon may possess widespread usefulness as a commercial elder cultivar and is worth further study.

According to the results of the two-way ANOVA analysis (Table S1), some of the differentially abundant proteins were assigned to N-treatment as the main effect. Phenotypically, three genotypes behave differently when cultivated side-by-side, suggesting that the different inherent genetic backgrounds play a role [31]. However, some differentially abundant proteins still had the same pattern of regulation across different genotypes with the same N-treatment indicating common pathways among genotypes.

Proteins tied to photosynthesis were differentially abundant in the study. These proteins have been previously highlighted in different species upon varied N-treatment [9,32,33,34,35,36,37]. Most of these photosynthetic proteins showed a lower expression level with increasing N-fertilizer concentrations, regardless of genotype (Figure 4). Our findings are consistent with earlier work which showed that the Rubisco activation state decreased at high N-concentrations and that Rubisco may serve as a storage protein in apple leaves [38]. A previous study on maize found lower intercellular CO2 concentrations in response to high N supply also consistent with our results [39]. A recent genetic study showed that the OsPIPs gene may be responsible for mediating insufficient CO2 assimilation [40]. We speculate one of the outcomes of decreased level of RCA in elderberry leaves with a high N-treatment could lower the CO2 concentration in the chloroplasts, that may be an indication of over-fertilization.

In addition to photosynthetic proteins, we could observe a considerable distribution of proteins involved in protein metabolism (synthesis, PTM, degradation, folding, assembly and cofactor ligation). Several components of the protein synthesis machinery (translation elongation and initiation factors and ribosomal proteins) except EF, were only expressed at a treatment of 112 kg N/ha. Elongation factors are known for regulating the translocation step in polypeptide chain elongation. The VAJ/GFA1/CLO gene encodes a translational elongation factor-2 family protein, which controls directional floral organ growth in Arabidopsis thaliana possibly through pre-mRNA splicing [41]. Similarly, eukaryotic initiation factors (EIF) are phosphorylated to modulate the rate of mRNA binding [42]. In tobacco, the expression of elongation factor EF-1 alpha increased in meristems, rapidly growing tissues, and developing gametophytes [43]. In our study, EF2 and EIF4A8 were found to be enriched in Bob Gordon and Wyldewood, suggesting enhanced protein synthesis to support the development and growth of plant. In Arabidopsis thaliana, mutations in genes encoding ribosomal proteins affect leaf developmental processes [44]. We suspect deficiency of the EF protein may play a role in perturbing the development of leaf under high N-treatment. Chaperone proteins (CPN21 and 60) are essential for facilitating protein folding, unfolding and preventing protein aggregation [45]. Others suggest it may regulate increased oxidative stress [46]. Consequently, the increased amounts of CPN21 may help to lower oxidative stress at a treatment of 112 kg N/ha. A comparative genetics study of mutants in Arabidopsis thaliana has implicated an inverse relationship between 26S proteasome (26SP) and 20S proteasome (20SP) abundance, with defectiveness 26SP level invariably causing the elevated abundances of free 20SP lead to increased tolerance to oxidative stress [47]. Our results agree well with the established correlation between 26SP and 20SP. We constantly observed enhanced expression of 20SP concurrently with the down-regulation of 26SP at different N levels. The enriched abundance of serine/threonine-protein phosphatase PP2A-4 catalytic subunit (PP2A-4) was presented for Adams II and Wyldewood at a treatment of 169 kg N/ha. Many functional roles of PP2A in stress signal pathways were discussed under biotic and abiotic conditions [48]. Its role in N stress is yet to be determined.

Proteins involved in amino acid metabolism, such as ketol-acid reductoisomerase (KARI), s-adenosylmethionine synthase 1, aspartate aminotransferase (AST), and GDH, were differently abundant in the three N groups. Amino acid biosynthesis and degradation tightly correlate with N-availability [49]. These proteins are involved in the biosynthesis of the branched-chain amino acids [50], glycine, serine, aspartate [51] and cysteine, respectively. KARI was also identified in maize ear with 2-DE gels under different N-conditions, which correlates with ear growth and grain yield [52]. In the case of AST, it is involved the signal transduction pathway to regulate asparagine to glutamine ratio in cobs and that ratio shows the plant N status [53]. Additionally, the ascending amount of AST likely coordinates with increasing maize productivity [54]. Taken together, the increased amount of these two proteins in Bob Gordon at a treatment of 169 kg N/ha may help explain the higher yield that is often observed under field conditions [31]. GDH and AMT, which are involved in amino acid degradation were found down-regulated for all three genotypes at a treatment of 112 kg N/ha.

Several proteins involved in redox homeostasis were found to exhibit altered expressions, including CAT, GDP-mannose 3,5-epimerase isoform X2 (GME), MDAR, 2-cys peroxiredoxin BAS1 (BAS1) and APX. Redox changes show well-organized interactions to prevent tissue damage and death under environmental disturbance [55]. Ascorbate-glutathione cycle plays a key role in H2O2 metabolism in plant cells [56]. H2O2 is reduced to water by APX using ascorbate as the electron donor and the reducing agent ascorbate is regenerated by MDAR and dehydroascorbate reductase using reduced glutathione [57]. The increased levels of APX and MDAR for Adams II were found at a treatment of 169 kg N/ha. Similarly, GME is involved in the regulation of ascorbate biosynthesis [58] and was found enriched for Bob Gordon at a treatment of 112 kg N/ha. APX was increased in the same genotype. Both BAS1 and CAT participate in the reduction and degradation of H2O2 into water [59,60]. There was an enhanced abundance of BAS1 and CAT in Wyldewood at a treatment of 112 kg N/ha. In summary, these results suggest a possible enhanced defense in the plant against oxidative stress at higher N-conditions.

One notable feature of this study was the identification of multiple protein isoforms associated with N-treatment [61]. About 62% of 101 non-redundant proteins were recognized as isoforms. These isoforms show different Mr and/or pI changes on the 2-DE gels in response to N-treatment. There are several reasons for this phenomena, including genetic variation, alternative splicing of RNA transcripts, and PTMs (ubiquitination, acetylation, glycosylation, and/or phosphorylation) [62]. Many of the protein isoforms reported in the present study would be omitted using the label free mass spectrometry method (i.e. geLC-MS), which provides overall changes in protein abundance [63]. Different protein isoforms have been found to play an important role in mediating N. In maize, three GS2 isoforms were identified in the root and the accumulation trend of phosphorylated isoform is consistent with plastidial GS activity. Moreover, nitrate and ammonium have been found to induce the accumulation of different isoforms [64]. Similarly, the physiological role of a new isoform of GS has been reported in the leaf of wheat, suggesting a complex and flexible regulation for GS isoforms underlying the wheat that is associated with N utilization and plant growth [65]. In the present study, two isoforms of GS were consistently found to have opposite expression levels for Adams II and Bob Gordon among the three N-treatments, while the same down-regulated changes were observed for Wyldewood. The modification of multiple protein isoforms may give us novel perspectives to understand plant response to N-treatment.

As described above, the comparison among different levels of N in all genotypes reveled 23 unique proteins representing 58 isoforms. Among these, 18 unique proteins illustrating 23 isoforms showed a significant genotype × N level interaction (p < 0.05). These conserved proteins were then classified into protein orthologous groups and summarized in the cellular critical pathways. As protein modification may play an important role in this study, the detection of PTMs by three-dimensional gel electrophoresis may warrant further investigation [66]. However, we note the stringent conditions used in our study for the identification of differential abundant proteins by 2DE should help to diminish issues with co-migration of proteins.

It is well known that under N starvation, accumulation of lipid degradation and toxic compounds in plants can induce ROS [67]. To protect the cell from oxidative damage, ROS-scavenging pathways are activated [68]. Our proteomics results revealed that a number of enzymes involved in protein metabolism, redox homeostasis and transportation (indicated as red squares in Figure 4) were up-regulated for the three genotypes at a treatment of 169 kg N/ha, suggesting that oxidative stress could be more efficiently reduced at high N-concentrations. As part of the plant defense system, detoxification of non-enzymatic components (i.e. GSH) is also activated to destroy ROS [69]. LGL, which is found up-regulated with N-treatment in our study, is involved in the GSH-based detoxification of methylglyoxal [70].

Changes in protein expression profiles in photosynthesis, the OPP pathway and glycolysis suggest that plants may require a certain amount of N to enable these proteins to produce energy for plant growth in its most efficient way. Similar results were found in maize leaf [37]. GS is a key enzyme in ammonium assimilation [71] and an overall decrease in abundance of GS among the three genotypes in response to increasing levels of N lower its utilization. SNAlm plays a role in carbohydrate recognition [30]. The abundance of SNAlm decreased at the lower treatment level and increased at the higher levels. Other proteins associated with minor CHO metabolism and development were down-regulated at highest fertilizer levels.

5. Conclusions

This is the first time that the critical response pathways in American elderberry leaves have been revealed at the proteome level as a function of fertilizer treatment. Determining the protein expression profiles of the three genotypes in response to N-treatment, suggests that Bob Gordon is a distinct genotype from the other two. The use of three commercial genotypes of elder allows us to find similar functional distributions underlying the N-treatment and to explore the conserved expression patterns of critical proteins. Moreover, identification of proteins specifically expressed in photosynthesis, protein metabolism, and redox homeostasis across three genotypes, indicates these pathways played major role in response to different N-treatments. The identification of proteins and critical pathways provides a novel perspective into the use of N-fertilizers.