We used B. napus Express 617, a dark-seeded winter cultivar (Lee et al., 2020). B. napus was first grown in the greenhouse under long day conditions and then transferred outside for natural vernalization, followed by additional growth outside. A. thaliana Columbia 0 (Col-0, NASC ID N1092) and Nössen-0 (Nö-0, NASC ID N3081) were used as wildtype controls. The f3h mutant (tt6-2, GK-292E08, NASC ID N2105575, Col-0 background) (Appelhagen et al., 2014) and the ans/fls1 double mutant (synonym ldox/fls1-2, ldox: SALK_028793, NASC ID N2105579, Col-0 background; fls1-2: RIKEN_PST16145, Nö-0 background) (Stracke et al., 2009) were used for the generation of transgenic lines. A. thaliana plants were grown in the greenhouse under a 16-h-light/8-h-dark cycle at 22°C before transformation.

BnaFLS homologs were identified with KIPEs v0.255 as described previously (Pucker et al., 2020). KIPEs was run with a minimal BLAST hit similarity of 40% to reduce the number of fragmented peptides derived from possible mis-annotations. As bait, peptide sequences from the sequence collection of functional F3H, FLS, and ANS sequences described in KIPEs were used. As subject species, the peptide sequence sets of several Brassica species were used (Supplementary Table 1). The alignment was constructed with MAFFT v.7 (Katoh and Standley, 2013) and trimmed to minimal alignment column occupancy of 10%. Next, a phylogenetic tree was built with FastTree v2.1.10 (Price et al., 2009) using 10,000 rounds of bootstrapping, including the bait sequences and 2-ODD-like sequences from A. thaliana derived from Kawai et al., 2014 (Supplementary File 1). The phylogenetic tree was visualized with FigTree v1.4.3¹ (Supplementary Figure 1). Classification of BnaFLS candidates was generated based on the corresponding A. thaliana orthologs.

A comprehensive summary about gene-specific features of BnaFLS candidates is summarized in Supplementary Table 2. GSDS 2.0 (Hu et al., 2015) was used to generate gene structure plots. Literature knowledge was used to identify MYB-recognition elements (MRE) within 1 kbp upstream of the translational start site of BnaFLS candidates (Supplementary Figure 2). The conserved MRE consensus sequence 5′-AcCTACCa-3′, identified as a SG7 recognition motif (Hartmann et al., 2005; Stracke et al., 2007) and the sequence motifs important for the binding of AthMYB21 (MYBPZM: 5′-CCWACC-3′) and AthMYB24 (MYBCORE: 5′-CNGTTR-3′) to AthFLS1 were used for screening (Battat et al., 2019; Shan et al., 2020).

Theoretical isoelectric points, as well as molecular weight values of the BnaFLS protein sequences were calculated with ExPASY V (Gasteiger et al., 2005; Supplementary Table 3). In addition, SignalP v. 5.0 (Almagro Armenteros et al., 2019b) and TargetP v. 2.0 (Almagro Armenteros et al., 2019a) were used to infer the presence of signal peptides and N-terminal presequences of BnaFLS candidates, respectively (Supplementary Tables 4, 5). TMHMM v. 2.0 (Krogh et al., 2001) was used to predict transmembrane regions within BnaFLS sequences (Supplementary Table 2). Finally, Plant-mPLoc v. 2.0 (Chou and Shen, 2010) was used to predict the subcellular localization of BnaFLS candidates (Supplementary Table 2). Amino acid sequence identities of BnaFLSs compared to FLS homologs of A. thaliana, B. rapa, and B. oleracea were calculated based on a MAFFT alignment (Supplementary Table 6).² Protein sequence alignments were visualized at http://espript.ibcp.fr/ESPript/ESPript/index.php v. 3.0 (Robert and Gouet, 2014) using the AthFLS1 pdb file derived from Pucker et al., 2020. Functionally relevant amino acid residues and motifs for FLS and F3H activity were highlighted.

In silico secondary structure models of relevant BnaFLS candidates were generated via I-TASSER (Roy et al., 2010) and visualized with Chimera v. 1.13.1 (Pettersen et al., 2004). The AthF3H PDB file derived from Pucker et al., 2020 was used for visualization. The generated PDB files of this work can be accessed via Supplementary File 2.

