Plant age-dependent dynamics of annatto pigment (bixin) biosynthesis in Bixa orellana

Abstract

Age affects the production of secondary metabolites, but how developmental cues regulate secondary metabolism remains poorly understood. The achiote tree (Bixa orellana L.) is a source of bixin, an apocarotenoid used in diverse industries worldwide. Understanding how age-dependent mechanisms control bixin biosynthesis is of great interest for plant biology and for economic reasons. Here we overexpressed miRNA156 (miR156) in B. orellana to comprehensively study the effects of the miR156–SQUAMOSA PROMOTER BINDING PROTEIN-LIKE (SPL) module on age-dependent bixin biosynthesis in leaves. Overexpression of miR156 in annatto plants (miR156ox) reduced BoSPL transcript levels, impacted leaf ontogeny, lessened bixin production, and increased abscisic acid levels. Modulation of expression of BoCCD4-4 and BoCCD1, key genes in carotenoid biosynthesis, was associated with diverting the carbon flux from bixin to abscisic acid in miR156ox leaves. Proteomic analyses revealed an overall low accumulation of most secondary metabolite-related enzymes in miR156ox leaves, suggesting that miR156-targeted BoSPLs may be required to activate several secondary metabolic pathways. Our findings suggest that the conserved BomiR156–BoSPL module is deployed to regulate leaf dynamics of bixin biosynthesis, and may create novel opportunities to fine-tune bixin output in B. orellana breeding programs.

Achiote tree, age, bixin, microRNA156, secondary metabolism, vegetative phase change

Introduction

The vegetative phase change pathway controls the transition between juvenile and adult phases in plants, affecting several fundamental metabolisms, including the production of secondary metabolites (Li et al., 2020). It is known that secondary metabolites vary at different stages of the plant life cycle, but how endogenous developmental cues regulate secondary metabolism remains poorly understood. Leaf age and developmental stage affect the production of secondary metabolites (Li et al., 2016; Vazquez-Leon et al., 2017). For instance, the biosynthesis of some monoterpenes and sesquiterpenoids begins as early as in the first cotyledon in Melaleuca alternifolia (Russell and Southwell, 2002). Synthesis and accumulation of essential oil of Cinnamomum verum initiates in young leaves (Li et al., 2016). Nevertheless, compounds associated with the sabinene hydrate-terpinen-4-ol-γ-terpinene pathways are produced at later stages of leaf development (Russell and Southwell, 2002). In terms of defense, total phenolics, alkaloids, and cyanogenic glycosides increase in some species when the plant ages (Elger et al., 2009). Although leaf ontogenetic changes are genetically programmed, they are not necessarily synchronized with whole-plant ontogenetic changes. New leaves are continuously produced throughout the vegetative phase, and long-lasting plants such as trees may represent a complex mosaic or a gradient of leaf tissues at different ontogenetic phases of development. These ontogenetic differences may impact the production of secondary metabolites in the plant as a whole. Localized on the surface or inside the leaves, secretory structures such as nectaries, resin ducts, secretory vesicles, salt glands, oil cells, and secretory trichomes are frequently among the main sites for synthesis and accumulation of secondary metabolites (Figueiredo et al., 2008). Therefore, distinct ontogenetic stages of leaves may affect the concentration of various primary and secondary plant products (Verma and Shukla, 2015). However, little is known about the molecular mechanisms underlying the age-dependent control of secondary metabolism in non-model plants, such as trees.

Popularly known as the achiote tree or annatto, Bixa orellana L. (Bixaceae) belongs to the order Malvales, whose representatives include cocoa, hibiscus, and cotton (The Angiosperm Phylogeny Group, 2003, 2009). Annatto is a source of the secondary metabolite bixin, an apocarotenoid commonly used in several industries worldwide (Ramamoorthy et al., 2010; Vilar et al., 2014). Bixin does not alter the nutritional value of foods, and has broad applications that vary from sunscreens to coloring processed food (Vilar et al., 2014). In leaves and seeds of B. orellana, apocarotenoids such as bixin are derived from the oxidative cleavage of carotenoids, through the action of specific enzymes previously identified, such as LYCOPENE CLEAVAGE DIOXYGENASE (BoLCD) and CAROTENOID DIOXYGENASES (BoCCDs) (Soares et al., 2011), ALDEHYDE DEHYDROGENASES (BoALDHs), and SABATH family methyltransferase (BoSABATH) (Cárdenas-Conejo et al., 2015). The bixin pathway may share common enzymes and substrates with other carotenoid-derived compounds, such as abscisic acid (ABA) (Cazzonelli and Pogson, 2010). The annatto tree produces another hydrophilic natural dye, norbixin, which is mostly used to provide color to non-oil-based products (Giuliano et al., 2003; Rivera-Madrid et al., 2013). Additional secondary metabolites such as terpenes, sesquiterpenes, and alkaloids are produced by annatto and have broad-spectrum antimicrobial properties (Karmakar et al., 2018; Santos et al., 2021). Understanding the age-dependent underlying mechanisms controlling the biosynthesis of bixin and other secondary metabolites in B. orellana is of great interest not only for plant biology, but also for the pharmaceutical, cosmetic, and textile industries (Valério et al., 2015).

The main network controlling the transition between juvenile and adult phases includes the highly conserved microRNA miR156 and its transcription factor targets, members of the SQUAMOSA PROMOTER-BINDING PROTEIN LIKE (SPL) gene family (Cardon et al., 1999; Wang et al., 2009; Wu et al., 2009; Yu et al., 2010; Morea et al., 2016). The miR156–SPL module has been extensively studied in several species, including Oryza sativa (Xie et al., 2012), Nicotiana tabacum (Feng et al., 2016), Pyrus pyrifolia (Qian et al., 2017), Solanum lycopersicum (Silva et al., 2014), Passiflora edulis (Silva et al., 2019), and Populus tremula (Lawrence et al., 2021). The miR156–SPL module is associated with leaf maturation, trichome formation, male fertility (anther development), shoot branching, phytomer length and composition, glucose and trehalose metabolism, floral organ development, and stress responses (Wang et al., 2009; Xing et al., 2010; Yu et al., 2010; Poethig, 2013; Fouracre and Poethig, 2016; Silva et al., 2019; Lawrence et al., 2021; Barrera-Rojas et al., 2023; Ferigolo et al., 2023). Interestingly, several miRNA regulatory modules, including the miR156–SPL module, have been shown to be important for controlling the biosynthesis of secondary metabolites (Padhan et al., 2016; Samad et al., 2020). The miR156 target, SPL9, directly binds to the promoter of the sesquiterpene synthase gene TPS21 to activate terpene synthesis (Yu et al., 2014). In Arabidopsis, miR156-targeted SPL9 negatively regulates anthocyanin accumulation by directly preventing expression of anthocyanin biosynthetic genes (Gou et al., 2011). SPLs may also directly or indirectly influence carotenoid formation in plants. The tomato colorless non-ripening (Cnr) locus encodes an Arabidopsis SPL3 homologue, and the Cnr mutation leads to low concentrations of total carotenoids, resulting in tomato fruits with a white pericarp (Fraser et al., 2001; Manning et al., 2006). More recently, banana (Musa acuminata) MuSPL16 was shown to be a positive regulator of carotenoid production through direct activation of carotenoid biosynthetic genes (Zhu et al., 2020), whereas papaya (Carica papaya) SPL1/SBP1 is a direct transcriptional repressor of PHYTOENE DESATURASE4 (PDS4), which encodes a key enzyme involved in the conversion of phytoene into lycopene (Han et al., 2019).

