Three biological replicates of buds, leaves, stems, roots (non-nodule containing segment adjacent to root nodules) and nodules from 3 year-old Clone 639 red alder saplings, that were greenhouse grown at Pacific Northwest National Laboratory (PNNL), were harvested at 3 time-points; early April (Spring), mid-June (Summer), and late October (Fall). Pollen and catkins were obtained from mature field-grown Clone 639 in Summer. Tissues were individually flash frozen immediately after harvest in liquid nitrogen and ground into fine frozen powders using a Freezer/Mill (SPEX, Metuchen, NJ) with a program that contained 2 cycles, with a 1 min pre-cool, a 2 min run time at 10 cycles/min and a 1 min cooldown period. Total RNA was extracted from each tissue type using the Spectrum™ Plant Total RNA kit (Sigma Aldrich) according to the manufacturer’s suggested protocol. Extracted RNA samples were stored at –80°C until ready for use.

Total RNA integrity was individually assessed using Fragment Analyzer (Advanced Analytical Technologies, Ankeny, IA) with the High Sensitivity RNA Analysis Kit. RNA Quality Numbers (RQNs) from 1 to 10 were assigned to each sample to indicate its integrity or quality. “10” stands for a RNA sample without any degradation, whereas “1” marks a completely degraded sample. RNA samples with RQNs ranging from 8 to 10 were only used for RNA library preparations with the TruSeq Stranded total RNA Library Prep Kit with Ribo-Zero rRNA Removal Kit (Illumina, San Diego, CA). Sizes of RNA libraries were assessed by Fragment Analyzer with the High Sensitivity NGS Fragment Analysis Kit, and their concentrations were measured by StepOnePlus Real-Time PCR System (ThermoFisher Scientific, San Jose, CA) with the KAPA Library Quantification Kit (Kapabiosystems, Wilmington, MA). Libraries were diluted to 2 nM with Tris-HCl buffer (10 mM, pH 8.5) and denatured with 0.1 N NaOH. Eighteen pM libraries were clustered in a high-output flow cell using HiSeq Cluster Kit v4 on a cBot (Illumina). After cluster generation, the flow cell was loaded onto HiSeq 2500 for sequencing using HiSeq SBS kit v4 (Illumina). DNA was sequenced from both ends (paired end) with a read length of 100 bp at the WSU Spokane Genomics Service Center. The raw bcl files were converted to fastq files using software program bcl2fastq2.17.1.14. Adaptors were trimmed from the fastq files during the conversion.