Ribonucleic acid (RNA) samples were isolated from seeds and leaves using the NucleoSpin^® RNA Plant kit (Macherey-Nagel, Düren, Germany) according to manufacturer’s instructions. Seed samples of the B. napus cultivar Express 617 were collected 23 and 35 days after flowering (DAF), while leave samples were collected 35 DAF. Samples were collected in triplicates. The RNA quality was validated using NanoDrop and Agilent 2100 to confirm the purity, concentration, and integrity, respectively. Based on 1 μg of total RNA, sequencing libraries were constructed following the TruSeq v2 protocol. Three seed and leaf samples per genotype were processed. Single end sequencing of 82 nt was performed on an Illumina NextSeq 500 at the Sequencing Core Facility of the Center for Biotechnology (CeBiTec) at Bielefeld University.

Read quality was assessed by FastQC (Andrews, 2018), revealing reads of good quality reaching a phred score of 35 or above. Next, reads were mapped to the Express 617 reference genome sequence (Lee et al., 2020) using STAR v. 2.7.1a (Dobin et al., 2013). STAR was run in basic mode allowing maximal 5% mismatches per read length and using a minimum of 90% matches per read length. These read mappings were used to manually correct the functional annotation of the BnaFLS candidates (Supplementary File 3). The corresponding corrected annotation file was used for downstream analysis.

Beside the newly generated RNA-Seq data, publicly available RNA-Seq data sets were used and retrieved from the Sequence Read Archive³ via fastq-dump v. 2.9.6⁴ to analyze the expression of the candidate genes across various organs (Supplementary Table 7). Kallisto v. 0.44 (Bray et al., 2016) was used with default parameters to quantify transcripts abundance. The heatmap was constructed with a customized python script (see text footnote 2) using mean transcripts per millions (TPMs) per organ. Condition-independent co-expression analysis was performed to identify co-expressed genes using Spearman’s correlation coefficient (see text footnote 2) by incorporating 696 RNA-Seq data sets (Supplementary Table 8). To filter for strong co-expression the Spearman’s correlation coefficient threshold was set to 0.7 as suggested by Usadel et al. (2009).

Genes were functionally annotated by transferring the A. thaliana Araport11 (Cheng et al., 2017) functional annotation to the B. napus Express 617 gene models. The annotation was used for the co-expression analysis. OrthoFinder v. 2.3.7 (Emms and Kelly, 2019) was applied using default parameters to identify orthogroups between the representative peptide sequences of Araport11 and the B. napus Express 617 peptide sequences as previously defined (Pucker et al., 2017). Remaining non-annotated genes were functionally annotated by using reciprocal best blast hits (RBHs) and best blast hits (BBHs) as described previously (Pucker et al., 2016; Supplementary Table 9).

All constructs generated in this work were produced via Gateway cloning technique according to manufacturer’s instructions and verified by DNA sequencing (Supplementary Table 10). Total RNA from leaves and seeds of Express 617 was extracted as described above (see section “Gene Expression Analysis: Ribonucleic Acid Extraction, Library Construction, and Sequencing”). Complementary DNA (cDNA) was synthesized with the ProtoScript^TM Reverse Transcriptase kit (Invitrogen, Karlsruhe, Germany) using ∼1 μg of total RNA and 1 μl of oligo (dT) and 1 μl of random-hexamer primers. cDNA fragments corresponding to the full-length ORFs of the candidate genes were then amplified via PCR with Q5^® High-Fidelity Polymerase PCR kit (NEB, Frankfurt am Main, Germany) using gene-specific gateway primers (Supplementary Table 9). The sizes of the amplification products were analyzed by gel electrophoresis and visualized by ethidium bromide on a 1% agarose gel. The amplicons were purified from the PCR reagent tube via the NucleoSpin^® Gel and PCR Clean-up Kit (Macherey-Nagel, Düren, Germany).