Although the control of the biosynthetic pathway of apocarotenoids in B. orellana is still poorly understood, recent evidence indicates that bixin biosynthesis is modulated by developmental cues. For instance, bixin levels increased during seed maturation (Moreira et al., 2023). We hypothesized that the miR156–SPL module regulates bixin production in annatto leaves in an age-dependent manner. To test this conjecture, we previously generated transgenic annatto trees overexpressing miR156 (Faria et al., 2022). Here, we showed that lower levels of bixin correlated with an increase in the expression of B. orellana miR156. Enhanced miR156 expression prolonged the juvenile phase and impacted leaf ontogeny, leading to reduced bixin levels in leaves as late as 90 d after planting. Moreover, the expression patterns and accumulation of the main enzymes involved in bixin and ABA biosynthetic pathways are modified in miR156-overexpressing annatto plants. Our findings suggested that the carbon flux in miR156-overexpressing B. orellana leaves was somewhat redirected from bixin to ABA production, indicating an age-dependent leaf dynamic of bixin biosynthesis.

Materials and methods

Plant material and growth conditions

Two miR156 overexpression (miR156ox) lines of B. orellana cv. Piave Vermelha were selected based on the miR156 expression levels. These lines were generated by introducing the AtMIR156a gene under the control of the CaMV35S promoter as described by Faria et al. (2022). Eight cloned plants of each miR156ox line and four non-transformed or control lines of the T₀ generation were acclimatized as described by Faria et al. (2019). Plants were transferred to 12 litre polyethylene pots with horticultural soil conditioner substrate (Tropstrato HT Hortaliças, Vida Verde Indústria e Comércio de Insumos Orgânicos Ltda, Mogi Mirim, SP, Brazil), and maintained under greenhouse conditions. The fertilizer Osmocote Plus 16-08-12 5-6M (ICL, USA) was added (2 g per pot) on a bimonthly basis. The plants were watered daily and maintained on an 11/13-h light/dark photoperiod at 25/16 °C (day/night). Access to the genetic heritage adhered strictly to current Brazilian biodiversity legislation and was approved by the Brazilian National System for the Management of Genetic Heritage and Associated Traditional Knowledge (SISGEN) under permission number AF3BC47. All genetically modified organism manipulations were performed in the Laboratório de Cultura de Tecidos II, Bioagro and in a genetically modified organism-adapted greenhouse following the CTNBio and CIBio-UFV biosafety rules.

Phenotyping and microscopy

miR156ox lines were maintained under greenhouse conditions and evaluated at 15, 30, 45, 60, 75, and 90 days after planting (DAP). The following variables were assessed: total leaf number, plant height (cm), number of phytomers, diameter of the main axis, length of the third leaf (from petiole insertion to the apical portion of the leaf blade), leaf area, and ramification index (total ramification length and main plant axis length) (Morris et al. 2001). We measured the area of all non-transformed leaves and then sampled the leaf area of miR156ox plants using the same number of leaves as in non-transformed plants. The total miR156ox leaf area was estimated as sampled leaf area×15, because miR156ox plants had 15 times more leaves than non-transformed plants. All images were analysed using ImageJ software version 1.43u (National Institutes of Health, Bethesda, MD, USA) (Schneider et al., 2012).

For scanning electron microscopy analysis, stems of non-transformed and miR156ox were cut longitudinally in the third phytomer rootwards, because the extrafloral nectaries were secreting in non-transformed plants. Then, the sectioned fragments were transferred to 0.1 M Karnovsky solution (Karnovsky, 1965) under a −250 mmHg vacuum for 1 h. The samples were dehydrated in a crescent ethanol series (50% to 100%) under −250 mmHg vacuum in each dehydration step. Following this, the samples were critical point dried with CO₂ (CPD 030, Balzers), and metallized with a thin layer of gold (20 nm) in a sputter-coater (Quorum Q150RS). Samples were examined in a scanning electron microscopy (LEO 1430 VP) at an acceleration voltage of 10.6 kV.

RNA extraction and cDNA synthesis

Expression profiles of bixin- and carotenoid-related genes, BomiR156, BomiR172, and BoSPLs in annatto were assessed by quantitative real-time PCR (qRT-PCR) using the CFX96 Touch Real-Time PCR Detection System (Bio-Rad Laboratories, Hercules, CA, USA). Approximately 10 mg of freeze-dried third fully expanded leaves were ground in liquid nitrogen and used for this analysis. The final reaction volume of 10 µl included 40 ng cDNA, 400 nM primers, SYBR-Green Supermix (Bio-Rad), and deionized water. Expression data were derived from four biological samples, with at least two technical replicates. Oligonucleotide sequences used in this study are listed in Supplementary Table S1. Expression data were normalized using the 40S RIBOSOMAL PROTEIN S9 (BoRPS9) (Moreira et al., 2018; Faria et al., 2022) (Supplementary Table S1), and expression levels of each gene and miRNA were calculated using the method (Livak and Schmittgen, 2001).

Identification, annotation, gene structure, and phylogenetic analysis