For each sample, RNA-seq reads were assembled into putative transcripts using the Trinity pipeline (Grabherr et al., 2011). Read coverage was normalized using bbnorm (https://sourceforge.net/projects/bbmap/), with normalized paired end reads then aligned to the red alder reference genome assembly using the default parameters of HISAT2 (Kim et al., 2019). The resulting SAM file was converted to BAM format and sorted, then the individual transcriptomes were assembled with Trinity (Grabherr et al., 2011) using the reference assembly (Hixson et al., 2023) as a guide, with parameters: –max_memory 50G –no_normalize_reads –genome_guided_bam –genome_guided_max_intron 10000. Assembly statistics were extracted with the TrinityStats.pl script provided as part of the Trinity package.

The bioinformatics component of this study built upon the red alder Clone 639 sequencing and assembly (Hixson et al., 2023). This study herein focused upon genes involved (or potentially involved) in generation of the red alder plant defense, plant dye, and medicinal natural products, i.e. shikimate-chorismate-phenylalanine-β-glucogallin, and phenylpropanoid/phenylpropanoid-acetate pathway derived metabolites, including diarylheptanoids, ellagitannins (hydrolysable tannins), proanthocyanidins (condensed tannins), and other plant phenols. Sequences of previously characterized proteins/enzymes were retrieved from NCBI and used as queries to search the Clone 639 predicted proteins using BLASTP.

Protein sequences for dehydroquinate dehydratase/shikimate dehydrogenase (DHQD-SDH), UDP glucosyltransferase (UGT), arogenate dehydratase (ADT), cinnamate 4-hydroxylase (C4H), p-coumarate 3-hydroxylase (C3H), ferulate 5-hydroxylase (F5H), flavonoid 3′-hydroxylase (F3′H), flavonoid 3′,5′-hydroxylase (F3′,5′H) and flavone synthase II (FNSII) homologs were also retrieved from databases of seven Fagales species: Black alder (Alnus glutinosa), silver birch (Betula pendula), grey oak (Casuarina glauca), European hazel (Corylus avellana), European beech (Fagus sylvatica), American black walnut (Juglans nigra) and pedunculate oak (Quercus robur) (see Supplementary Table S1 for database source).

Unrooted phylogenetic trees were generated using Clustal Omega (Madeira et al., 2024) and rendered using iTol (Letunic and Bork, 2021).

Gene models of red alder Clone 639 genes of interest were curated in the Apollo annotation system (Dunn et al., 2019), assisted by RNA-seq data, to correct errors and determine the most likely gene models. GFF files of the Apollo gene models were exported and Transcripts Per Million (TPM) were calculated using TPMCalculator (Vera Alvarez et al., 2019). TPM values were plotted using the GD::Graph perl module to arrive at an expression atlas for each gene.

Aliquots of each of the same cryogenically pulverized samples (each done in triplicate) harvested for RNA-seq analysis from roots, root nodules, stems, leaves, buds, pollen and catkins were individually used for metabolomics analyses. For metabolite extraction, each previously pulverized sample (as described above) was initially freeze-dried. Each freeze-dried sample (5 – 10 mg) was then extracted with MeOH–H₂O (70:30, v/v) containing 0.1 mM naringenin (Aldrich) internal standard (IS) at a ratio of 1 mL to 50 mg tissue. After adding extraction solvent, each solution was vortexed for 30 s, sonicated in cold H₂O for 15 min, vortexed again and centrifuged (20,000 g × 15 min, at 4°C). Each supernatant was transferred individually to vials and subjected to ultra-high performance liquid chromatography – quadrupole time-of-flight – mass spectrometry (UPLC-QToF-MS) analyses.

This employed an ACQUITY™ UPLC system (Waters, Milford, MA, USA) equipped with a photodiode array (PDA) eλ detector (Waters) coupled to a Xevo™ G2 QToF mass spectrometer (Waters MS Technologies, Manchester, UK) using MassLynx (V4.1) software. Separations of metabolites used an BEH C₁₈ column (Waters, 2.1 × 150 mm, 1.7 µm particle size), with a linear gradient for separation: 100% A (0.1% HCO₂H in H₂O) over 0.5 min, then sequentially to 45% B (0.1% HCO₂H in CH₃CN) over 25.0 min, to 100% B over 1.0 min, to 100% A in 2 min and held at 100% A for an additional 11.5 min. Flow rate was 0.2 mL min^–1, with column and sample temperatures kept at 25 and 10°C, respectively. Injection volumes were 2 μL for LC–MS analysis, with UV-visible spectra recorded between 200 and 500 nm (1.2 nm resolution). An electrospray ionization (ESI) source was used to detect masses of eluted compounds (m/z range: 100 – 1000 Da) and Ar was the collision gas. Detection settings were as follows: Negative ion mode (capillary voltage at 2.0 kV; cone voltage at 30 eV; collision energy at 6 eV and at 27 eV) and positive ion mode (capillary voltage at 3.0 kV; cone voltage at 30 eV; collision energy at 6 eV and at 18 eV). Sodium formate (5 mM in 2-propanol-H₂O, 90:10, v/v) was used for calibrating the mass spectrometer, with leucine enkephalin (2 ng µL^–1 in A:B, 50:50, v/v) employed as lock-mass. Metabolomics data processing was carried out as in Lu et al. (2017).

All metabolite peak information (i.e. retention time, UV spectrum, and MS fragmentation analysis) from LC–MS data was used for targeted and non-targeted metabolite annotation. For targeted metabolite identification and annotation of compounds in red alder tissue extracts, authentic standards were used to establish coincident retention times, comparable molecular ions, and corresponding MS fragmentation patterns. For some metabolites (i.e. p-hydroxycinnamoyl quinic acid isomers), the previously reported order of their elution and fragmentation patterns allowed for the basis of their identifications (Gutiérrez Ortiz et al., 2018). Identification of metabolic features were performed using RStudio (Smith et al., 2006) as in Lu et al. (2017). Metabolomics data (i.e. meta-table of red alder tissue extract UPLC/EIMS analyses) included 673 features with 69 identified and/or putatively annotated (see Supplementary Data Set and Supplementary Table S2 ), with the total ion chromatograms (TIC) of the various red alder Clone 639 tissues individually used to help identify metabolites. The 69 metabolites above belonged to 8 distinct plant metabolite classes including diarylheptanoids, ellagitannins and precursors, proanthocyanidins, flavonoids, organic acids, gallic acid derivatives, phenolic acids, and polyamines, these being either fully identified or provisionally annotated.

Authentic alder diarylheptanoid standards included metabolites, oregonin, rubranol, hirsutanonol, rubranoside A, rubranoside B, alnuside A, alnuside B, platyphylloside, and 2′′′-O-p-coumaroyloregonin, 5-hydroxy-1,7-bis(4-hydroxyphenyl)heptan-3yl-β-D-apiofuranosyl-(1→6)-β-D-glucopyranoside, 1,7-bis-(4-hydroxyphenyl)heptan-3yl-β-D-apiofuranosyl-(1→6)-β-D-glucopyranoside, 7-(3,4-dihydroxyphenyl)-1-(4-hydroxyphenyl)-3(R)-β-D-glucosyloxy-heptane, and 5-O-methylhirsutanonol (Lai et al., 2012; Novaković et al., 2013), with each of these being used to analyze MS fragmentation patterns of all diarylheptanoids present in red alder Clone 639. Characteristic ion fragments of diarylheptanoids (m/z values of 327.12, 331.16, 313.14 311.13, 295.13, 205.09, and 189.09) were also used to obtain extracted ion content (EIC) from red alder tissues for additional diarylheptanoid identification/tentative annotation.

Other authentic standards used for comparison purposes included various flavonoids, organic acids, phenolic acids, spermidine derivatives, ellagitannins and precursors thereof, and proanthocyanidins:

Calculated and experimentally determined masses, together with ppm differences of each feature were used to assess annotation accuracy (see ppm error values listed in Supplementary Table S2 ). Relative abundances (normalized to internal standard naringenin equivalents) of annotated and tentatively annotated features from metabolic classes are shown using Excel. The features were classified as to occurrence(s) in nodules, roots, stems, buds, leaves, pollen, and catkins, respectively, with normalized data as obtained above.

Supplementary Figure S1 depicts a simplified biochemical shikimate-chorismate pathway to phenylalanine (Phe) via shikimate, as well as to quinic acid and β-glucogallin via pathway offshoots from 3-dehydroquinate. The main carbon flux overall to these different metabolites is well known to be directed towards the phenylpropanoid pathway as this primarily results in structural lignin deposition (Nature’s second most abundant terrestrial vascular plant biopolymers), as well as to various lower molecular weight phenylpropanoids (particularly diarylheptanoids in red alder) and related mixed biochemical pathway intermediates. Also metabolically important are the quinic acid and β-glucogallin metabolic pathway branches, albeit representing relatively minor carbon flux offshoots to various quinate derivatives and to gallo- and ellagitannins, respectively.

These different metabolic branches, however, all utilize the 3-dehydroquinate intermediate, with the latter being differentially converted into shikimic, quinic, and gallic acids, respectively ( Supplementary Figure S1 ). As described below, while only some of the overall possible branch point pathway intermediates are detected in the metabolomics analyses of different plant tissues at different time-points, it was considered instructive nevertheless to compare/contrast these trends with transcript abundances at the different time points. [Note: while these analyses were taken at three distinct time-points in the growing season or were harvested (catkins and pollen) in the summer for mature field grown Clone 639, such analyses give no indication of either any potential metabolic turnover or transcriptional abundance changes and/or protein turnover].

A BLAST search of our red alder predicted proteins using previously known bona fide shikimate-chorismate-phenylalanine and related pathway genes from Arabidopsis and other species gave 27 putative genes ( Supplementary Table S3 ; Figure 2A ); these are shown placed in the order of their biochemical pathway transformations.

Their corresponding 27 proteins had plastid targeted signal peptides as determined using DeepLoc 2.0 (Thumuluri et al., 2022) and LOCALIZER 1.0.4 (Sperschneider et al., 2017), with the exceptions of the presumed cytosolic bifunctional 3-dehydroquinate dehydratase/shikimate dehydrogenase (DHQD-SDH) sub-family, a putative chorismate mutase (ArCM2), and the putative β-glucogallin-forming UDP glucosyltransferases (UGTs) ( Supplementary Table S3 ). The existence of a cytosolic chorismate mutase (CM2) has long been known (d’Amato et al., 1984), with its gene widely conserved across plant species (Westfall et al., 2014). However, its physiological function was until recently unknown. In petunia (Petunia hybrida), Qian et al. (2019) reported that its cytosolic PhCM2 was able to form prephenate (from chorismate), which can next be converted first to phenylpyruvate by phenylpyruvate dehydratase (PDT), and then to Phe by action of a cytosolic phenylpyruvate aminotransferase (PPY-AT) (discussed later).

Amino acid sequence alignments for each enzymatic step leading to Phe, relative to established functional proteins in other plant species, are in Supplementary Figures S2 – S11 . Provisional red alder amino acid identity prioritization was based on comparison to known bona fide proteins ( Supplementary Table S3 ). Some of the putative red alder proteins though had relatively low identities (ca 50% or lower) with respect to corresponding known homologs from Solanum lycopersicum, Corydalis sempervirens and Arabidopsis, i.e. including members of the ArDHQD-SDH, shikimate kinase (ArSK), ArCM2, and arogenate dehydratase (ArADT2 and ArADT3) families.

In red alder, dehydroquinate synthase (ArDHQS), ArSK, 5-enolpyruvylshikimate-3-phosphate synthase (ArEPSPS), chorismate synthase (ArCS) and prephenate aminotransferase (ArPPA-AT) are apparently single gene families. The others are encoded by small (3 – 6) multi-gene families, i.e. 3-deoxy-D-arabinoheptulosonate 7-phosphate synthase (ArDAHPS, 5 sub-family genes), ArDHQD-SDH (6), ArCM (3) and ArADT (3), respectively.

Relative transcript abundances for genes putatively encoding each enzymatic step in each tissue at different time-points are in Figure 2A and Supplementary Figure S12 . Highest expression levels were noted across essentially all tissues and time-points for ArDAHPS3 and ArDAHPS4 (highest expression levels overall), ArDHQS, ArDHQD-SDH1, ArEPSPS, ArCS, ArCM1, ArPPA-AT and, to a lesser extent, ArADT3, respectively. These specific transcript abundances however presumably represent genes encoding each of the actual biochemical steps involved in the predominant phenylpropanoid carbon flux to lignin, this being evident for nearly all tissues and time-points, with the exception of the lower overall levels in pollen. Conversely, ArSK transcript levels were notably different (low throughout) for all tissue types and time-points.

Interestingly, about half of the 27 potential shikimate-chorismate pathway genes are reasonably strongly expressed in the 7 distinct tissue types and time-points and are thus likely to encode the main sub-family isoforms involved overall.

As regards red alder’s other ArDAHPS isoforms, it was noted that both ArDAHPS1 and ArDAHPS2 marginally had the highest amino acid identity to bona fide DAHPS from Arabidopsis. Yet, their encoding genes only had moderate to low expression levels in catkins, these being even lower in all other tissues ( Supplementary Figure S12A ). Furthermore, ArDAHPS5 had low gene expression across all tissue types (highest in bud, leaf, and stem tissues), the significance of which is currently unknown.

Of the 5 single gene family members, four (ArDHQS, ArEPSPS, ArCS, and ArPPA-AT) were expressed in all tissues at different time-points, with several being highest in Fall collected tissue ( Supplementary Figure S12B ). Somewhat comparable expression patterns for these 4 single gene family members were also noted for catkins, but not for pollen. Conversely, ArSK was the only single gene family member that consistently displayed low expression levels, i.e. weakly in nodules, roots, stems, buds, while barely detectable in leaf, catkin and pollen tissues.

Shikimate is formed from 3-dehydroquinate via 3-dehydroshikimate. In bacteria and fungi, these conversions are catalyzed by two monofunctional enzymes, dehydroquinate dehydratase (DHQD) and shikimate dehydrogenase (SDH), whereas in plants a bifunctional enzyme, DHQD-SDH, is used ( Supplementary Figure S1 ).

Interestingly, A. thaliana has only one DHQD-SDH (AT3G06350, AtDHQD-SDH), whereas other plants have multiple gene copies such as red alder here. The AtDHQD-SDH reportedly shows a strict requirement for shikimate as substrate [in vitro assays were carried in the reverse direction, as no activity could be measured with quinate as substrate even at very high concentrations (Gritsunov et al., 2018)]. DHQD-SDH genes, whose corresponding proteins had quinate activity, have also been characterized in Populus trichocarpa (PoptrDHQD-SDH2 and 3) (Guo et al., 2014), whereas DHQD-SDHs with the capacity of producing gallic acid were identified in walnut (Juglans regia) and Vitis vinifera (Muir et al., 2011; Bontpart et al., 2016).

Phylogenetic analysis of DHQD-SDHs from different species within the dicotyledons showed that they cluster into five groups (I–V) (Bontpart et al., 2016). We thus constructed an unrooted phylogenetic tree using the sequences from Bontpart et al. (2016), the six red alder proteins (ArDHQD-SDH1 – ArDHQD-SDH6), as well as DHQD-SDHs from other Fagales genomes available [A. glutinosa, B. pendula, C. glauca, C. avellana, F. sylvatica, J. nigra, and Q. robur ( Supplementary Table S1 ; Figure 2B )]. The 6-red alder DHQD-SDHs encompass all five sub-families, as do all Fagales DHQD-SDHs.

The bifunctional ArDHQD-SDH1 in Group I presumably converts 3-dehydroquinate into shikimate, with this envisaged to occur via an initial dehydration catalyzed by the DHQD domain to afford 3-dehydroshikimate, followed by a NADPH dependent reaction employing the SDH domain to afford shikimate. By contrast, ArDHQD-SDH2 in Group II can be envisaged to catalyze the NADH dependent conversion of 3-dehydroquinate into quinate. The Group III ArDHQD-SDH3 and Group IV ArDHQD-SDH4 likely produce gallic acid via initial formation of 3-dehydroshikimate (via the DHQD domain as above) followed by NADP⁺-dependent dehydrogenation using the SDH domain to produce gallic acid.

Group V homologs (ArDHQD-SDH5 and ArDHQD-SDH6) are, however, presently of unknown biochemical function.

Subsequent structural, bioinformatics, and biochemical approaches identified a primary sequence motif (³⁷⁹SX[TG]³⁸¹; AtDHQD-SDH numbering) that could be used to predict substrate specificity between shikimate- and quinate-forming DHQD-SDHs (Gritsunov et al., 2018). That is, the Arabidopsis AtDHQD-SDH wild-type (WT) protein, harboring a Thr amino acid residue at position 381, only uses shikimate. However, when this residue was mutated to Gly, the mutant protein was able to oxidize quinate to 3-dehydroquinate (in the reverse reaction) (Gritsunov et al., 2018) and all known DHQD-SDHs in Group II have a Gly residue at this relative position ( Supplementary Figure S4 ), i.e. perhaps suggesting that ArDHQD-SDH2 (and its Fagales homologs) produce quinate from 3-dehydroquinate. Interestingly, all group II DHQD-SDHs have a conserved Ser-Val-Gly (SVG) motif, as do Group V enzymes (although the latter’s biochemical function is as yet unknown). By contrast, Group III enzymes have a Ser-Cys-Thr (SCT) motif at these positions, whereas Group IV has a corresponding Ser-Tyr-Thr (SYT) motif, both of these being considered to be involved in gallate formation.

A cofactor key binding residue variation was also observed in all Group II sequences, e.g. the aspartate-isoleucine-aspartate (⁴⁸³DID⁴⁸⁵) motif is usually associated with NAD⁺ binding, in place of asparagine-arginine-threonine (⁴⁸³NRT⁴⁸⁵) which uses NADP⁺ in most other DHQD-SDHs (including in A. thaliana).

Relative ArDHQD-SDH sub-family gene expression (heatmap) and transcript abundances are in Figure 2A and Supplementary Figure S12C . As indicated above, the most highly expressed overall was the Group I ArDHQD-SDH1, particularly in Fall collected aerial tissues, versus a lower extent in root and nodules tissues ( Figure 2A ; Supplementary Figure S12C ), i.e. in accordance with its presumed role in producing shikimic acid for phenylpropanoid/phenylpropanoid-acetate pathways. The putative quinate Group II forming ArDHQD-SDH2 was second most highly expressed during Summer in stems, buds, and catkins, as well as to a much lower extent in nodules, roots, leaves, and pollen. By contrast, the presumed gallic acid-forming Group III ArDHQD-SDH3 and Group IV ArDHQD-SDH4 genes were only mildly expressed in root and bud tissues and, to an even lesser extent, in all other tissue types, except for pollen where they were not detected (discussed later below). The Group V ArDHQD-SDH5 and ArDHQD-SDH6, of unknown biochemical functions, had faint expression levels in stem, bud, and leaf tissues collected in Summer. Their expression patterns though may give some clues as to their functions, i.e. in terms of seeking correlations in the future with unique stem, bud, and leaf metabolites perhaps only found in these tissues.

ArCM (3 genes) and ArADT (3 genes) expression was also examined ( Figure 2A ; Supplementary Figures S12D, E ). Of the 3 ArCM genes, ArCM1 was most highly expressed, particularly in Fall collected tissues and catkins, except for pollen ( Figure 2A ; Supplementary Figure S12D ). ArCM2 was also relatively highly expressed, but only in pollen and catkin tissues. By contrast, ArCM3 was barely – if at all – detectable in bud and pollen tissues.

As regards the 3 red alder ADTs, ArADT3 was generally the most highly expressed followed by ArADT2 and ArADT1 ( Figure 2A ; Supplementary Figure S12E ). To better understand their possible physiological roles, we constructed an unrooted phylogenetic tree using ADT sequences from gymnosperms and angiosperms including the 7 Fagales species mentioned above ( Supplementary Figure S13 ). ArADT2 clusters in the group with AtADT2 in ADT subgroup II (Cho et al., 2007), this enzyme being considered involved in protein synthesis in vascular plants, as knockouts of the gene were lethal (Corea et al., 2012). By contrast, ArADT3 is in the same cluster as AtADT3, AtADT4, AtADT5 and AtADT6 in subgroup III. Of the latter, double knockouts of AtADT4 and AtADT5 gave substantial reductions in lignin levels and related metabolites in Arabidopsis, with the corresponding quadruple knockout (adt3/4/5/6) having even larger effects (Corea et al., 2012; Höhner et al., 2018). ArADT1 (in subgroup I) which was expressed to the lowest level has not yet been directly correlated with a downstream metabolic process.

Interestingly, in petunia, an alternative transcription start site in PhADT3 results in a cytosolic PhADT3 which can use prephenate generated by the cytosolic PhCM2 (as indicated earlier) to form phenylpyruvate. The latter can then be converted to small amounts of Phe through action of a cytosolic phenylpyruvate aminotransferase (PPY-AT) (Qian et al., 2019).

Best potential candidate genes for gallic acid metabolism to β-glucogallin were ArUGT1 and ArUGT2, based on sequence comparison with bona fide β-glucogallin-forming UDP glucosyltransferases (UGTs) in other plant species (Mittasch et al., 2014; Cui et al., 2016; Ono et al., 2016). ArUGT1 was relatively mildly expressed in all tissue types (except for pollen), whereas ArUGT2 was essentially not detected, except for a very faint level in catkins and nodules (Spring/Summer) ( Figure 2A ; Supplementary Figure S12F ). Neighbor-Joining analysis also showed that both ArUGT1 and ArUGT2 cluster within the phylogenetic group L of UGTs (Ross et al., 2001), together with ester-forming glucosyltransferases and more particularly with all characterized β-glucogallin-forming UGTs ( Supplementary Figure S14 ). The Fagales homologs also cluster within this group. (Other red alder genes in group L, ArUGT3, ArUGT4 and ArUGT5, are discussed later below).

It was next instructive to ascertain which metabolites in the biochemical pathways leading to shikimic, quinic, and gallic acids were detectable in the seven tissue types, i.e. in an attempt to correlate gene expression with their metabolite occurrence in the different tissues. As previously indicated, some 673 features were detected of which 69 distinct metabolites were identified and/or putatively annotated. Numerous metabolites were identified based on either comparison with authentic standards, or via analyses of elution profile behavior and mass spectrometric fragmentation profiles relative to literature reports (discussed later). Provisional annotations do not though represent full unambiguous identification. Peak intensities in all chromatograms were normalized to the internal standard naringenin.

In this study, only quinic, shikimic, and gallic (as its β-glucogallin derivative) acids ( Figure 2C ) were unambiguously detected in all 7 tissue types to varying amounts ( Figure 2D ). Quinate was present to the largest amount in leaves, and also detected in all seven tissue types.

Mixed metabolic pathway derivatives were also detected with quinate esters being the most abundant.

In terms of correlations with transcript levels, shikimate (derived from ArDHQD-SDH1) and quinate (derived from ArDHQD-SDH2) metabolite levels were generally low in roots, nodules, and pollen ( Figures 2A, D ). Transcript level profiles for ArDHQD-SDH1 and ArDHQD-SDH2 showed comparable profiles. By contrast, β-glucogallin was at low levels in essentially all tissues examined (discussed later).

Both metabolomics and transcriptome data were thus in excellent agreement as regards each tissue type gene expression and the corresponding shikimic/quinic/gallic acid derived pathway constituents present.

As indicated above, phenylpropanoids, largely derived from Phe, are mainly ultimately converted into abundant cell wall lignin biopolymers, as well as to various low molecular weight compounds ( Supplementary Figure S15 ). In red alder, the latter include phenylpropanoid-derived esters/glucosides ( Figure 2D ), amides ( Figure 3B ), and diarylheptanoids ( Figure 4 ; Supplementary Figure S16 ), with the diarylheptanoids vastly dominating in relative amounts.

A BLAST search of the red alder predicted proteins yielded 36 putative genes encoding all phenylpropanoid pathway transformations from Phe to lignin/lignan precursor monolignols, p-coumaroyl, coniferyl and sinapyl alcohols (Figure 3C; Supplementary Figures S17 – S26 ; Supplementary Table S4 ), with gene families ranging from 1 – 10 members. (The first three enzymatic steps in the phenylpropanoid pathway, namely phenylalanine ammonia lyase (PAL), cinnamate 4-hydroxylase (C4H), and 4-coumarate CoA ligase (4CL) are also shared with the flavonoid (phenylpropanoid-acetate) pathway as described later).

Red alder homologs were generally ranked in decreasing sequence identity to bona fide homologs in other plant species, with these mainly prioritized again in terms of highest sequence identity, i.e. ca 70% or higher to those of known function. However, some red alder phenylpropanoid pathway multi-gene family members had low homologies to bona fide enzymatic steps, i.e. those annotated as ArC4H2 and ArC4H3 (ca 62% identity), Ar4CL4 (60.3% identity), caffeic acid O-methyltransferase (ArCOMT2, 62.3% identity), caffeoyl CoA O-methyltransferase (ArCOMT3, 58% identity), cinnamoyl CoA reductase (ArCCR2, 50.3% identity), and 8 cinnamyl alcohol dehydrogenases (ArCAD3 to ArCAD10, with 50.1 – 47.5% identities), respectively. In some of these annotated proteins, such low amino acid level identities may reflect different catalytic functions as regards substrate(s).

An unrooted phylogenetic tree was constructed that included sequences of characterized hydroxylases in the phenylpropanoid pathway (i.e. C4H, p-coumarate 3-hydroxylase, C3H, and ferulate-5 hydroxylase, F5H families), these belonging to the CYP73A, CYP98A, and CYP84A families, respectively, including those from the seven Fagales species mentioned above ( Supplementary Figure S27 ). ArC4H1 clustered with CYP73A family members Helianthus tuberosus HtC4H (Teutsch et al., 1993), Arabidopsis AtC4H (Mizutani et al., 1997), and Medicago sativa MsC4H (Fahrendorf and Dixon, 1993), whereas ArC4H2 and ArC4H3 clustered with a French bean (Phaseolus vulgaris) PvC4H (Nedelkina et al., 1999) as a separate CYP73A subgroup. Additionally, ArF5H clustered with the CYP84A AtF5H (Meyer et al., 1996). Finally, the four alder C3Hs (ArC3H1 – 4) clustered with the CYP98A, AtC3H (Schoch et al., 2001), ObC3H from basil (Ocimum basilicum) (Gang et al., 2002) and CcC3H1/2 from coffee (Coffea canephora) (Mahesh et al., 2007).

Transcript abundances of the aforementioned 36 genes in different tissues and organs are provided in a cumulative (relative) heatmap ( Figure 3A ), and for each enzyme class individually in Supplementary Figure S28 . Within specific gene families, the isoforms generally with highest expression levels ( Figure 3A ; Supplementary Figures S28A–I ) in the seven tissue types included ArPAL2, ArC4H1, hydroxycinnamoyl CoA:shikimate hydroxycinnamoyl transferase ArHCT1, ArC3H1/2, Ar4CL4, ArF5H, ArCOMT1/2, cinnamoyl CoA reductase ArCCR1, as well as ArCAD1, ArCAD3, and ArCAD5. These mainly had high to highest sequence identities at the amino acid level relative to bona fide enzymatic steps in other plant species; however, Ar4CL4, ArCAD3, and ArCAD5 had lower identity levels. Taken together, these expression levels are indicative of an intact pathway to the monolignols with the corresponding expected isoforms.

Catkin tissue also had moderate to relatively high expression of ArPAL2/ArPAL1 ( Figure 3A ; Supplementary Figure S28A ), ArC4H1/ArC4H2/ArC4H3 and ArF5H ( Figure 3A ; Supplementary Figure S28B ), ArHCT1 ( Figure 3A ; Supplementary Figure S28C ), ArC3H1/ArC3H2 ( Figure 3A ; Supplementary Figure S28D ), Ar4CL1/Ar4CL2 ( Figure 3A ; Supplementary Figure S28E ), ArCOMT1 ( Figure 3A ; Supplementary Figure S28F ), ArCCOMT1/ArCCOMT2 ( Figure 3A ; Supplementary Figure S28G ), ArCCR1 ( Figure 3A ; Supplementary Figure S28H ) and ArCAD1/ArCAD3/ArCAD4 ( Figure 3A ; Supplementary Figure S28I ).

ArC3H2 also displayed notable expression levels in bud tissues, whereas ArC3H3 and ArC3H4 ( Supplementary Figure S28D ) had essentially low to no expression in leaf, pollen, and catkin tissues.

Of the CADs, ArCAD1 was most highly expressed in nodules, roots, and catkins, and to a lesser extent in all other tissues ( Supplementary Figure S28I ), this being a homolog most closely related to phenylpropanoid/monolignol/lignin formation proper (Kim et al., 2004; Jourdes et al., 2007).

ArCAD3, however, was highly expressed in nodules, roots, stems, buds, and leaves, particularly in the Fall collected tissues, as well as in catkins, whereas ArCAD5 had highest expression levels in lignifying stems and in buds. By contrast, ArCAD7 expression was highest in both nodule and root tissues, the physiological significance of the latter being currently unknown given its low sequence identity.

Three additional ArUGTs (ArUGT3, ArUGT4, and ArUGT5) were present in the red alder genome, albeit of lower homology (~70%) to the afore-mentioned β-glucogallin-forming UGTs ( Figure 2A ; Supplementary Table S3 ). These, and other Fagales UGTs, clustered between the β-glucogallin-forming UGTs and hydroxycinnamic acid 1-O-glucosyltransferases from Arabidopsis (UGT84A1–UGT84A4) and Brassica napus (BnUGT84A9) (Milkowski et al., 2000a, b), perhaps suggesting a role in downstream plant phenol monomer glucosylation ( Supplementary Figure S14 ). The putative ArUGT3 and ArUGT4 genes were generally moderately to weakly expressed in nodules, roots, and stems. ArUGT4 was also expressed in buds, and, to a lesser extent, in leaves and catkins ( Supplementary Figure S12F ), whereas ArUGT5 had lower levels of expression in stem, bud, leaf, pollen, and catkin tissues, but was not detected in nodules or roots. ArUGT3 – ArUGT5 expression was not detected in pollen tissue.

Under the conditions employed, our metabolomics analysis did not detect all possible monomeric phenylpropanoid pathway products per se, but instead various derivatives thereof. This is because most of the carbon flux went directly into lignin biopolymer deposition (not shown) and diarylheptanoids ( Figures 4A, B ; Supplementary Figures S16A–D ), with the latter discussed later. Much less abundant low molecular weight phenylpropanoid derivatives detected included the aforementioned mixed metabolic pathway hydroxycinnamoyl quinate ester derivatives, of which trans-4-O-p-coumaroyl quinic acid, cryptochlorogenic (4-O-caffeoyl quinic acid), and chlorogenic acid had highest (albeit variable) accumulation levels across stems, buds, leaves, and catkins and, to a lesser extent, in roots ( Figures 2C, D ; Supplementary Table S2 ). Pollen tissue, by contrast, lacked any of these esters under the conditions employed.

In addition, five other phenylpropanoid derivatives were found in the various tissues, with spermidine derivatives highest in pollen, i.e. tri-p-coumaroyl spermidine, N¹,N¹⁰ -bis(p-coumaroyl)spermidine (confirmed with authentic standard), and a N¹,N¹⁰ -bis(p-coumaroyl)spermidine like metabolite, as well as two phenylpropanoid diarylheptanoid esters (i.e. 2′′′-O-p-coumaroyl oregonin and cinnamoyl hirsutanol-5-O-β-glucopyranoside) in other tissues ( Figure 3B ; Supplementary Table S2 ).

Of these, the N¹,N¹⁰ -bis(p-coumaroyl)spermidine like metabolite differs from the authentic N¹,N¹⁰ -bis(p-coumaroyl)spermidine standard chromatographically. This could perhaps be due to either different cis/trans p-coumaroyl isomers or different p-coumaroyl positional isomers.

Hydroxycinnamic acid amide-linked spermidine metabolites had previously been reported in various Alnus species, and are reportedly specific to the Fagales in pollen exine (Meurer et al., 1988).

While the main metabolic flux is of course lignin and diarylheptanoid targeted, the genes that are significantly more highly expressed (ArPAL1/2, ArC4H1, ArHCT1, ArC3H1/2, Ar4CL4, ArCOMT1/2, ArCCR1, and ArCAD1; Figure 3A ) can be considered as the main bona fide phenylpropanoid pathway genes controlling flux into lignins, diarylheptanoids ( Figure 4 ; Supplementary Figure S16 ), phenylpropanoid esters ( Figure 2D ), and phenylpropanoid amides ( Figure 3B ).

Previously discussed expression levels ( Supplementary Figure S12C ) of ArDHQD-SDH2 (encoding the envisaged quinic acid forming isoform), and ArDHQD-SDH1 (for shikimate formation and leading mainly to lignin, hydroxycinnamic acid intermediates, and diarylheptanoids), were thus also in good agreement with metabolite occurrences ( Figures 2D , 3B , 4A, B ; Supplementary Figure S16A–D ; Supplementary Table S2 ). ArDHQD-SDH2 expression was highest in stem, bud, and catkin tissues reflecting the higher levels of hydroxycinnamic acid quinate esters present. ArDHQD-SDH1 was more highly expressed in all tissue types, as to be anticipated for shikimate formation for hydroxycinnamic acid, diarylheptanoid, and lignin biosynthesis. Pollen, by contrast, only had traces of quinic acid and diarylheptanoids.

A simplified flavonoid biochemical pathway, such as to catechin and epicatechin, is summarized in Supplementary Figure S29 , and involves 11 additional specific enzymatic steps/enzyme classes beyond those already described above for PAL, C4H, and 4CL ( Supplementary Figure S15 ).

Entry point substrates for flavonoid and PA formation are phenylpropanoid-derived p-coumaroyl CoA and acetate pathway derived malonyl CoA as shown ( Supplementary Figure S29 ). Metabolomics analyses of the various red alder tissues resulted in detection of 17 different flavonoids and proanthocyanidins ( Supplementary Table S2 ), of which apigenin, 7-O-methyl apigenin (genkwanin) [not shown], luteolin, isoquercitrin [not shown], quercitrin [not shown], (+)-catechin, (–)-epicatechin, and procyanidin B2 ( Figure 5B ) were unambiguously identified based on comparison with their authentic standards. Relative abundances of known and those of other provisionally annotated flavonoids (examples in Figure 5B ) are shown in Figures 5C–E for the different tissues and timepoints in the growing season. While several of these flavonoids were quercetin derivatives, together with catechin/epicatechin, the relative amounts and compositions differed substantially depending upon tissue type.

Leaf, pollen, and catkin tissues had the highest concentrations of flavonoids (in naringenin equivalents), relative to the other tissues ( Figure 5E ). Most abundant metabolites in leaves were quercetin derivatives, these mainly being quercitrin (quercetin-3-O-rhamnoside), quercetin glucuronide, and other quercetin derivatives. By contrast, dihexosylquercetin was present as the most abundant flavonoid in pollen, while this and other quercetin derivatives were in lower abundance in catkins.

Stem and bud tissues displayed somewhat different metabolic profiles ( Figure 5D ). Of their quercetin derivatives, both isoquercitrin (quercetin 3-O-glucoside) and quercetin glucuronide were the most abundant of this sub-class of metabolites. In those same tissues, catechin and epicatechin were either of highest overall abundance or near highest level of abundance with smaller levels of the related PAs, procyanidin B2 and a PA trimer, at all three time points in the growing season. Buds also accumulated both a putative viscidulin III like metabolite and dimethoxyluteolin at all time points. Interestingly, levels of apigenin were highest in the Fall collected bud tissue.

Root and nodule tissues ( Figure 5C ) displayed very distinct metabolomic profiles, both in relative amounts (these being overall lowest in the two tissue types) and in compositions. Their tissues mainly accumulated differing amounts of catechin/epicatechin and PAs with overall amounts being decreased in Fall tissue. Small amounts of quercetin glucuronide and quercitrin were also present with the latter at near constant levels for all time points.

A red alder genome BLAST search gave 33 putative genes in the flavonoid pathway (Figure 5A; Supplementary Tables S4, S5 ; Supplementary Figures S30 – S40 ). ArPAL2 (as well as ArPAL1 in catkins), ArC4H1, and Ar4CL4 had been discussed in the preceding phenylpropanoid pathway section.

Following these three enzymatic steps, the further downstream flavonoid pathway gene families ranged from one (chalcone isomerase ArCHI; flavanone-3 hydroxylase ArF3H; flavonoid 3′, 5′-hydroxylase ArF3′,5′H; and leucoanthocyanidin dioxygenase ArLDOX); two (flavonoid 3′-hydroxylase ArF3′H1 and ArF3′H2; flavanol synthases ArFLS1 and ArFLS2; dihydroflavanol 4-reductases ArDFR1 and ArDFR2; leucoanthocyanidin reductases ArLAR1 and ArLAR2; and anthocyanidin reductases ArANR1 and ArANR2); four (chalcone synthases ArCHS1, ArCHS2, ArCHS3, and ArCHS4, as well as flavone synthases II ArFNSII-I, ArFNSII-2, ArFNSII-3, and ArFNSII-4). The most likely candidate genes generally had sequence identities ca ≥70%, relative to bona fide genes in other plant species, with exception of putative FNSII genes which had much lower sequence identities (ca 48 – 57%).

The corresponding heatmap of differing expression levels of ArPAL, ArC4H, Ar4CL, ArCHS, ArCHI, ArF3H, ArF3′H, ArF3′,5′H, ArFLS, ArDFR, ArLAR, ArLDOX, ArANR and ArFNSII gene family members in the different tissues is shown in Figure 5A . [Gene expression patterns for different isoforms of the individual enzyme classes are in Supplementary Figure S41 , with ca 27 of the 33 possible gene candidates having readily detectable gene expression levels albeit sometimes in a tissue specific manner].

For the enzymatic steps shared with the phenylpropanoid pathway, it had already been discussed that ArPAL2, ArC4H1, and Ar4CL4 were the most highly expressed genes. For the flavonoid specific pathway, CHS genes ArCHS1, ArCHS2, and ArCHS3 (with amino acid sequence identities of ~ 86, 74, and 86%, respectively, to bona fide CHSs) were highly expressed in nodules, roots, and catkins collected in Fall, and in Summer/Fall for the catkins, buds and stems ( Figure 5A ; Supplementary Figure S41A ). ArCHS2 expression was also highly expressed in leaves in Fall, whereas ArCHS1 and ArCHS3 were both expressed to lower levels. In contrast, ArCHS4 was essentially not detected. Both ArCHS2 and ArCHS4 were, however, of somewhat lower sequence identity (ca 72–74%) to bona fide CHSs from other plant species ( Supplementary Table S5 ). Furthermore, although CHS homologs suggest gateway entry points to flavonoids, related CHS-like polyketide synthase (PKS) proteins have been implicated as a gateway entry point to diarylheptanoids (Katsuyama et al., 2009) as discussed later below. ArCHI gene expression, by contrast, was variable, but highest in nodules and roots (Fall) and in buds (Summer and Fall) ( Figure 5A ; Supplementary Figure S41B ).

The combined enzymatic activities of ArCHS and ArCHI result in naringenin formation ( Supplementary Figure S29 ), which can then be metabolized to differentially afford quercetin, luteolin, kaempferol, apigenin, catechin, and epicatechin skeleta.

Other genes with overall high to low expression levels in different tissue types were ArF3H, ArF3′,5′H ( Supplementary Figure S41B ), ArF3′H1, ArF3´H2 ( Supplementary Figure S41C ), ArDFR2, ArDFR1 ( Supplementary Figure S41E ), ArLAR2, ArLAR1 ( Supplementary Figure S41F ), and ArLDOX, ArANR1, ArANR2 ( Supplementary Figure S41G ). Leaf, pollen and catkin tissues generally had much lower expression levels of these genes, with the exception of ArFLS1 and ArLAR1 in pollen and catkin tissues ( Supplementary Figure S41D, F ).

From our metabolomics analyses, red alder tissues have at least 2 PAs, procyanidin B2 and a PA trimer ( Figures 5B–E ), 6 ellagitannins or “hydrolysable tannins” ( Figure 6C ), and 25 diarylheptanoids ( Figures 4A, B ; Supplementary Figures S16A–D ; Supplementary Table S2 ). While biochemical pathways for formation of substrates prior to entry into these three metabolic classes are well-established, the proteinaceous biochemical entry points to PAs, ellagitannins, and potentially diarylheptanoids, represent a major gap in our scientific knowledge, as do their downstream (post-gateway entry) metabolic processes. That is, genes encoding proteins affording entry points to both ellagitannins and PAs, as well as potentially to diarylheptanoids, are either unknown or poorly understood, in spite of their widespread and enormous chemical structural diversity across the plant kingdom.