The purified cDNA fragments corresponding to the full-length ORFs of the candidate genes were then recombined into pDONR^TM /Zeo (Invitrogen, Karlsruhe, Germany) using the Gateway BP Clonase II Enzyme Mix (Invitrogen, Karlsruhe, Germany) and the attB recombination sites of the respective gateway primers (Supplementary Table 10). Each entry clone was then used to transfer the CDS into the destination vector pLEELA (Jakoby et al., 2004) or pDEST17 (Invitrogen) via the Gateway LR Clonase II Enzyme Mix (Invitrogen, Karlsruhe, Germany). In pLEELA, the rapeseed coding sequences are under control of a double 35S promoter. pDEST17 was used for heterologous protein expression during the in vivo E. coli bioconversion assay under the control of the T7 promotor. The following constructs were available from previous studies: pDEST17-AthF3H, pDEST17-AthFLS1 (Busche et al., 2021), pDONR-AthFLS3, pDONR-AthFLS5, pDONR-AthANS (Preuss et al., 2009). The respective BnaFLS CDS sequences are listed in Supplementary File 4.

The bioconversion assay in E. coli subsequent HPTLC analysis of the methanolic extracts were performed as described in Busche et al. (2021). Successful heterologous expression of the recombinant proteins via SDS-PAGE was shown (Supplementary Figure 3).

The generated pLEELA-BnaFLSX constructs were used to transform the A. thaliana f3h knock out mutant, as well as the ans/fls1 double mutant using the A. tumefaciens strain GV3101:pM90RK (Koncz and Schell, 1986) according to the floral dip protocol (Clough and Bent, 1998). Selection of T1 plants was carried out by BASTA selection. Surviving plants were genotyped for the respective wildtype and mutant alleles, as well as the insertion of the transgene into the genome and its expression via PCR and RT-PCR (Supplementary Table 10). The genotyping for the presence of the transgene was repeated with T2 plants. T2 plants were used for the generation of flavonol-containing methanolic extracts as described below (see section “Flavonol Content Analysis by High-Performance Thin-Layer Chromatography”). T2 plants of the transformed ans/fls1 mutants and T3 plants of the transformed f3h mutants were used for DPBA-staining of young seedlings (see section “In situ Flavonoid Staining of Whole Seedlings”).

The flavonol glycosides were extracted and analyzed as previously described (Stracke et al., 2009). A. thaliana stems were homogenized in 80% methanol and incubated for 15 min at 70°C and then centrifuged for 10 min at 16,100 × g. The supernatants were vacuum-dried at 60°C and sediments were dissolved in 1 μl of 80% methanol mg^–1 starting material for HPTLC analysis. In total, 3 μl of each sample were spotted on silica-60 HPTLC-plates. The A. thaliana accessions Col-0 and Nössen-0, as well as the ans/fls1 double mutant were used as controls for the ans/fls1 A. thaliana complementation lines. For the f3h complementation lines, Col-0 and the f3h mutant were used as controls. The mobile phase consisted of a mixture of 66.7% ethyl acetate, 8% formic acid, 8% acetic acid, and 17.3% water. Flavonoid compounds were detected as described before (Stracke et al., 2009).

The visualization of flavonoids via DPBA-staining with whole seedlings was performed as described (Stracke et al., 2007), with the following minor adaptations: the bleached seedlings were stained to saturation in a freshly prepared aqueous solution of 0.25% (w/v) DPBA, 0.01% (v/v) Triton X-100, and 20% ethanol (v/v).

We identified a monophyletic group of 13 BnaFLS candidates through phylogenetic analysis using F3H, ANS, and 2-ODD-like protein sequences as outgroup to classify members of the 2-ODD family (Figure 2, Supplementary Figure 1, and Supplementary Table 2) of B. napus. The BnaFLS candidates were further classified within the FLS gene family based on their phylogenetic relationship to their most likely A. thaliana orthologs (Figure 2). Thereby, we identified two BnaFLS1, two BnaFLS2, five BnaFLS3, and four BnaFLS4 candidates in the B. napus cultivar Express 617. BnaFLS1-1 was identified on chromosome C09, while its homeolog BnaFLS1-2 is located on chromosome A09.

The genomic structure of the BnaFLS candidate genes comprises 3-4 exons and the encoded proteins display a length range from 270 to 336 amino acids (aa) (Table 1, Supplementary Figure 4, and Supplementary Table 2). Considering the chromosomal rearrangements as described for the cultivar Darmor-bzh (Chalhoub et al., 2014), homeologs were identified (Table 1).