To identify annatto sequences containing SBP/SPL domains (PF03110), all 16 members of the Arabidopsis and 19 of the rice (Oryza sativa) SPL gene family (Cardon et al., 1999; Xie et al., 2006) were downloaded from TAIR (https://www.arabidopsis.org/) and Phytozome v13 (https://phytozome-next.jgi.doe.gov/) and used as a query to conduct local BLASTn against all 73,381 transcripts generated by RNA-Seq from Bixa orellana (Moreira et al., 2023). All retrieved coding sequences were translated in all six frames using EMBOSS Transeq (https://www.ebi.ac.uk/Tools/st/emboss_transeq/) and analysed by PFAM v35.0 (Salazar et al., 2020) to maintain only those containing the SBP/SPL domain (Supplementary Table S2). The coding region of each BoSPL gene was predicted through a Gene Structure Display Server (GSDS 2.0) program (http://gsds.gao-lab.org/, Hu et al., 2014).

All unique and primary coding sequences from tomato (S. lycopersicum, ITAG4.0 assembly), cacao (Theobroma cacao, v2.1 assembly), poplar (Populus trichocarpa, v4.1 assembly), grape (Vitis vinifera, v2.1 assembly), and switchgrass (Panicum virgatum, v4.1 assembly) were downloaded at Phytozome v13 and used to recovery sequences containing the PF03110 domain. The SBP/SPL gene names were used as described previously for tomato (Salinas et al., 2012), poplar (Hou et al., 2013; Guo et al., 2021), grape (Díaz-Riquelme et al., 2012; Hou et al., 2013), and switchgrass (Wu et al., 2016). All protein sequences for each multigene family were aligned using the default settings of Multiple Sequence Comparison by Log-Expectation (MUSCLE) (http://www.ebi.ac.uk/Tools/msa/muscle/; Madeira et al., 2022). A maximum likelihood (ML) phylogenetic analysis was performed using PhyML3.0 (Guindon et al., 2010), using Smart Model Selection (SMS) to select the best model (Lefort et al., 2017) (like VT+R+F substitution model) and approximate likelihood-ratio test branch support testing (Anisimova and Gascuel, 2006). The phylogenetic tree was visualized using iTOL software (https://itol.embl.de; Letunic and Bork, 2021).

To identify potential BomiR156 and BomiR172 mature sequences, microRNA databases available for Arabidopsis, tomato, cocoa, grape, poplar, and rice were used, such as miRBase (https://mirbase.org/, Kozomara and Griffiths-Jones, 2010) and PMRD: plant microRNA database (http://bioinformatics.cau.edu.cn/PMRD/; Zhang et al., 2010). The identification of miR156 target regions in silico and SPL gene sequences in Bixa orellana was performed through multiple alignments in MUSCLE and psRNATarget (https://www.zhaolab.org/psRNATarget/; Dai et al., 2018), using 20–21 nt miR156 sequences. Given that most sequences used in this study contain only the coding sequence, it is possible that the unidentified miR156 target regions may be localized in the 3ʹ-untranslated region of particular BoSPLs.

Bixin content and bixin channel quantification

Bixin content was quantified as described (Napoleão et al., 2017; Faria et al., 2019). Approximately 10 mg of freeze-dried third fully expanded leaves from 3-month-old plants was ground in liquid nitrogen, followed by repeated extraction with 300 μl methanol: isopropanol: acetic acid (20:79:1; v/v/v). Samples were agitated four times for 20 s, sonicated for 10 min, kept on ice for 30 min, and centrifuged at 20 000 g for 10 min at 4 °C. The supernatant was collected, filtered (Econofltr PVDF 13 mm, 0.2 μm; Agilent Technologies, Santa Clara, CA, USA), and analysed by liquid chromatography–tandem mass spectrometry (LC-MS/MS). Samples of the final supernatant (5 μl) were automatically injected into a triple quadrupole (QqQ) LC-MS instrument (6430; Agilent Technologies, Waldbronn, Germany) equipped with an Eclipse Plus C18 column (2.1 × 50 mm, 1.8 μm) and a Zorbax SB-C18 guard column (1.8 μm, Agilent Technologies) at 26 °C. The mobile phase consisted of 0.02% acetic acid in water (solvent A) and 0.02% acetic acid in acetonitrile (solvent B), at a constant flow rate of 300 μl min⁻¹. A linear gradient was applied as follows: 0–3 min, 2% to 60% B; 3–8 min, 60% to 99% B; 8–11 min, 99% B; 11–12 min, 99% to 2% B; and 12–15 min, 2% B. An electrospray ionization source was used under the following conditions: gas temperature of 300 °C, nitrogen flow rate of 10 l min⁻¹, nebulizer pressure of 35 psi, and capillary voltage of 4000 V. Bixin was analysed via multiple reaction monitoring of ion pairs with mass transitions (395/157) and quantified via calibration curves using pure standards (1–200 μg and 0.1–500 ng; Sigma-Aldrich, St Louis, MO, USA). Data were analysed using Mass Hunter Workstation software (Agilent Technologies).

Bixin channels are visible secretory structures that develop during organ maturation (Lima et al., 2013). For bixin channel counting, images of the abaxial leaf surfaces of 3-month-old non-transformed, miR156ox_1, and miR156ox_4 plants were captured using a stereomicroscope (SZX7; Olympus, Tokyo, Japan) coupled to a Moticam 580 5.0 MP camera (Nikon, Tokyo, Japan). Images were analysed using ImageJ version 1.49 (Schneider et al., 2012). Bixin secretory channels were quantified by sampling the basal, medial, and apical leaf portions, with three different areas per image, using three leaves per genotype. The number of bixin channels was divided per mm².

Abscisic acid profile

ABA levels were quantified in the third (from apex to base) fully expanded leaves collected from 3-month-old greenhouse-grown non-transformed and miR156ox plants. Leaf samples (approximately 110 mg) were powdered in liquid nitrogen. A 300 μl aliquot of extracting solution (methanol: isopropanol: acetic acid 20:79:1) was added to 1.5 ml microtubes. The samples were vortexed four times for 20 s each, sonicated for 5 min, kept at 4 °C for 30 min, and centrifuged at 13 000 g for 10 min and 4 °C. Subsequently, 350 μl of supernatant was collected in a new microtube. The pellet was run through the same steps, and the supernatant was collected. ABA identification and quantification was carried out by LC-MS/MS using an Agilent 1200 Infinity Series chromatograph coupled to a 6430 triple quadrupole mass spectrometer as described by Napoleão et al. (2017).

Proteomic analysis

For proteomic quantification, leaf samples were collected from 3-month-old non-transformed and miR156ox_4 plants. Total protein extraction was performed according to Damerval et al. (1986) with some modifications. Briefly, 10 mg of freeze-dried pulverized leaves in liquid nitrogen were mixed with 1 ml of chilled solution containing 10% (w/v) trichloroacetic acid (Sigma-Aldrich) in acetone (Merck, Darmstadt, Germany) and 20 mM dithiothreitol (GE Healthcare, Piscataway, NJ, USA), and the mixture was vortexed and centrifuged (Reis et al., 2021). The resulting pellets were washed three times with cold acetone containing 20 mM dithiothreitol and centrifuged for 5 min at each wash step. The pellets were air dried and resuspended in solubilization buffer containing 7 M urea, 2 M thiourea, 2% Triton X-100, 1% dithiothreitol (all GE Healthcare), 1 mM phenylmethanesulfonyl fluoride (Sigma-Aldrich), and complete protease inhibitor cocktail (Roche Diagnostics, Mannheim, Germany). The mixture was vortexed and centrifuged. Supernatants were collected and the protein content measured.

The protein samples were precipitated using methanol–chloroform to remove any contaminants (Nanjo et al., 2012), resuspended in a solution consisting of 7 M urea and 2 M thiourea, and digested using Microcon-30 kDa filter units (Millipore, Billerica, MA, USA) according to the filter-aided sample preparation methodology described by Reis et al. (2021). After digestion, the resulting peptides were quantified via absorbance measurements at 205 nm using a NanoDrop 2000c spectrophotometer (Thermo Fisher Scientific) and then injected into a NanoAcquity ultra-performance LC coupled to a Q-TOF SYNAPT G2-Si mass spectrometer (Waters, Wilmslow, UK).

Spectral processing and database searches were performed using Protein Lynx Global SERVER software version 3.02 (Waters), and label-free quantification analyses were performed using ISO Quant software version 1.7 (Distler et al., 2014). Protein identification was performed against a non-redundant protein databank for B. orellana generated by transcriptome sequencing and de novo assembly (Moreira et al., 2023). The mass spectrometry proteomics data have been deposited to the ProteomeXchange Consortium via the PRIDE (Perez-Riverol et al., 2022) partner repository with the dataset identifier PXD036944. Only proteins that were present or absent (for unique proteins) in all three runs of biological replicates were considered in the differential accumulation analysis. Proteins with significant Student’s t-test (two-tailed, P≤0.05) results were considered differentially accumulated proteins, as up-accumulated if the log₂ fold change (FC) was greater than 0.6 and down-accumulated if the log₂ FC was less than −0.6. For functional annotations, we used OmicsBox software (https://www.biobam.com/omicsbox) and UniProtKB (www.uniprot.org). Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment pathway analysis between differentially accumulated proteins (P≤0.01) was performed in Metascape (Zhou et al., 2019) after a BLAST search of NCBI (https://www.ncbi.nlm.nih.gov) to obtain Arabidopsis reference sequences. The predicted protein interaction networks were constructed using Arabidopsis homologs in B. orellana, which were identified through a STRING search followed by downstream analysis in Cytoscape version 3.9.1 (Shannon et al., 2003).

Statistical analysis

The software IBM SPSS Statistics v. 22 was used for all statistical analyses. For outlier detection and removal, the boxplot method was used. In the phenotyping parameters, qRT-PCR analyses, bixin and ABA measurements, and bixin channel number, the mean differences of each miR156ox line were compared with the control non-transformed, using Student’s t-test (P≤0.05). For bixin content in Fig. 1A, ANOVA followed by Tukey’s honestly significant difference (HSD) test (P<0.05) was used.

Fig. 1.

Variations in transcript levels of miR156 and SPLs are correlated with bixin content in Bixa orellana leaves. (A) Average total bixin content. Values are means ±SE (n=4). Distinct letters represent statistical differences at 5% level evaluated by Tukey’s HSD test. (B, C) Relative expression levels by qRT-PCR of B. orellana CAROTENOID CLEAVAGE DIOXYGENASE4 (BoCCD4-4) (B) and miR156 (BomiR156) (C). Values are means ±SE (n=4). (D) BoSPL open reading frames (light grey, coding sequence (CDS)), the conserved SBP domain (white) and in silico prediction of miR156 recognition site in BoSPL transcripts (black). Sequences were sized based on their nucleotide number (below the boxes). (E) Relative expression levels by qRT-PCR of BoSPL6, BoSPL13a, and BoSPL13b. Values are means ±SE (n=4). *P≤0.05 (Student’s t-test). DAP, days after planting.

Results

High bixin content in B. orellana leaves is correlated with a decline in miR156 levels

A recent report has provided evidence that bixin biosynthesis in seeds is developmentally regulated (Moreira et al., 2023). To determine whether bixin biosynthesis may be also developmentally regulated in annatto leaves, we quantified bixin levels in the fully expanded third leaf from plants at three distinct stages (45, 90, and 300 DAP). We observed an increase in bixin content as the plant ages, concomitant with an increase in the transcript levels of B. orellana CAROTENOID CLEAVAGE DIOXYGENASE4 (BoCCD4-4) (Fig. 1A, B), which is responsible for lycopene cleavage to generate bixin aldehyde (Carballo-Uicab et al., 2019). To test whether age-dependent accumulation of bixin in leaves might be regulated by the miR156–SPL node, we initially measured the abundance of BomiR156, and found that BomiR156 transcript levels dropped ~2 to 9-fold in leaves over the same period (Fig. 1C). This observation indicates that miR156 accumulation and bixin content are inversely correlated in leaves.

To characterize the molecular basis of this phenomenon in more detail, we identified the genes that encode SPL transcription factors in B. orellana. We took advantage of the recently published B. orellana transcriptomes in several organs, including leaves and seeds (Cárdenas-Conejo et al., 2015; Moreira et al., 2023). Our data revealed a set of 13 annatto SPL family members, designated as BoSPL with their respective number considering the closest related phylogenetic SPL/SBP genes (Fig. 1D; Supplementary Fig. S1; Supplementary Table S2). With the exception of the VIII clade (Supplementary Fig. S1), all 13 BoSPL genes were found in eight known clades (Salinas et al., 2012). Five annatto SPL genes (BoSPL2, 6, 10, 13a, and 13b) showed the predicted recognition site for miR156, whilst eight did not (Fig. 1D), probably due to the incomplete transcript sequence available or because these eight SPL genes are not targeted by miR156. We selected BoSPL6, -13a, and -13b as their homologs have been associated with the control of secondary metabolism in distinct species (Feyissa et al., 2019; Cao et al., 2021), and measured their transcript accumulation by qRT-PCR. PCR primers were designed to a unique sequence within the predicted coding region of each BoSPL gene containing in silico the miR156-predicted recognition site (Supplementary Table S1; Fig. 1D). We observed higher transcript levels of BoSPL6, -13a, and -13b mainly at 300 DAP (Fig. 1E), which coincided with the higher content of bixin in leaves (Fig. 1A). These data demonstrate that miR156-targeted BoSPLs are developmentally regulated and their expression patterns correlate positively with bixin levels in leaves.

Higher levels of miR156 affect growth and leaf shape of B. orellana

To evaluate whether high levels of miR156 may impact bixin biosynthesis, we employed B. orellana plants overexpressing miR156 (miR156ox; Faria et al., 2022). For phenotypic analyses, we selected two in vitro miR156ox lines with similar vegetative architecture and ramification patterns to avoid any bias that would affect subsequent growth ex vitro. To our knowledge, to date there are no studies characterizing in detail the transition between juvenile and adult phases in B. orellana trees. Therefore, we initially compared vegetative architecture and leaf development between miR156ox and control or non-transformed plants growing under greenhouse conditions. Three-month-old miR156ox plants displayed drastic modifications of canopy growth patterns, with multiple shoots arising from axillary buds (Fig. 2A). In addition, miR156ox plants had smaller leaves compared with non-transformed/control plants (Fig. 2B). These observations suggest that larger heart-shaped leaves may be a morphological marker for adult phase in B. orellana.

Fig. 2.

miR156 overexpression impacts proper canopy architecture and leaf development in B. orellana. (A) Vegetative architecture and size of non-transformed (Nt) plants and plants overexpressing miR156 (miR156ox) in two distinct lines (miR156ox_1 and miR156ox_4) at 90 days after planting (DAP) in a greenhouse. Bar: 10 cm. (B) Leaf shape and size of non-transformed (Nt) and miR156ox_1 and miR156ox_4 lines. Bar: 10 cm. From left to right: first to sixth leaf. (C) Total leaf number until 60 DAP. (D) Average phytomer number assessed up to 90 DAP in non-transformed, miR156ox_1, and miR156ox_4 plants. (E) Branching index assessed up to 60 DAP in Nt, miR156ox_1, and miR156ox_4 plants. (F) Average values for main stem height assessed up to 90 DAP in non-transformed, miR156ox_1, and miR156ox_4 plants. (G) Average values for main stem diameter until 90 DAP in non-transformed, miR156ox_1, and miR156ox_4 plants. (H) Third leaf length (from petiole insertion to the end of the leaf blade) at 75 DAP in non-transformed, miR156ox_1, and miR156ox_4 plants. (I) Leaf area at 90 DAP (n=20) in non-transformed, miR156ox_1, and miR156ox_4 plants. (J) Total leaf area (leaf area×number of non-transformed leaves (15)) at 90 DAP in non-transformed, miR156ox_1, and miR156ox_4 plants. Values are means ±SE (n=3). *P≤0.05 (Student’s t-test).

Next, we characterized in detail greenhouse-grown plants every 15 d until 90 DAP. As expected, both miR156ox lines produced more leaves than their non-transformed counterparts at all evaluated time points (Fig. 2C), as well as a higher phytomer number after 30 DAP (Fig. 2D). Quantification of the branching ratio revealed 11-fold more lateral branches in miR156ox than in non-transformed plants (Fig. 2E), indicating a strong loss of apical dominance similar to other species (Barrera-Rojas et al., 2023). Stem height was higher in miR156ox plants than in non-transformed plants until 75 DAP, although the opposite was observed at 90 DAP (Fig. 2F). Although not statistically significant at all time points evaluated, our data suggest that miR156ox plants displayed less secondary growth compared with non-transformed plants, as indicated by their reduced main stem diameter compared with non-transformed/control plants (Fig. 2G). We also evaluated leaf growth over time in non-transformed and miR156ox plants, analysing the first fully expanded leaf of B. orellana (third leaf). At 15 DAP, the third leaf of miR156ox plants grew faster than those of non-transformed plants, but the opposite was observed from 45 DAP onwards (Fig. 2H). Finally, we quantified leaf area by photographing the leaves after 90 DAP. Because of the high number of leaves in miR156ox plants (Fig. 2A–C), we used an equivalent number of leaves in miR156ox and non-transformed plants. Therefore, the leaf area for both genotypes was comparable. miR156ox leaf area was 15 times smaller than in non-transformed plants (Fig. 2I). No significant differences were observed between total miR156ox and non-transformed leaf areas (Fig. 2J), when the former was adjusted by multiplying the measured miR156ox leaf area, which is the same as for non-transformed leaves, by a factor of 15, because miR156ox plants exhibited 15 times more leaves than non-transformed plants. These findings implied that miR156ox lines produced more leaves with distinct ontogenetic growth, but with overall leaf area index similar to control plants.

Leaf shape is commonly used to distinguish juvenile versus adult identity, but light intensity can sometimes affect leaf shape without necessarily operating via age-dependent mechanisms (Jones, 1995). Therefore, we tested whether miR156ox plants are delayed in terms of leaf development by evaluating additional traits such as the presence of sugar-secreting extrafloral nectaries (EFNs). The appearance of EFNs produced by acacias (genus Vachellia) is tightly correlated with the decline in the expression of miR156 and increasing expression of SPLs, suggesting that these structures evolved by co-opting a preexisting age-dependent program (Leichty and Poethig, 2019). While non-transformed plants produced EFNs in several nodes, miR156-overexpressing B. orellana plants did not produce any visible EFNs and displayed apparently fewer trichomes in the stems (Supplementary Fig. S2, S3), indicating that high levels of miR156 impact the timing of vegetative development in achiote trees.

Annatto miR156-overexpressing leaves have reduced bixin content

Given that leaf development and ontogeny in miR156-overexpressing B. orellana plants are distinct from those of non-transformed plants (Fig. 2A; Supplementary Figs S2, S3), we investigated whether the bixin production in leaves would be affected as well. miR156ox leaves displayed ~2.5-fold lower bixin levels than control plants (Fig. 3A). Bixin is produced in all plant tissues and accumulated in secretory structures known as bixin channels (Almeida et al., 2021). We then hypothesized that modifications in leaf structures, such as secretory structures, would be specifically impacted in miR156ox plants. In the basal region of the leaf, the number of bixin channels did not differ between miR156ox and non-transformed counterparts. In contrast, there were fewer bixin channels in the medial and apical parts of miR156ox leaves (Fig. 3B, C), which may be another fundamental adult trait of B. orellana. The fewer bixin channels observed in miR156ox leaves may be a result of the lack of structural maturation of these tissues (Almeida et al., 2021). Indeed, we observed a less complex leaf venation pattern in miR156ox leaves compared with non-transformed ones (Fig. 3D; gray arrows), which reinforced the negative correlation between leaf maturation and high levels of miR156.

Fig. 3.

Density of bixin channels and the amount of bixin per leaf area are reduced when miR156 is overexpressed. (A) Average total bixin content in B. orellana leaves from non-transformed/control (Nt) plants and plants overexpressing miR156 (miR156ox_1 and miR156ox_4) at 90 days after planting (DAP). The results were adjusted to the mass of fresh leaves from four plants of each genotype. (B) Average number of bixin channels in different B. orellana abaxial leaf parts at 180 DAP. Values are means ±SE. *P≤0.05 (Student’s t-test). (C) Images of representative 90 DAP B. orellana leaves. Bars: 3 cm. (D) Detail of bixin channels of B. orellana leaves at 180 DAP. The black arrows point to bixin channels and the gray arrows indicate veins on the abaxial face of the leaf blade. Bar: 1 cm.

To determine whether the observed B. orellana leaf phenotypes were consistent with high levels of miR156, we quantified miR156 transcript levels by qRT-PCR in the greenhouse-grown plants. As expected, mature BomiR156 transcript abundance was significantly higher in both miR156ox lines compared with non-transformed plants (Fig. 4A), whereas the opposite was observed for BomiR172 (Fig. 4B). This observation is consistent with the opposing miR156 and miR172 expression patterns observed during Arabidopsis floral transition (Wang et al., 2009). We recently showed that miR172 transcript levels significantly increased in adult leaves of Passiflora edulis, whereas an opposite expression pattern was observed for miR156 (Silva et al., 2019). Collectively, these observations implied that miR156 and miR172 may have opposing roles associated with the leaf ontogenetic process in both B. orellana and P. edulis. We then evaluated transcript levels of miR156-targeted BoSPL6, -13a, and -13b (Fig. 1), of which all were repressed in miR156ox lines (Fig. 4C–E). These observations indicate that these BoSPLs may be required for proper leaf ontogeny and likely bixin biosynthesis in annatto leaves.

Fig. 4.

High expression of miR156 leads to low levels of miR172 and BoSPL transcripts in B. orellana leaves. Relative expression profiles by qRT-PCR of B. orellana miR156 (BomiR156) (A), BomiR172 (B), BoSPL6 (C), BoSPL13a (D), and BoSPL13b (E) in leaves from non-transformed (Nt) plants and plants overexpressing miR156 (miR156ox_1 and miR156ox_4) at 90 d after planting in the greenhouse. Values are means ±SE (n=3). *P≤0.05 (Student’s t-test).

Carotenoid- and bixin-associated biosynthetic genes are misexpressed in miR156-overexpressing annatto leaves

Transcriptomic profiling of distinct developmental B. orellana seed stages showed differential gene expression associated with carotenoid and bixin biosynthesis (Moreira et al., 2023). The decline in bixin production observed in the miR156ox leaves is greater than the decline in the number of bixin channels and might occur earlier in leaf development (Fig. 3). These observations suggest that, in addition to leaf ontogeny, the miR156 overexpression affects additional mechanisms that might include the control of the expression of carotenoid and bixin biosynthesis-associated genes, such as BoCCD4-4 (Fig. 1). To test this conjecture, we analysed the expression of several genes associated with carotenoid and bixin biosynthesis. Bixa orellana 1-DEOXY-D-XYLOSE-5-PHOSPHATE SYNTHASE 2A (BoDXS2a) is responsible for converting pyruvate and glyceraldehyde-3-phosphate into 1-deoxy-d-xylulose-5-phosphate, one of the first steps of lycopene (one of the bixin precursors) biosynthesis (Vaccaro et al., 2014). BoDXS2a was repressed in miR156ox_1 compared with non-transformed leaves (Fig. 5A). Transcripts of annatto PHYTOENE SYNTHASE 1 (BoPSY1), coding for an enzyme that converts dimethylallyl pyrophosphate into phytoene (Fu et al., 2014), were ~2.0-fold less abundant in miR156ox than in non-transformed counterparts (Fig. 5B). Although annatto ZETA-CAROTENE DESATURASE (BoZDS) mRNA levels did not significantly differ between miR156ox and non-transformed leaves, BoCCD4-4 transcript levels were strongly reduced in miR156ox leaves (Fig. 5C, D). The SABATH methyltransferase gene family appear to be responsible for converting methylate norbixin to bixin, and represent the end point of bixin production (Cárdenas-Conejo et al., 2015; Faria et al., 2020). SABATH1, -3, and -4 are expressed in leaves (Moreira et al., 2023), and all three SABATH genes were down-regulated in miR156ox leaves (Fig. 5E–G). In particular, BoSABATH1 expression showed ~8-fold reduction in both transgenic lines. Annatto ALDEHYDE DEHYDROGENASE FAMILY 3 MEMBER I1 (BoALDH3I1) is also responsible for converting bixin aldehyde into norbixin, and BoALDH3H1-1 transcripts accumulated in leaves and seeds (Cárdenas-Conejo et al., 2015). Low levels of BoALDH3H1-1 were observed in miR156ox compared with non-transformed leaves (Fig 5H). However, given that miR156ox_4 plants did not differ from non-transformed, we tested another BoALDH gene, BoALDH3H1, which was repressed in both transgenic lines (Fig. 5I).

Fig. 5.

miR156 overexpression impacts expression of genes associated with carotenoid, bixin, and abscisic acid (ABA) biosynthesis in B. orellana leaves. (A–M) Relative expression levels by qRT-PCR of 1-DEOXY-D-XYLOSE-5-PHOSPHATE SYNTHASE 2A (BoDXS2a) (A), PHYTOENE SYNTHASE 1 (BoPSY1) (B), ZETA-CAROTENE DESATURASE (BoZDS) (C), CAROTENOID CLEAVAGE DIOXYGENASE4 (BoCCD4-4) (D), SABATH methyltransferase genes BoSABATH1 (E), BoSABATH3 (F) and BoSABATH4 (G), ALDEHYDE DEHYDROGENASE FAMILY 3 MEMBER I1 (BoALDH3I1) (H), ALDEHYDE DEHYDROGENASE FAMILY 3 MEMBER H1 (BoALDH3H1) (I), CAROTENOID CLEAVAGE DIOXYGENASE1 (BoCCD1) (J), LYCOPENE ε-CYCLASE2 (Boɛ-Lyc2) (K), LYCOPENE β-CYCLASE 1 (Boβ-Lyc1) (L) and -2 (Boβ-Lyc2) (M). (N) ABA content. Leaves were collected from non-transformed (Nt) plants and plants overexpressing miR156 (miR156ox_1 and miR156ox_4) at 90 d after planting in the greenhouse. Values are means ±SE (n=3). *P≤0.05 (Student’s t-test).

Bixa orellana CAROTENOID CLEAVAGE DIOXYGENASE1 (BoCCD1) may divert the carbon flux away from bixin production (Cárdenas-Conejo et al., 2015). miR156ox leaves exhibited more than 2-fold increase in BoCCD1 transcripts compared with non-transformed counterparts (Fig. 5J). Carotenoid biosynthesis proceeds through two different pathways, depending on the activity of two enzymes: LYCOPENE ε-CYCLASE and LYCOPENE β-CYCLASE. LYCOPENE ε-CYCLASE constitutes a key control point for the continuation of lycopene to lutein or β-carotene biosynthesis (Cazzonelli and Pogson, 2010; Fu et al., 2019). Bixa orellana LYCOPENE ε-CYCLASE2 (Boɛ-Lyc2) mRNA levels did not differ between miR156ox and non-transformed leaves, although miR156ox showed slightly fewer transcripts (Fig. 5K). In contrast, B. orellana LYCOPENE Β-CYCLASES 1 and 2 (Boβ-Lyc1 and Boβ-Lyc2) showed higher expression in miR156ox leaves (Fig. 5L, M). Collectively, these observations suggested that low levels of miR156-targeted SPLs may interfere with β-carotene biosynthesis in annatto leaves. Given that β-carotene-derived xanthophylls serve as biosynthetic precursors of ABA (Nambara and Marion-Poll, 2005), we assessed ABA levels in non-transformed and miR156ox leaves. Both transgenic lines presented higher ABA levels when compared with non-transformed (Fig. 5N), suggesting a possible deviation of carbon flux from bixin to the ABA biosynthetic pathway.

Overall secondary metabolism is affected in miR156ox leaves as revealed by proteomic analysis

To determine how miR156 overexpression impacts the overall production of secondary metabolite-related proteins/enzymes in B. orellana leaves, we employed proteomic analysis (Supplementary Table S3). A total of 255 differentially accumulated proteins were detected, of which 89 were up-regulated and 159 were down-regulated in miR156ox leaves. Gene ontology (GO) analysis (Supplementary Table S3) revealed that approximately 6% of the low accumulated proteins were associated with cell cycle (e.g. tubulin β-1 chain). A similar trend was observed for other GO-related categories, such as cellular component biogenesis, developmental processes, and cellular component organization. Classifying differentially expressed proteins into their respective KEGG pathways (Supplementary Fig. S4) revealed that most proteins associated with sugar and nucleotide metabolisms, amino acid biosynthesis, secondary metabolite biogenesis, and ribosome biosynthesis were less abundant in miR156ox leaves (Fig. 6A). Proteins associated with flavonoid, anthocyanin, and lignin biosynthesis were the most altered by miR156 overexpression. A volcano plot analysis of the enzymes related to secondary metabolism (Supplementary Table S4) highlighted five key candidates linked to pigments and other secondary metabolite pathways (Fig. 6B). Among the high-abundance enzymes in miR156ox leaves, we found ZEAXANTHIN EPOXIDASE (ZEP), which catalyses the conversion of zeaxanthin to violaxanthin and is a key enzyme in ABA biosynthesis (Hieber et al., 2000). Low-abundance enzymes associated with secondary metabolism in miR156ox leaves were FLAVONOL SULFOTRANSFERASE-LIKE, PHENYLALANINE AMMONIA-LYASE, ANTHOCYANIDIN 3-O-GLUCOSYLTRANSFERASE 5-LIKE, and ANTHOCYANIDIN SYNTHASE, all involved in defense against abiotic and biotic stresses (Landi et al., 2015; Zhang and Liu, 2015; Qiu et al., 2016; Jin et al., 2019).

Fig. 6.

Proteins associated with secondary metabolism are differentially accumulated in miR156-overexpressing (miR156ox) leaves. (A) Functional annotation of differentially accumulated proteins. Biological processes were shown to be enriched by Fisher’s exact test at P=0.01 comparing proteins in miR156-overexpressing (miR156ox_4) and non-transformed (Nt) leaves. (B) Volcano plot for identifying secondary metabolism-related proteins. Each point represents the difference in expression (fold-change) between miR156ox_4 and non-transformed annatto plants plotted against the level of statistical significance. Red circles: down-regulated proteins; gray triangles: proteins whose levels have not changed; blue diamonds: up-regulated proteins; A: ZEAXANTHIN EPOXIDASE; B: FLAVONOL SULFOTRANSFERASE-LIKE; C: PHENYLALANINE AMMONIA-LYASE; D: ANTHOCYANIDIN 3-O-GLUCOSYLTRANSFERASE 5-LIKE; E: ANTHOCYANIDIN SYNTHASE. (C) STRING protein–protein interaction diagram. Protein–protein interaction network from miR156ox and non-transformed plants comparing with putative roles in secondary metabolism. The network summarizes predicted associations for a particular group of proteins based on Arabidopsis orthologs identified by STRING using the B. orellana amino acid sequence. Only first neighbors are shown in the interactions. All the predicted interactions in the network showed high confidence (>0.7) according to the combined STRING score. Different colors indicate different KEGG pathways. The straight line shows the interaction between proteins. Yellow line: literature evidence; green line: neighborhood evidence; purple line: experimental evidence; light blue line: database evidence; and black line: co-expression evidence.

Proteins shown in the Supplementary Table S4 were used to determine mutual functional relationship, which identified four different groups (Fig. 6C): group I (TT5, DFR, F3H, LDOX, FLS1, BAN, CHIL, and TT4) are associated with the biosynthesis of amino acids, and with general secondary metabolism; group II (JAT, PAL1, CCoa, CAD5, OMT1, ELI3-1, and GAD9) are linked to the biosynthesis of phenylpropanoid and lignin; group III (LYC, ZDS, DXR, HDS, CDPMEK, and CCD1) are related to carotenoid and ABA biosynthesis, classified also as terpenoid backbone biogenesis; and group IV (DHS2, DHSI, AT2G, MEE32, EMB1144, AT3G, and ASA2) are enzymes implicated in the biosynthesis of specific amino acids. As revealed by the composition of group III, the main enzymes associated with carotenoid biosynthesis were clustered together (Fig. 6C). Collectively, our data indicated that the biosynthesis of ABA and carotenoids was impacted in the miR156ox leaves, likely by modifying the accumulation of several enzymes, including those from the BoCCD1 and/or BoCCD4-4 pathways (Fig. 5). Future studies using miR156-resistant BoSPL versions and knock-out mutants may shed light on the precise role of the miR156–SPL pathway in modulating bixin and ABA contents in B. orellana leaves.

Discussion

miR156-targeted BoSPL activity is likely required for B. orellana vegetative phase change

Bixa orellana plants show few obvious vegetative modifications from juvenile to adult developmental stage (Baliane, 1982), which makes it challenging to determine modifications during leaf and shoot development. Determining how phase-specific traits contribute to B. orellana fitness is fundamental for better management of this economically important tree. For instance, traits associated with phase change such as leaf shape, leaf tissue maturation, and the emergence of stem structures, are closely related to the improvement in the chances for plant survival and reproduction (Lawrence et al., 2023). Here, we characterized the role of the miR156–SPL module, and found that modifications during vegetative phase change are not very different from those observed in Arabidopsis, an annual herb. miR156-overexpressing B. orellana plants exhibited smaller leaf blades shaped as triangles rather than hearts (as commonly displayed by adult B. orellana plants; Baliane, 1982), lack of extrafloral nectaries, and a rather thinner bark layer compared with control plants (Fig. 2). In poplar (Populus spp.), other woody species with subtle vegetative phase changes, the overexpression of miR156 also leads to more evident shoot modifications (Rubinelli et al., 2013; Lawrence et al., 2021). Collectively, these observations reinforce that vegetative phase change is an ancient and fundamental aspect of shoot development, which might respond to natural selection (Leichty and Poethig, 2019).

The higher number of leaves observed in miR156ox annatto plants is a common phenotype observed in SPL-deficient plants, such as Petunia and tomato (Silva et al., 2014; Zhou et al., 2021). Thus, it is possible that the low expression levels of SPLs are the main causative factor underlying the miR156ox annatto plant phenotypes. Smaller leaf size, which indicates modifications in leaf ontogeny, might be a result of differential modulation of genes associated with organ size that are regulated by the miR156–SPL module (Zhang et al., 2015; Fouracre and Poethig, 2019; Barrera-Rojas et al., 2020; Fang et al., 2021; Ferigolo et al., 2023). Even though miR156ox lines have more leaves, they exhibited similar photosynthetic area to their non-transformed counterparts, a phenomenon also observed in Eucalyptus grandis (Levy et al., 2014), Solanum tuberosum (Bhogale et al., 2014), and Arabidopsis (Xie et al., 2017). High levels of miR156, together with the reduced expression levels of BoSPLs and BomiR172 (Fig. 4), suggested that even under greenhouse conditions, miR156ox annatto plants displayed similar overall gene regulation as in vitro-raised plants (Faria et al., 2022). Importantly, BomiR172 transcripts were less abundant in miR156ox leaves, which might correlate with their delay in leaf maturation, as observed for rice and Arabidopsis (Zhu et al., 2009; Jung et al., 2011). However, a recent report showed that AtmiR172 is more associated with reproductive maturation than vegetative phase change, which indicates that modifications in Arabidopsis leaf ontogeny may be independent of the miR172 regulation by the miR156–SPL module (Zhao et al., 2023).

miR156-targeted BoSPL activity is likely required for proper bixin and ABA biosynthesis in B. orellana leaves

The levels of bixin in fully expanded leaves (similar developmental stage) correlated negatively with miR156 accumulation when the plant ages (Fig. 1). Accordingly, miR156ox leaves exhibited fewer bixin channels, which likely impacted bixin content (Fig. 3). Although bixin production is apparently regulated independently in annatto seeds and leaves, recent evidence demonstrates the presence of reddish latex rich in carotenoids secreted by pigment glands/anastomosed articulated laticifers (Almeida et al., 2021) or carotenoid storage cells in different organs (Louro and Santiago, 2016; Moreira et al., 2023). In seeds, these secretory structures are particularly numerous (Lima et al., 2013), and therefore seeds produce significantly higher amounts of the pigment. The secretion of bixin by ‘pigment glands’ or carotenoid storage cells in seeds, leaves, and other organs suggests a conserved genetic pathway modulating bixin production (Moreira et al., 2023). A similar pattern is observed in Eucaluptus brevistylis, in which the secondary metabolite content, oil quantity, and gland maturation are closely related to leaf age (Goodger et al., 2018). In a similar way, levels of anthocyanin, flavones, and flavonols were modified in miR156-overexpressing poplar plants (Wang et al., 2020). Bixin gland ontogenesis depends on the gradual coalescence of cells with shared cell walls, which become degraded as seen in laticifer channels (Almeida et al., 2021). Such process depends on organ maturation (Lima et al., 2013; Lee and Ding, 2016), which is delayed in miR156ox plants.

In addition to changes in leaf ontogeny, the miR156–SPL module may also function in the regulation of bixin biosynthesis-associated enzymes. Although this study lacks transcriptome data from annatto leaves, by combining specific gene expression and proteomic analyses, we shed light on how the miR156-targeted BoSPLs might regulate carotenoid- and bixin-associated pathways at the molecular level. The main results are summarized in the Fig. 7. For instance, the overexpression of miR156 attenuated the transcription of BoCCD4-4 and SABATH genes, which encode key enzymes in bixin biosynthesis (Brandi et al., 2011; Varghese et al., 2021). Moreover, genes coding for lycopene-associated downstream enzymes were down-regulated in miR156ox plants, except for BoALDH3I1 (Fig. 5), whose up-regulation may be explained by the BoALDH3 expression being controlled via the ABA stress response pathway rather than directly by the miR156–SPL module (Brocker et al., 2013). The ALDH gene family has additional functions in plants, as shown in Arabidopsis, where ALDHs play a central role in salt tolerance through putrescine and γ-aminobutyrate production (Zarei et al., 2016). It will be interesting in future studies to determine how distinct bixin biosynthetic pathway-associated genes are temporally expressed in fully expanded leaves to confirm their age-dependent patterns.

Fig. 7.

The miR156–SPL module may control bixin and ABA biosynthesis. Summary of carotenoid-related transcripts, proteins/enzymes, and metabolites modified by the overexpression of miR156 in B. orellana leaves. Dotted line: indirectly related to the route; continuous line: directly related in the metabolic pathway. Proteins/enzymes names as described in Fig. 5.

Reduction of BoSPL expression in miR156ox leaves may have led to the up-regulation of enzymes involved in β-carotenoid production, including BoCCD1 (Fig. 5), which might be responsible for diverting carbon from bixin (Cazzonelli and Pogson, 2010; Cárdenas-Conejo et al., 2015). MuSPL16 modulates secondary metabolism in banana by directly promoting LYCOPENE Β-CYCLASE expression and thereby increasing β-carotene content in the fruits (Zhu et al., 2020). Here, we showed that high levels of miR156 seem to increase LYCOPENE Β-CYCLASE gene expression in leaf tissue (Fig. 5). It is possible that BoSPLs, such as BoSPL6 and/or BoSPL13, might directly control annatto β-Lyc1 and β-Lyc2 expression to modulate carotenoid and ABA biosynthesis in leaves. This conjecture deserves further investigation in future studies.

Given that bixin and ABA contents were altered in miR156ox leaves (Figs 3, 5), we investigated whether additional secondary metabolic pathways were affected by the overexpression of miR156. LC-MS analysis revealed that the abundance of several secondary metabolism-associated proteins was reduced in miR156ox leaves (Supplementary Fig. S4B), with a few exceptions such as SUPEROXIDE DISMUTASE (SOD), which responds to high ABA levels (Zhang et al., 2008; Salazar-Chavarría et al., 2020). ZEP, which is crucial for ABA biosynthesis (Marin et al., 1996; Hieber et al., 2000), also accumulated at high levels in miR156ox leaves (Fig. 6). Thus, it is possible that the accumulation of ZEP and SOD in miR156ox leaves reflects the rerouting of carbon flux from bixin production towards ABA biosynthesis during the juvenile phase. There are relatively few examples of juvenile and adult molecular traits, such as secondary metabolites, that might contribute to plant fitness or be expected to be selectively advantageous (Leichty and Poethig, 2019; Lawrence et al., 2021). Our observation that B. orellana adult leaves produce more bixin and less ABA than juvenile ones suggests that it may be advantageous to accumulate stress-related metabolites, such as ABA, early in development, when annatto trees are more susceptible to biotic and abiotic stresses. In contrast, increased apocarotenoids such as bixin in leaves may be advantageous once adult trees are established and less susceptible to stresses.

The miR156–SPL module has been well characterized in numerous species and is highly conserved in angiosperms (Morea et al., 2016). However, how different metabolites are developmentally regulated by this miRNA regulatory module is still unclear. Here, we illustrate one of these fine-tuned controls of secondary metabolic pathways in B. orellana, a worldwide economically important species whose molecular characterization is still in its infancy.

Supplementary data

The following supplementary data are available at JXB online.

Fig. S1. Maximum-likelihood phylogeny of the SPL/SBP proteins of Bixa orellana, Arabidopsis, tomato, cacao, grape, rice, and switchgrass.

Fig. S2. Bixa orellana extrafloral nectary development is negatively impacted by the miR156 overexpression.

Fig. S3. Bixa orellana extrafloral nectary and trichome formation are negatively impacted by the miR156 overexpression.

Fig. S4. Proteomic analysis of Bixa orellana miR156-overexpressing (miR156ox_4) and non-transformed (Nt) leaves.

Table S1. Oligonucleotides used in this study.

Table S2. List of SPL coding genes and their respective gene name for each plant species analysed.

Table S3. List of proteins identified in miR156ox_4 and non-transformed/control (Nt) plants, and catalogued by Gene Ontology function analysis.

Table S4. List of secondary metabolism-associated proteins extracted from the proteome list.

Acknowledgements

We thank the Chr. Hansen Ind. Com. Ltda. (Valinhos, SP, Brazil) for kindly donating the seeds of Bixa orellana Piave Vermelha, and Dr A. R. Leichty for kindly making available the miR156-overexpressing vector. We also thank the Universidade Federal de Viçosa (UFV) for providing all relevant infrastructure, the Núcleo de Análise de Biomoléculas (NUBIOMOL) for ABA and Bixin quantification, Gilmar Valente from Núcleo de Microscopia e Microanálises (NMM) for microscopy analyses, and the Genética Molecular de Bactérias Laboratory (DMB—BIOAGRO) for RT-qPCR equipment.

Author contributions

KLGM, MBSD, DVF, FTSN, and WCO designed the study. KLGM, DVF, and MBSD, performed the experiments. KLGM, DVF, TRO, TCAF, ER, and WCO analysed the data. LASS, TRO, DVF, WCO, and FTSN provided technical assistance. KLGM, LASS, DSB, DVF, TRO, MGCC, CSC, VS, ER, WCO, and FTSN wrote the manuscript.

Conflict of interest

There are no conflicts of interest to declare.

Funding

This work was supported by The Fundação de Amparo à Pesquisa do Estado de Minas Gerais—FAPEMIG (Grant numbers APQ-02372-17 and APQ-00772-19), the Conselho Nacional de Pesquisa e Desenvolvimento Científico e Tecnológico (CNPq), Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (CAPES) (Finance Code 0001), and partially by FAPESP (Grant no. 18/17441-3). KLGM received a scholarship from CNPq (146061/2019-5). DVF received a postdoctoral scholarship from CNPq (155615/2018-1).

Data availability

Data on the mass spectrometry proteomics are available at the Proteome X Change Consortium via PRIDE (https://www.ebi.ac.uk/pride/archive), under the dataset identifier PXD036944. All primary data to support the findings of this study are openly available upon request.