No FLS5 and FLS6 homologs were identified in B. rapa, B. oleracea, and B. napus (Figure 2). Additionally screened B. napus cultivars (Gangan, No2127, Quinta, Shenglii, Tapidor, Westar, ZS11, Zheyou7) were in line with these results. As a FLS6 homolog is present in Raphanus sativus, a very close relative to B. rapa, B. oleracea and B. napus, the latter three might have lost FLS6 very recently. FLS5 was not found in the analyzed species of Brassiceae, Arabideae, Eutremeae, and Coluteocarpeae, while at least one copy was present in Camelineae and Boechereae indicating that FLS5 might have recently emerged in the latter tribes.

The expression of all BnaFLS candidate genes was analyzed by newly generated and publicly available RNA-Seq data (Table 2 and Supplementary Table 7). As seeds are the major organ for agronomical relevance, we screened for BnaFLS candidates expressed in seeds. In total, five genes were found to be expressed in seeds: BnaFLS1-1, BnaFLS1-2, BnaFLS2-1, BnaFLS3-3, and BnaFLS3-4. These five BnaFLS candidate genes revealed organ- and seed developmental-specific expression patterns (Table 2).

Both BnaFLS1 candidates revealed similar expression patters, showing the highest expression in late anther development, petals, and seeds. The expression of both BnaFLS1s tend to increase in siliques from 10 to 40 days after flowering (DAF). A similar expression pattern was observed in the seed coat revealing a development dependent expression. The biggest differences in BnaFLS1-1 and BnaFLS1-2 expression were observed in the embryo, where BnaFLS1-2 is higher expressed compared to BnaFLS1-1 indicating organ-specific transcriptional regulation at least for this organ. In contrast to BnaFLS1s, both BnaFLS3s are only marginally expressed in anthers and petals. While the expression of BnaFLS1s peaks during late seed and silique development, the expression of both BnaFLS3s peak in the early developmental stages. BnaFLS3-4 is highly expressed during seed coat development. Contrasting expression patterns of BnaFLS3-3 and BnaFLS3-4 were identified in e.g., seed coat samples indicating again organ-specific transcriptional regulation. BnaFLS2-1 was only marginally expressed in all analyzed organs, showing the highest expression in seed coat and roots. In summary, these findings indicate a role of BnaFLS1-1, BnaFLS1-2, BnaFLS2-1, BnaFLS3-3, and BnaFLS3-4 in seeds.

The five BnaFLS candidates expressed in seeds were used for downstream in-depth sequence- and functional analysis of the encoded proteins. The candidates revealed similar genomic structures and an alternative splice variant of BnaFLS2-1 was detected (Figure 3, Supplementary Table 2, and Supplementary Figure 5).

To get first insights into which biological pathways the five BnaFLS candidates expressed in seeds might be involved, we identified co-expressed genes (Supplementary Tables 11–15). Interestingly, the genes with the most similar expression pattern to BnaFLS1-1 are part of the flavonoid biosynthesis or the general phenylpropanoid pathway, including 4CL, CHS, CHI, F3H, F3′H, FLS1-2, UGT84A2, GSTF12, and MYB111. Similar results were obtained for BnaFLS1-2, which is co-expressed with homolog(s) of 4CL, CHS, CHI, F3H, FLS1-1, UGT84A2, and MYB111. Both BnaFLS1 genes contain the conserved subgroup 7 MYB-recognition element (MRE) motif in their putative promotor sequences (Supplementary Figure 2).

BnaFLS3-4 was identified to be co-expressed with genes which mostly lack a functional annotation. However, BnaFLS3-4 is strongly co-expressed with a MYB61 homolog. AthMYB61 is a known regulator of seed coat development. For BnaFLS2-1 (Spearman’s correlation coefficient < 0.59) and BnaFLS3-3 (Spearman’s correlation coefficient < 0.69) no genes with strong co-expression could be identified. This is likely due to the very weak expression of BnaFLS2-1 and the broad expression pattern of BnaFLS3-3 (Table 2).

To shed light on the potential functionalities of the BnaFLS candidates, the encoded proteins were compared to the well-characterized 2-ODD-members FLS, F3H, and ANS from A. thaliana (Table 3). BnaFLS1-1 and BnaFLS1-2 share > 91% sequence identity to AthFLS1, while BnaFLS2-1 has 57.4% sequence identity to AthFLS2. BnaFLS3-2 and BnaFLS3-3 revealed a sequence identity of 66.8% to AthFLS3. When comparing all BnaFLS candidates to AthF3H and AthANS, the protein identity ranged from 26.7-31 to 33.6-38.4%, respectively. The two BnaFLS1 candidates share 98.2% sequence identity, differing in 6 aa positions, while both BnaFLS3 candidates have 97.4% sequence identity, differing in 8 aa positions. The high sequence similarity between the BnaFLS candidates and their respective AthFLS orthologs implies close structural relationships and related functions.

The five BnaFLS candidates expressed in seeds were analyzed with respect to conserved amino acids and motifs important for FLS functionality (Figure 4). Both BnaFLS1 candidates contain all conserved amino acids and motifs. All remaining candidates lack the motifs potentially important for FLS activity, namely “SxxTxLVP”-, “CPQ/RPxLAL”-, and the N-terminal “PxxxIRxxxEQP,” in parts or completely. However, all BnaFLS candidates possess the conserved residues for ferrous iron- and 2-oxoglutarate-binding. Only BnaFLS2-1 revealed three amino acid exchanges in the five substrate binding residues analyzed, which are H103N, K173R, and E266D. BnaFLS3-3 and BnaFLS3-4 carry a G235A (G261 in AthFLS1) amino acid exchange. As some FLSs are bifunctional showing F3H-side activity, BnaFLS candidates were additionally screened for residues important for F3H activity (Figure 4). Besides the previously described G235A exchange of both BnaFLS3 candidates, all five BnaFLS candidates possess the residues described to play a role for F3H activity. The high conservation of relevant motifs and amino acids suggested both FLS1 candidates to be bifunctional. Due to the incomplete motifs and exchanges in conserved amino acids of BnaFLS3-3, BnaFLS3-4, and BnaFLS2-1 the FLS and/or F3H activity of these candidates might be affected.

Moreover, all BnaFLS candidates were predicted to contain no transmembrane helices, signal peptides or N-terminal presequences (mitochondrial-, chloroplast-, thylakoid luminal transfer peptide) and are therefore assumed and predicted to be located in the cytoplasm (Supplementary Tables 2, 4, 5).

For the functional characterization of BnaFLS1-1, BnaFLS1-2, BnaFLS2-1, BnaFLS3-3, and BnaFLS3-4 in vivo bioconversion assays in E. coli as well as analysis of stablely transformed A. thaliana knock out mutants were performed. The reproducibility of the bioconversion assay was ensured by showing that the observed functionalities of the well-known 2-ODD members AthFLS1, AthFLS3, AthFLS5, AthF3H, and AthANS match literature-based knowledge (Supplementary Figure 6). As expected, AthF3H showed clear F3H activity. In line with previous reports, AthFLS1 was identified as bifunctional possessing FLS activity and F3H side activity and AthANS showed FLS and F3H side activity. None of these activities could be detected for AthFLS5. Although AthFLS3 was reported to have FLS activity under extended assay conditions in E. coli, we could not detect FLS or F3H activity.

To investigate whether the bifunctionality of both BnaFLS1s compared to both BnaFLS3s, which showed only F3H activity, might be based on structural differences in silico, 3D models were generated (Figures 7A–F). The BnaFLS1s showed three major differences compared to both BnaFLS3s, which offer insights into the potential mechanisms of bifunctionality: (i) Both BnaFLS3 models revealed a shorter N-terminus compared to BnaFLS1s, resulting in the loss of the presumably FLS-specific “PxxxIRxxxEQP”-motif and α-helices (Figures 4, 7). (ii) The amino acid G261 proposed to be involved in proper folding is only present in both BnaFLS1s, while BnaFLS3s carry an alanine at this position. This residue is located between the transition of a beta-sheet from the jellyroll core structure to an α-helix. The hydrophobic side chain of alanine likely reduces the space in the catalytic center. (iii) Both BnaFLS3s show only partial overlaps with the “SxxTxLVP”- and “CPQ/RPxLAL”-FLS-specific sequence motifs (Figure 4). However, these mismatches do not have a substantial effect on the overall secondary structure in these regions (Figures 7E,F). Moreover, an extended N-terminus is not essential for F3H activity since it is absent in BnaFLS3-3 and BnaFLS3-4 (Figures 7D–F).

Although flavonols are of agronomical, ornamental, nutritional, and health importance, the major players of the flavonol biosynthesis in the oil and protein crop B. napus have not been investigated in great detail yet. So far, only one FLS gene was identified via transient expression in tobacco (Vu et al., 2015). However, as there are several members of the BnaFLS gene family expressed in seeds, it is necessary to characterize the encoding enzymes to infer which genes contribute to flavonol biosynthesis in B. napus seeds.

The members of the BnaFLS gene family are more closely related to each other than to any of the other 2-ODDs, which is in line with the results for the AthFLS gene family (Owens et al., 2008). In contrast to the AthFLS gene family, which is located on chromosome 5 in close proximity (Owens et al., 2008), the members of the BnaFLS gene family are distributed across seven chromosomes. Considering the chromosomal rearrangements described for B. napus cultivar Darmor-bzh (Chalhoub et al., 2014) and also the chromosomal positions of the B. rapa (Guo et al., 2014) and B. oleracea (Parkin et al., 2014) FLS genes, high local synteny of the FLS loci to those of B. napus Express 617 was identified. This syntenic relation allowed the assignment of 6 homeologous pairs of the B. napus FLS gene family. The homeolog pair BnaFLS3-3 and BnaFLS3-4 is located on the pseudochromosome C02 and clusters together with one additional unassigned BnaFLS3-5 homolog. The position of BnaFLS3-4 and BnaFLS3-5 on C02 in the Express 617 assembly likely derives from a mis-assembly as inferred by manual curation of the locus and the frequent assignment of the respective homologs to A02 in other long-read B. napus cultivar assemblies like Westar and Shengli (Song et al., 2020). Moreover, the respective B. rapa homologs BraFLS3-4 (Bra029212) and BraFLS3-5 (Bra029211) are located on A02. BraFLS3-4 and BraFLS3-5 are assumed to originate from duplication of the syntenic AthFLS2 to AthFLS5 tandem array. This duplication is part of the whole genome duplication (WGD), but only some genes of the array were retained in B. rapa (Guo et al., 2014). BraFLS3-5 is assumed to be derived from a gene duplication event of BraFLS3-4 (Guo et al., 2014), which is underlined by the close proximity of the two homologs BnaFLS3-4 and BnaFLS3-5 (only 2.9 kbp apart in the Express 617 assembly, see Table 1). Moreover, BraFLS2-2 (Bra038647) and BraFLS3-2 (Bra038648) are assumed to have emerged by WGD events as described before (Guo et al., 2014). A similar originating mechanism is assumed for BolFLS2-1 (Bo3g103270) and BolFLS3-1 (Bo3g103260) and thus for their respective homologs BnaFLS3-1 and BnaFLS2-1. Therefore, these ancient duplication events shaped the B. napus FLS gene family.

In A. thaliana FLS5 encodes a full-length protein, which contains amino acid exchanges important for hydrogen bonding of the substrate most likely resulting in a non-functional polypeptide (Owens et al., 2008; Preuss et al., 2009). In line with our results, no FLS5 homolog was identified in B. rapa (Guo et al., 2014). However, we could also not detect FLS5 in Brassiceae, Arabideae, Eutremeae, and Coluteocarpeae, but FLS5 was detected in the Camelineae, which include A. thaliana, as well as in Boechereae. Thus, we postulate that FLS5 emerged after the divergence of the common ancestor of the parental species of B. napus (B. rapa and B. oleracea) and A. thaliana rather than that it was frequently lost after the WGD events of the tandem array as postulated by Guo et al., 2014.

FLS6 was characterized as a pseudogene in A. thaliana (Owens et al., 2008; Stracke et al., 2009) and no FLS6 homolog was identified in B. rapa (Guo et al., 2014). These findings are in line with our results showing that FLS6 was lost very recently in B. rapa and B. oleracea and consequently is not present in B. napus, since FLS6 is still present in Raphanus sativus. As FLS6 was identified as pseudogene and FLS5 is known to encode a non-functional protein in A. thaliana (Owens et al., 2008; Preuss et al., 2009; Stracke et al., 2009), the parental species B. oleracea and B. rapa have already eliminated these “unnecessary” genes.

However, some BnaFLS genes are retained as they still encode functional proteins like BnaFLS3-3 and BnaFLS3-4, which encode for proteins with F3H activity. Importantly, both BnaFLS3s show a higher sequence identity with functional FLSs compared to F3H homologs, although exhibiting only F3H activity. This fact provides clear evidence that a classification solely based on amino acid sequences is not sufficient to infer functionalities of FLS family members and very likely 2-ODDs in general.

Bifunctionality has so far not been reported for a FLS from B. napus. By using two independent methods, we demonstrated bifunctionality of the two BnaFLS1 homeologs, which exhibit F3H and FLS activity. Thus, BnaFLS1-1 and BnaFLS1-2 are responsible for flavonol production in planta. We hypothesize that the respective orthologs of B. oleracea (Bo9g174290) and B. rapa (Bra009358) are bifunctional enzymes as well (Supplementary Table 6). Moreover, two additional members of the BnaFLS gene family have been functionally characterized. Interestingly, BnaFLS3-3 and BnaFLS3-4 revealed only F3H activity, while no FLS activity was detected. By incorporating sequence and structural analyses of 3D secondary structure models of BnaFLS1s vs. BnaFLS3s, we proposed a set of evolutionary events underlying the mechanisms of bifunctionality. Both BnaFLS3s lack several amino acids at the beginning of the N-terminus, which could cause the loss of FLS activity as it harbors the “PxxxIRxxxEQP” motif. This motif was previously proposed to be important for FLS activity as it distinguishes FLS from 2-oxoglutarate-/FeII-dependent dioxygenases with other substrate specificities (Owens et al., 2008; Stracke et al., 2009). Additional support for the relevance of this N-terminal region is provided by an AthFLS1 protein lacking the first 21 amino acids which showed no FLS activity (Pelletier et al., 1999; Owens et al., 2008). Moreover, amino acid exchanges in the “CPQ/RPxLAL”- and “SxxTxLVP”-motif in both BnaFLS3s possibly impact FLS activity. In addition, both BnaFLS3s carry a G235A (G261 in AthFLS1) amino acid exchange in comparison to BnaFLS1s, which might be relevant for bifunctionality as this exchange reduced the activity of a mutated Citrus unshiu FLS by 90% (Wellmann et al., 2002). This glycine is conserved across 2-ODDs and is suggested to play a role in proper folding (Wellmann et al., 2002). In accordance, we identified the A235 of BnaFLS3s and G261 of BnaFLS1 located between the transition of a beta-sheet from the jellyroll core structure to an α-helix, thereby the hydrophobic side chain of the alanine might reduce the space in the catalytic center. We propose that FLS bifunctionality is likely influenced by a combination of the identified motifs and residues rather than a single causative change as observed before for other flavonoid enzymes (Gebhardt et al., 2007; Seitz et al., 2007). The impact of each motif or amino acid on FLS bifunctionality needs further investigations that go beyond this study. As these sequence differences of BnaFLS3s do not abolish F3H activity, we uncovered that a truncated N-terminus and G261 are not essential for F3H activity. This is of importance as G261 was reported to be important for F3H activity (Britsch et al., 1993) while it may only play a minor role in conservation of F3H activity.

In addition to the FLS activity of BnaFLS1s, the 2-ODD member ANS might be able to contribute to flavonol production, as AthANS exhibit FLS and F3H side activities in vitro (Turnbull et al., 2004). In planta, FLS is the major enzyme in flavonol production as AthANS was not able to fully substitute AthFLS1 in vivo which is visible in the flavonol deficient fls1-2 mutant (Owens et al., 2008; Stracke et al., 2009).

BnaFLS2-1 is most likely a pseudogene. Although BnaFLS2-1 is still marginally expressed as shown by RNA-Seq data, it carries amino acid exchanges within 3/5 substrate binding residues in addition to a truncated N-terminus, which render the protein non-functional. In addition, an alternative transcript of BnaFLS2-1 was discovered (Supplementary Figure 5) that leads to a frameshift and thus likely encodes a non-functional protein as well. In A. thaliana, a heterologous expressed mutated FLS carrying one of the identified amino acid exchanges, namely K202R (K173R in BnaFLS2-1) is described to possess only 12% of the wild type FLS activity (Chua et al., 2008). In accordance, AthFLS2 encodes a most likely non-functional protein, which also harbors a truncated N-terminus (Owens et al., 2008). We assume that BnaFLS2-1 might be derived from a gene duplication event, losing its original function over time due to a pseudogenization process similar to that proposed for the AthFLS gene family members (Preuss et al., 2009; Stracke et al., 2009). The rather low expression of BnaFLS2-1 across various organs supports this hypothesis.

The spatio-temporal patterns of flavonol accumulation in B. napus are characterized by the activity of multiple BnaFLS genes. Both BnaFLS3s are expressed in early seed development while BnaFLS1s are expressed during late seed development (Table 2). The similar expression patterns of both BnaFLS1s are expected because they are homeologs. Thus, their expression patterns in the parental species B. rapa and B. oleraceae were likely to be very similar as they fulfill similar functions. In line with these results, BnaFLS1-1 and BnaFLS1-2 share co-expressed genes of the flavonoid and phenylpropanoid pathway. Both BnaFLS1s are co-expressed with MYB111, a regulator of flavonol biosynthesis (Stracke et al., 2007) and contain SG7 MRE in their putative promoter regions. Additionally, genes important for flavonoid transport into the vacuole and anthocyanidin/flavonol glycosylation like GSTF12 (TT19) and UGT84A2 (Kitamura et al., 2004; Yonekura-Sakakibara et al., 2012) were identified to be co-expressed with BnaFLS1s. These results further support the role of BnaFLS1-1 and BnaFLS1-2 as major players of flavonol biosynthesis in B. napus seeds. Moreover, transcriptomic and functional analysis of BnaFLS1s indicate gene redundancy.

Both BnaFLS1s were mainly expressed in reproductive organs as observed for AthFLS1 (Owens et al., 2008). BnaFLS3-4 was identified to be co-expressed with the well-known transcription factors MYB61, MYB123, and MYB5 which play a role in flavonoid biosynthesis and seed coat development in A. thaliana (Penfield et al., 2001; Li et al., 2009; Xu et al., 2014). This indicates a likely conserved transcriptional regulation between these two closely related species and supports the importance of flavonols during reproductive processes, e.g., pollen tube growth (Muhlemann et al., 2018).

In line with metabolomic studies showing that phenolic and flavonoid seed content maximized 35 days after flowering (DAF) (Wang et al., 2018), the expression of BnaFLS1s was higher at 35 DAF compared to 23 DAF. In accordance, most kaempferol and quercetin derivates reach their abundance peak at 35 DAF (Wang et al., 2018). Thus the expression pattern of BnaFLS1s fit well with the flavonol accumulation pattern of developing seeds, where flavonols contribute to seed quality (Wang et al., 2018).

Finally, the expression of BnaFLSs family members is not restricted to seeds. Some BnaFLSs were identified to be expressed in roots including BnaFLS3-3 and BnaFLS3-4 indicating a role of those BnaFLS family members in flavonoid biosynthesis in roots.

Engineering and breeding of flavonol content is of agronomical, economical, and ornamental importance (Takahashi et al., 2007; Cook et al., 2013; Yin et al., 2019). Flavonol content in petals influences pollinator attraction and drives microevolution of pollinators (Sheehan et al., 2015; Grotewold, 2016). Moreover, flavonols possess ROS scavenging activities and provide protection against UV-B radiation (Harborne and Williams, 2000). Besides the potential of engineering flavonol biosynthesis, anthocyanin and proanthocyanindin production can be engineered as FLS and DFR compete for substrates (Figure 8), thereby influencing important agronomical traits e.g., seed color (Luo et al., 2016). This study identified two bifunctional BnaFLS1s which are highly expressed in seeds and can thus be harnessed to engineer the metabolic flux of seed flavonol biosynthesis in the future (Figure 8). For example, the main bitter off-taste component in rapeseed protein isolates is kaempferol 3-O-(2″′-O-Sinapoyl-β-sophoroside) (Hald et al., 2019). Thus, the results of this study provide the basis for breeding low-phenolics lines with focus on the reduction of e.g., kaempferols in seeds, thereby supporting the use of rapeseed protein in human consumption.