Transcriptome-Wide Analyses Provide Insights into Development of the Hedychium coronarium Flower, Revealing Potential Roles of PTL

The flower of Hedychium coronarium possesses highly specialized floral organs: a synsepalous calyx, petaloid staminodes and a labellum. The formation of these organs is controlled by two gene categories: floral organ identity genes and organ boundary genes, which may function individually or jointly during flower development. Although the floral organogenesis of H. coronarium has been studied at the morphological level, the underlying molecular mechanisms involved in particular organ morphologies that emerge in flower development still remain poorly understood. Here, we used comparative transcriptomics combined with Real-time quantitative PCR to investigate gene expression patterns of ABC-class genes in H. coronarium flowers, as well as the homolog of PETAL LOSS (HcPTL). Our studies found that stamen/petal identity or stamen fertility in H. coronarium was not necessarily correlated with the differential expression of HcAP3 and HcAG. We also found a novel spatio-temporal expression pattern for HcPTL mRNA, suggesting it may have evolved a lineage-specific role in the morphogenesis of the Hedychium flower. Our study provides a new transcriptome reference and a functional hypothesis regarding the role of a boundary gene in organ fusion that should be further addressed through phylogenetic analyzes of this gene, as well as functional studies.


Introduction
The Zingiberales is a monophyletic order that includes the so-called "ginger group" containing the Marantaceae, Cannaceae, Costaceae, and Zingiberaceae, and the "banana group" containing the remaining four families (Musaceae, Strelitziaceae, Heliconiaceae, and Lowiaceae) (Dahlgren and Rasmussen 1983;Kirchoff 1988a;Kress 1990). The four families in the ginger group share a synapomorphy: petaloid staminodes take the place of most of the stamens in mature flowers (Kirchoff 1997;Kirchoff et al. 2009). These petaloid staminodes, sharing positional homology with stamens in the banana families but sharing some properties commonly associated with the perianth rather than with the androecium, were considered to take on the function of perianth whorls, such as pollinator attraction (Kirchoff 1991).
The Zingiberaceae encompass species with high ornamental, culinary and medicinal value, which commonly possess colorful and morphologically diverse flowers. A typical Zingiberaceae flower possesses four concentric whorls of floral organs: the outer perianth whorl (whorl 1) which consists of three sepals; the inner perianth whorl (whorl 2) composed of three petals; a single fertile stamen and five staminodes comprise the androecial whorl (whorl 3), in which 2-4 petaloid staminodes normally fuse to form a distinctive labellum; and the innermost whorl (whorl 4), which bears a tricarpellary pistil. In contrast with the calyx and corolla that are visually inconspicuous in the mature flower of Zingiberaceae, the petaloid staminodes and the resulting labellum usually constitute the bulk of the overall 1 3 floral display. The single fertile stamen and the five staminodes are derived from two concentric androecial whorls, with the inner androecial whorl organs sharing a common primordium with the petal primordia, and the outer androecial members then developing in a space between the petal/ inner androecial whorl (Kirchoff 1988a). Hedychium coronarium, a Zingiberaceae species with a showy flower shape and inviting fragrance, is widely cultivated in tropical and subtropical regions due to its high ornamental and medicinal value. H. coronarium was considered as an ideal model to study the origin of particular floral forms bearing highly specialized floral organs. While the floral organogenesis of H. coronarium has been studied in depth at the morphological level (Kirchoff 1997), the underlying molecular mechanisms involved in its floral organ development are still unknown.
Previous research documented the importance of MADSbox transcription factors in the control of floral organ identity specification during flower development in the model plant Arabidopsis thaliana (Bowman et al. 1991;Coen and Meyerowitz 1991;Weigel and Meyerowitz 1994;Theiβen 2001). Insights from this body of research gave rise to the ABC model of floral organ identity, which posits that the differentiation of the canonical floral organs in a flower is due to combinatorial interactions among the products of a set of transcription factors, all of which are MADS-box genes (except for APETALA2 [AP2]) (Bowman et al. 1991). Thus, from the outmost part of the floral primordium to its center, sepals are specified through the expression of A-class genes APETALA1 [AP1] and AP2; petals are specified through the concerted action of A plus B-class genes PISTILLATA [PI] and APETALA3 [AP3]; stamens are differentiated through the action of B plus the C-class gene AGAMOUS [AG] and the carpel is specified by the sole expression of AG. Furthermore, A-class and C-class genes repress each other (Bowman et al. 1991). The ABC model provided a conceptual framework for comparative developmental genetics, thus, the role of orthologs of ABC-class MADS-box genes has been investigated in many non-model angiosperms, showing that some variants have evolved, notably: A-class genes are not necessary for calyx formation in many plants and rather, AP2-like genes seem to have conserved a cadastral role (restricting the domain of expression of the C-class gene in the outer whorls), while the perianth and stamen whorls can have nuanced interactions among B and C-class genes that have enabled variable morphologies (Kramer 2019). In the Zingiberales, several MADS-box genes that control the specification of floral organ identities have been characterized to explore their potential roles in the formation of petaloid structures including staminodes and the resulting labellum (Bartlett and Specht 2010;Song et al. 2010;Yockteng et al. 2013;Almeida et al. 2013Almeida et al. , 2015Fu et al. 2014). Recent research suggests that particular events of gene duplication and differential losses of paralogs could underlie the tendency towards staminode laminarization across the different families that comprise the order (Almeida et al. 2015;Piñeyro-Nelson et al. 2017). However, this body of research has focused on analyzing the functional evolution of floral MADS-box genes and its relation with floral organ inception and diversification in members of the Zingiberales, and as such cannot uncover the molecular mechanisms underlying particular organ morphologies that emerge during late flower development, such as the synsepalous calyx fused by the margins of adjacent sepal primordia through intercalary growth or the variable fusion patterns among petaloid staminodes that will form the labellum in Zingiberaceae species (Zhao et al. 2020).
Floral organ fusion is a pervasive macroevolutionary trend in angiosperms and provides great potential for flower diversity (Endress 2011;Phillips et al. 2020). In eudicots, floral organ fusion has been mapped to mutations on different genes with the function of boundary delimitation, such as CUP-SHAPED COTYLEDON 1-3 (CUC1-3), LATERAL ORGAN BOUNDARIES (LOB) and PETAL LOSS (PTL) (Aida et al. 1997;Shuai et al. 2002;Takeda et al. 2011;Lampugnani et al. 2012). However, there are only a few functional studies focusing on the role of boundary gene orthologs in monocot floral development, mainly focused in the grasses. While a recent study in Zingiber zerumbet has identified several boundary genes which may contribute to the formation of the synsepalous calyx in this species (Zhao et al. 2020), how these gene exert their function during floral organ differentiation and growth still needs further exploration. Notably, PTL is a plant-specific trihelix transcription factor involved in regulating perianth architecture in the Arabidopsis flower (Brewer et al. 2004). Mutations in the PTL gene have a phenotype of petal defects and sometimes fused sepals, indicating that it plays a role in promoting petal initiation and boundary delimitation between adjacent sepal primordia (Griffith et al. 1999;Brewer et al. 2004). In Arabidopsis, the PTL transcripts were detected in four inter-sepal zones in the early-developing flower, where they repressed organ growth in regions between newly initiating sepals, ensuring their separation (Brewer et al. 2004). Until now, the expression pattern and potential function of PTL gene during flower development of Zingiberaceae plants are still unclear. In addition, the existence of common primordia that emerge during early floral development in many Zingiberales (e.g., species in Zingiberaceae, Costaceae and Heliconiaceae) (Kirchoff 1988b(Kirchoff , 1997Kirchoff et al. 2009), is a characteristic of many floral lineages such as species in the Caryophyllales, Ranunculales, etc. The gene expression patterns that induce floral organ differentiation from common primordia remain poorly understood. Genes such as CUC and PTL may play potential roles in the differentiation of petals and stamens from these common primordia, as well as in the fusion of laminar structures such as staminodes.
H. coronarium, as a non-model plant and a suitable system to study the origin of particular floral forms (e.g., having a common primordium, petaloid staminodes and fused floral organs), thus, was selected in this study to investigate the molecular and genetic mechanisms involved in floral development and organ fusion. While a reference genome is lacking and few candidate-gene studies have been performed in H. coronarium, the molecular mechanisms involved in its flower development can be investigated through the use of comparative genomic methods such as transcriptome analysis. RNA-seq, a cost-effective and highly efficient technology built on small-read sequencing, has accelerated gene discovery across a diverse set of species, while extending our comprehension of complex gene regulatory networks to a wide breath of organisms from the tree of life. In recent years, there has been an explosion in the use of plant transcriptomes as part of a comparative analysis pipeline to isolate homeotic gene orthologs, as well as to study differential gene expression among flowers at different developmental stages in a variety of non-model species, including some within the Zingiberales (Zhang et al. 2013;Huang et al. 2015;Li et al. 2017;Almeida et al. 2018;Su et al. 2018;Xiao et al. 2019;Zhao et al. 2020). Here, we conducted de novo transcriptome sequencing of two floral developmental stages for H. coronarium and used comparative transcriptomics combined with Real-time quantitative PCR (qRT-PCR) to investigate gene expression patterns of orthologs of selected ABC-class homeotic genes. We further assayed the spatio-temporal gene expression pattern of HcPTL (the homolog of PTL in H. coronarium), across different stages of floral organ differentiation through mRNA in situ hybridization assays. The results presented in this paper provide important bioinformatic resources for further investigation of genes and gene networks in this species and evidence for the potential functional diversification of a boundary gene that could underlie some of the floral diversity observed across the Zingiberales.

Plant Materials and RNA Extraction
Young inflorescences of H. coronarium were collected from June to August in 2017 (for transcriptome sequencing and qRT-PCR) and 2018 (for scanning electron microscopy, SEM), respectively, from South China Botanical Garden, Chinese Academy of Sciences (SCBG, CAS) at Guangzhou, China. For transcriptome sequencing, samples of flower primordia (Hp) and flowers showing organ differentiation (Hd) were removed from inflorescences, flash frozen in liquid nitrogen and maintained at -80 ℃ for subsequent RNA extraction. Each sample was combined from three inflorescences to reduce the effect of variation between individuals. Three biological replicates were performed for each Hp and Hd sample, respectively. Dissection of different floral organs for qRT-PCR analyses was performed on one-centimeter-long inflorescences. Total RNA used for transcriptome sequencing was extracted with a mirVana™ miRNA Isolation Kit (Thermo Fisher Scientific, USA) following the manufacturer's protocol and treated with DNase I (TianGen, China).

Scanning Electron Microscopy
Inflorescences with multiple floral buds were collected and dissected under a stereomicroscope. Bracts and larger floral organs were removed, and the floral buds were preserved in a fixative with 2.5% glutaraldehyde and 2% paraformaldehyde made with 0.1 M phosphate buffer overnight. All floral buds were washed in a 0.1 M phosphate buffer three times in 2 h and dehydrated in an alcohol series (30%, 50%, 70%, 80%, 90%, 100%, 100%, 100%). The materials were critical point freeze-dried with a Leica EM CPD 300, mounted on stubs, gold-coated in a Leica EM ACE600, and observed under a JSM-6360LV SEM (JEOL) operated at 10 kV.

Sequencing and De Novo Assembly
H. coronarium mRNA was enriched by Oligo(dT) beads and fragmented using divalent cations under elevated temperature prior to sequencing. The first-strand cDNA was synthesized with random primers, followed by second-strand cDNA synthesis with DNA polymerase I, RNase H, dNTPs and buffer. The double-strand cDNA fragments were purified with a QiaQuick PCR extraction kit, followed by end repair and poly(A) addition (He and Jiao 2014). Illumina sequencing adapters were ligated to the fragments, size selected by agarose gel electrophoresis and PCR amplified (He and Jiao 2014). The libraries were constructed using the Illumina paired-end sequencing technique in the size of 150nt on an Illumina HiSeq™ 4000 platform by Gene Denovo Biotechnology Co. (Guangzhou, China). Clean reads were obtained by removing adapters, unknown nucleotides, and low-quality sequences. The de novo assembly was performed using Trinity v2.1.0 as previously described for de novo transcriptome assembly in the absence of a reference genome (Grabherr et al. 2011). 1 3 protein annotations. GO annotations were analyzed with the Blast2GO 2.3.5 program based on the Nr annotation (Conesa et al. 2005). RPKMs (Reads per Kilobase per Million mapped reads) were calculated to reflect the relative expression level of each unigene (Mortazavi et al. 2008). Differential expression analyses were operated using the edgeR 3.12.1 (Robinson et al. 2009). Genes with a fold change ≥ 2 and a false discovery rate (FDR) < 0.05 were considered as differentially expressed genes (DEGs). GO enrichment and KEGG pathway analysis of DEGs were performed using an online website OmicShare tools (http:// www. omics hare. com/ tools/). The DEGs were presented as a heatmap, constructed using the program TBtools (Chen et al. 2018).

Real-Time Quantitative PCR
Twelve unigenes (unigene0006394, unigene0036304, uni-gene0032128, unigene0066414, unigene0032127, uni-gene0000995, unigene0043696, unigene0048812, uni-gene0033259, unigene0054002, unigene0052504 and unigene0038064) were selected for validation using Realtime quantitative PCR (qRT-PCR). The primers were designed with the Integrated DNA Technologies Prim-erQuest tool (https:// sg. idtdna. com/ prime rquest/ home/ index) and are listed in Supplementary Table S1. The qPCR reactions were performed using Top Green qPCR Super-Mix (TransGen Biotech, China) in a total volume of 10 μl, with a reaction mixture containing 5 μl of Top Green qPCR SuperMix, 0.2 μl (10 μmol/l) of each primer, 0.2 μl Passive Reference Dye II (50 ×), 4 μl of water, and 0.4 μl of cDNA template. Each reaction was performed in three technical replicates. The qPCR assays were performed in an ABI Prism 7500 Fast Real-time PCR Detection system under the program of 95 °C for 30 s, followed by 40 cycles of 95 °C for 15 s and 60 °C for 35 s. The dissociation stage condition was set at 95 °C for 15 s, 60 °C for 60 s, 95 °C for 15 s, and 60 °C for 15 s to test for primer specificity. All the data were normalized against the expression of the reference gene β-Actin. The relative quantities of transcripts were calculated by the comparative Ct method ( 2 −ΔΔC t ) (Livak and Schmittgen 2001). Graphs were constructed using Origin9.0.

Construction of Phylogenetic Tree
The evolutionary history was inferred using the Neighbor-Joining method. The bootstrap consensus tree inferred from 1000 replicates is taken to represent the evolutionary history of the taxa analyzed. Evolutionary analyses and phylogenetic tree were conducted in MEGA6. The CDS sequence of A. thaliana is from TAIR (https:// www. arabi dopsis. org/). The CDS sequences of Canna indica, Zingiber zerumbet, Ravenala madagascariensis, Musa ornata, Heliconia bourgaeana, Thalia dealbata, Cheilocostus speciosus, Alpinia galanga and Globba schomburgkii are from transcriptome (some are unpublished).

Description of Developmental Events
To better elucidate the early stages of floral organogenesis in H. coronarium, we divided the continuous process of its flower development into 10 stages (as shown in Fig. 2), using a series of landmark events. We fully considered the uniqueness of the developmental process of H. coronarium flower while referring to the stages of early flower development proposed for Arabidopsis thaliana (Smyth et al. 1990), as well as the previous flower development descriptions made for H. coronarium by other authors (Kirchoff 1997). Preliminary investigations have demonstrated considerable variability in the pattern of flower development in H. coronarium, much of which is correlated with the position of the individual flower on the inflorescence (Kirchoff 1997). Given this variability, the detailed description put forth is based on the development of the first flower of the cincinnus, developing in the central portion of the inflorescence (see Fig. 1c; f1). Stage 1 begins with the formation of the floral primordium (fp). Stage 2 commences when sepal primordia (sep1, sep2, sep3) sequentially arise. A ring primordium (rp) forms interior to the sepals (Stage 3). Three common primordia (cp) then separate from the ring primordium and the posterior common primordium first separate into a posterior petal (p) and a fertile stamen (s) (Stage 4). Two abaxial common primordia (cp2, cp3) then differentiate into a petal to the exterior and an inner androecial member to the interior. The synsepalous calyx (ca) and the flower cup form (Stage 5). When the calyx gradually encloses all internal floral organs, the carpel primordia (c) initiates style formation from the flower cup region. The outer androecial primordia initiate (Stage 6). The abaxial outer androecial member ceases growth and the two inner androecial members (ia1, ia2) later form the labellum (Stage 7). The adaxial outer androecial members form the petaloid staminodes (Stage 8). During Stage 9, the style grows rapidly and becomes longer than the stamen and labellum. When flower development is almost complete, the mature flower becomes resupinate, with the labellum oriented adaxially (Stage 10). The complete developmental stages described above can be divided into two phases: the early developmental phase, encompassing 7 stages that end with the formation of the labellum, by this time all floral organs have initiated.
The late phase of flower development (Stage 8-10), refers to the growth, cellular expansion and maturation of the differentiated floral organs.

Transcriptome Sequencing and Assembly
To study the molecular mechanisms of development in the flower of H. coronarium at a transcriptome-wide scale, two developmental stages were analyzed: undifferentiated flower primordia (Hp, stage 1) and differentiated flowers at different developmental stages (Hd, encompassing stages 2-9). In all, 287,535,456 raw reads were generated, and a total of 276,075,968 clean reads were recovered after low-quality reads and adapters were removed. The reads of six RNA-seq libraries (Hp and Hd, with three biological replicates for each) were pooled together for subsequent de novo assembly. A total of 80,909 unigenes with a mean size of 983 bp and N50 length of 1770 bp were obtained, the length of which ranged between 201 and 16,715 bp. Sixteen-thousand six-hundred and fifteen (16,615) unigenes were ≥ 1.5 Kb in length and 4815 unigenes were ≥ 3 Kb in length. The length distribution of unigenes can be found in Supplementary Table S2.

Gene Annotation and Differential Expression Analysis
Sequences of assembled unigenes were searched against public protein databases, using an E-value threshold of 1E-5. As shown in Fig. 3a, 38 ,492, 27,169, 22,809 and 14,768 unigenes were matched to the Nr, Swiss-prot, KOG and KEGG databases, respectively. A total of 38,743 unigenes were found in at least one database, while 11,278 unigenes shared similarity to sequences deposited in all four reference databases. Species similarity based on Nr annotation demonstrated that the closest match was with the Musa acuminata transcriptome (51.5%, 19,825 unigenes, Fig. 3b), a species from the same order (Zingiberales) as Hedychium. The next-closest matches were to other monocot taxa: Elaeis guineensis (4.8%,1839 unigenes) and Phoenix dactylifera (4.4%,1695 unigenes) (Fig. 3b). The unigenes assigned to the KOG database were classified into 25 functional categories, of which the largest number was assigned to the general function Fig. 2 Ten developmental stages of the H. coronarium flower. Stage 1: floral primordium forms; Stage 2: sepal primordia arise; Stage 3: ring primordium forms; Stage 4: adaxial petal and fertile stamen arise/abaxial common primordia form; Stage 5: abaxial petals and inner androecial primordia arise. Flower cup and ring calyx form; Stage 6: carpel primordium and outer androecial members initiate; Stage 7: abaxial outer androecial member ceases growth and labellum forms; Stage 8: adaxial outer androecial members form petaloid staminodes; Stage 9: style extends rapidly; Stage 10: floral organs mature and resupinate. ab abaxial side, ad adaxial side, ax inflorescence axis, c carpel, ca calyx, CM cincinnus meristem primordium, cp common primordium, fc flower cup, fp floral primordium, ia inner androecial primordium, l labellum, oa outer androecial primordium, p petal, ps petaloid staminode, rp ring primordium, s stamen, sb secondary bract, se sepal, sti stigma, sty style. * indicates aborted stamen. The color of carpel represents its longitudinal depth. The darker the color, the longer the style. The colors of SEM are false colors added to highlight certain floral organs or organ parts category, followed by the category of transcripts related with signal transduction mechanisms (Fig. 3c).
Comparison of relative gene expression levels revealed that 215 unigenes demonstrated significant expression differences between the floral primordium and differentiated flowers (Hp vs Hd). Among all these DEGs, 187 were significantly up-regulated and 28 were significantly downregulated. A GO term enrichment analysis was conducted for all of these DEGs (Fig. S1). In the category pertaining biological processes, the analyzed DEGs were largely enriched in the metabolic process, cellular process and single-organism process. In the molecular function category, binding and catalytic activity were the two dominant enriched terms. In the cellular component category, the most represented items were cell, organelle and cell part.

RNA-seq Expression Validation by qRT-PCR
To validate the relative gene expression levels detected for different transcripts in our bioinformatic analysis, twelve unigenes were selected for real-time quantitative PCR analysis. The detailed information of the selected unigenes and primers are shown in Supplementary Table S1. All of the twelve unigenes assayed were successfully amplified, showing single bands of the expected sizes. Furthermore, the expression patterns detected by qRT-PCR were consistent with the overall trends detected in the RNA-seq analysis (Fig. S2). These results demonstrate the robustness and credibility of the quantitative RNA-seq data generated in this study.

Differential Expression of ABC-Class Genes in Developing Floral Organs of H. coronarium
We performed additional qRT-PCR assays to investigate the organ-specific expression patterns for several ABC-class genes on RNA purified from floral organs dissected from one-centimeter long inflorescences (Fig. 5). The relative expression levels of HcAP1, HcAP3 and HcAG unigenes were assayed in the six different organs that develop in four concentric floral whorls, namely the sepal (whorl 1), petal (whorl 2), petaloid staminode (whorl 3), labellum (whorl 3), stamen (whorl 3) and carpel (whorl 4). The A-class gene HcAP1 was expressed in all six floral organs, showing the highest level of expression in the petal. The B-class gene HcAP3 was expressed in the second and third floral whorls including petals, petaloid staminode, labellum and stamen; while, its expression in sepal (whorl 1) and carpel (whorl 4) was almost undetectable. The C-class gene HcAG was expressed in the third and fourth floral whorls, including the petaloid staminode, labellum, stamen and carpel, while its expression was hardly detected in the first and second floral whorls. The qRT-PCR patterns of expression in the four floral whorls are consistent with the patterns depicted in the classic ABC model (Bowman et al. 1991) as well as those documented in other angiosperms (Kramer 2019).

Phylogenetic Analyses
Phylogenetic analyses (Fig. 6) indicate that our recovered sequences fall into two clades-banana group and ginger group (Except A. thaliana). CDS sequence analysis shows close phylogenetic relationship of four Zingiberaceae species including A. galanga, Z. zerumbet, G. schomburgkii and H.coronarium, which are clustered to one branch in neighbor-joining tree. C. speciosus (Costaceae) is clustered with four Zingiberaceae species. T. dealbata and C. indica, respectively, from Marantaceae and Cannaceae (two sister families), are clustered together.

HcPTL Expression During Early Development of H. coronarium Flowers
The spatio-temporal expression pattern of HcPTL during early floral development (Stage 1-Stage 7) was assessed through an mRNA in situ hybridization assay for this transcript (Fig. 7). In Stage 1 (Fig. 7a, b), HcPTL was expressed throughout the floral primordium, as well as in the primordium of the cincinnus meristem (Fig. 7b). At an intermediate stage between Stage 2 and Stage 3 (Fig. 7c, d), the expression of HcPTL was maintained in the cincinnus meristem primordium (Fig. 7d). Within the developing flower meristem HcPTL was expressed in new arising sepals, and it was also detected in the domain where the petal and fertile stamen will later emerge from the ring primordium (Fig. 7d). At Stage 4, when petal and the fertile stamen were already clearly visible (Fig. 7e, f), HcPTL continued to be expressed in the sepals and the inter-sepal boundary. At this stage, the expression pattern of HcPTL in the common primordium is seemingly correlated with the emerging dorsoventrality in the floral meristem, as the expression of this transcript was higher in the dorsal part of the common primordium (Fig. 7f). Similarly, the expression level of HcPTL in the petal derived from the dorsal part of the common primordium was relatively higher than in the fertile stamen, which derived from the ventral part of this primordium (Fig. 7f). At stage 5 (Fig. 7g, h), HcPTL was continually expressed throughout the petal and fertile stamen, and the expression domain also included the synsepalous calyx (Fig. 7h). Notably, the signal was almost undetectable in the inner androecial whorl (Fig. 7h). By stage 6 (Fig. 7i, j), HcPTL expression decreased in the petal and calyx; at this stage, the expression of HcPTL was confined to the ventral part of Fig. 4 Predicted transcription factors. The inner whorl represents the relative expression level of unigenes in Hp while the outer whorl represents the relative expression level of unigenes in Hd. Genes were clustered by RPKM value with color scale representing the log-transformed RPKM value. a Expression pattern of MADS-box transcription factors. M-type and MIKC-type categories represented in orange and blue bars (innermost whorl), respectively. * indicates A-class (orange), B-class (blue), C-class (green) and E-class (purple) genes. b Expression pattern of differentially expressed transcription factors during flower development. Color bars in the innermost whorl represent different transcription factors family and the corresponding colors refer to the TF family identification (legend right) ◂ the fertile stamen, with weak signals detected in its dorsal part (Fig. 7j). A faint expression of HcPTL was detected in the initiating carpel. By stage 7, when all floral organs had already initiated (Fig. 7k, l), HcPTL expression could only be detected in the fertile stamen, covering an area that is 3-4 cells wide. Due to the lack of longitudinal sections of the carpel at this stage, we have no further knowledge of the expression of HcPTL in this organ.

Analysis of the Floral Transcriptome of H. coronarium Yields Insights into its Flower Development
RNA-Seq has become an important tool to unravel the molecular players underlying floral diversification in non-model plants that lack a reference genome, such as the Zingiberales. In recent years, whole flower transcriptome libraries of exemplar species from the six families (Musaceae, Lowiaceae, Zingiberaceae, Costaceae, Marantaceae and Cannaceae) comprising the Zingiberales have been generated and analyzed in previous studies (Almeida et al. 2018;Zhao et al. 2020). Our transcriptome data are consistent in terms of unigene numbers, unigene length distributions, and quality scores with the aforementioned studies. Moreover, given that we sequenced two different developmental floral stages in H. coronarium, we were able to assess differential gene expression, finding that 215 unigenes showed differences between the Hp and Hd stages, with 187 significantly up-regulated and 28 significantly down-regulated. The GO analysis for these DEGs annotated many of them as related to cellular processes and metabolism. Another notable finding is that while the annotated 2142 transcription factor sequences were a small fraction of the overall unigenes recovered (2.65%), we found 56 gene families represented. Interestingly, we also found some of  (25) stage. We also detected 57 transcription factors that were differentially expressed between stages (Fig. 4b), and many of these are MADS-box genes (discussed further in next section).
Comparing our unigenes with previous studies where lineage-specific genes had been retrieved only in certain lineages within the Zingiberales, some of the most interesting genes were homologs to LOB-domain genes, a plant-specific transcription factor family which has been implicated in defining organ boundaries between Arabidopsis flower organs through the negative regulation of brassinosteroid accumulation (Shuai et al. 2002;Bell et al. 2012), as well as root development, auxin signaling and stress response (Zhang et al. 2020). The shared floral transcription factors recovered for all species across the exemplar Zingiberales species analyzed, such as LOB40, 41 and 6-like homologs (Almeida et al. 2018), were also identified in the H. coronarium floral transcriptome, suggesting that these genes could be conserved regulators of floral development across the Zingiberales. In contrast with previous findings, a LOB36like gene which had been formerly recovered only in the banana clade (Almeida et al. 2018), was recovered in our H. coronarium. We also recovered a LOB18-like sequence, which was previously recovered only in the Canna-Calathea lineage (Cannaceae-Marantaceae lineage) and in Zingiber officinale (Zingiberaceae) (Almeida et al. 2018). As expected given its close phylogenetic relationship, the LOB4-like gene, which had been previously recovered only in Z. officinale (Almeida et al. 2018) and Z. zerumbet (Zhao et al. 2020), was also retrieved in the H. coronarium floral transcriptome. A CUC2-like gene, previously recovered only in the Canna-Calathea and Zingiber lineages (Almeida et al. 2018), was also identified in our data for H. coronarium; while PTL, previously only recovered in Zingiber, was also retrieved in our Hedychium floral transcriptome, consistent with their close phylogenetic relationship. Most of the genes highlighted above, have been implicated in floral organ boundary formation in A. thaliana (Bell et al. 2012;Aida et al. 1997;Lampugnani et al. 2012), but knowledge of their roles outside of core eudicots is still scarce.
Although we did not sequence more developmental stages on account of the difficulty in obtaining plant materials, we believe the transcriptome analyses presented here add to the body of data that can provide a crucial bioinformatic basis for future functional studies on the molecular mechanisms of flower development and diversification in the Zingiberaceae.

MADS-Box Genes Involved in Floral Organ Identity Specification Recovered in H. coronarium
In this study, we recovered 51 MADS-box genes, including those involved in floral organ specification according to the ABC model (Bowman et al. 1991). Previous comparative developmental studies addressing the functional conservation of the orthologs of these MADS-box genes have been conducted in different non-model angiosperm lineages, including the Zingiberales, unraveling their role in floral meristem determinacy and floral organ identity specification in both eudicots and monocots (Yamaguchi et al. 2006;Song Fig. 7 HcPTL expression during early development of H. coronarium flowers. Green lines on scanning electron micrographs (a, c, e, g, i and k) indicate planes of section for in situ hybridizations. mRNA in situ hybridizations (b, d, f, h, j and l) were detected by HcPTL antisense probes in longitudinal sections of H. coronarium flowers in Stage 1 (b), an intermediate stage between Stage 2 and Stage 3 (d), Stage 4 (f), Stage 5 (h), Stage 6 (j) and Stage 7 (l). c carpel, ca calyx, CM cincinnus meristem, cp common primordium, fp floral primordium, ia inner androecium, la labellum, oa outer androecium, p petal, rp ring primordium, s stamen, sb secondary bract, se sepal, sty style. Bars = 100 μm et al. 2010;Heijmans et al. 2012;Almeida et al. 2013;Fu et al. 2014;Tian et al. 2016). In the Zingiberales, the families comprising the ginger group show a trend towards stamen modification into petaloid structures, as in all lineages petaloid staminodes replace several fertile stamens in the two concentric stamen whorls present in developing flowers (Kirchoff 1997;Kirchoff et al. 2009). This kind of homeotic conversion was first considered to be possibly caused by the differential expression of B-and C-class genes in the petaloid staminodes, thereby allowing a shift in organ identity (Wake et al. 2011). However, evidence provided by later studies in Alpinia hainanensis (Zingiberaceae) and Canna indica (Cannaceae) demonstrated that the C-class gene is expressed in both the petaloid staminode and the fertile stamen (Song et al. 2010;Almeida et al. 2013;Fu et al. 2014;Tian et al. 2016), indicating that differential gene expression may not be the main reason for the shift in organ identity. In our work, the expression patterns of the homologs of ABC-class genes analyzed in H. coronarium is in line with previous observations. Based on our quantitative transcriptomic data, the A-class gene HcAP1 displays a higher expression in the undifferentiated floral primordium; while the B and C-class genes HcAP3 and HcAG, respectively, are significantly up-regulated in flowers with differentiated floral organs. The relatively high expression of HcAP1 in the floral primordium might be due to the observation in other taxa that this gene plays an important role in maintaining floral meristem identity and promoting sepal differentiation, functioning at earlier stages before sepal initiation, or, alternatively, that it plays a cadastral role, limiting the domain of expression of the C-class gene (Kramer, 2019). In contrast, B-and C-class genes exert their function in the inner three whorls (petal, stamen and carpel whorls), which was evidenced by a higher expression in organ-differentiated flowers. qRT-PCR detection in different floral organs showed that the expression of the A-class gene HcAP1 was detected not only in sepals, but also in all four floral whorls, with a higher expression in petals (Fig. 5a). The concerted expression of HcAP3 (B-class) and HcAG (C-class) genes was detected in floral organs corresponding to whorl 3, both in the fertile stamen and the sterile stamen-derived petaloid structures (petaloid staminode and the resulting labellum), with no significant differences in their expression levels across these structures (Fig. 5b, c); HcAG, as expected, had a higher expression in the developing carpel (Fig. 5c), where HcAP3 expression was absent.
Our results support the suggestion of A-class gene function as discussed in the (A)BC model put forth by Causier et al. (2010), and further indicate that the formation of petaloid structures and stamen fertility are not necessarily correlated with the differential expression of B-and C-class genes. Instead, independent gene duplication events in the B-and C-class AG lineage, which have been documented in several eudicot and monocot taxa, including species within the Zingiberales (Kramer et al. 2004;Zahn et al. 2006;Almeida et al. 2015), may potentially lead to functional divergence among paralogs or dissociation among genetic modules thereby increasing the genetic substrate that can lead to morphological diversity (Kramer et al. 2004). For instance, the duplication of the AG gene followed by functional divergence of the duplicated AG copies during the diversification of the ginger families, may explain some of the observed morphological changes in Zingiberales flowers, indicating possible mechanisms for the evolution of androecial petaloidy within this order (Almeida et al. 2015;Piñeyro-Nelson et al. 2017). Future research should address if paralogs of the B-class genes exist in H. coronarium as well as investigate the function of a potential AG paralog detected in the transcriptome analyses (Fig. 4a).
Additionally, while MADS-box transcription factors expressed in developing flowers have been studied in multiple species, these genes have been ignored for their potential roles in floral organ pigmentation. A recent study in orchids found that the expression levels of AP3 and AGL6 displayed more differences between lip segments (the yellow hypochile v. the purple-red epichile), than between the sepal/petal segments which have a similar coloration during floral development, suggesting that AP3 and AGL6 might participate in color differentiation in the perianth and lip segments during flower development ). Similar to orchids, the labellum of different genera in the Zingiberaceae displays a wide variation in color and shape, as well as ontogenetic identity (Zhao et al. 2020). In this context, the potential involvement of MADS-box genes in flower color differentiation across different floral organs has yet to be studied. The MADS-box transcription factors retrieved in this study provide a global candidate-gene list for future studies that can address not only the molecular mechanisms underlying floral organ identity specification, but also assess if they play a role in pigmentation in the different petaloid organs present in Zingiberaceae species.

HcPTL may Promote Growth in the Developing Floral Meristem and Floral Organs
Although PTL loss of function could influence the number of petals, PTL expression could not be detected in early-developing petals, suggesting that the observed petal defects associated with PTL loss of function may be indirect, perhaps influenced by the overgrowth of nearby inter-sepal zones (Brewer et al. 2004). Thus, while PTL has been regarded as a boundary gene modulating connation between adjacent organs within the same floral whorl, its activity has not been widely investigated in monocots. If PTL functions as a conserved boundary-specification regulator, it should not be expressed in inter-sepal boundaries during the formation of the synsepalous calyx in the Zingiberaceae. However, the expression pattern of a PTL-like homolog in developing Hedychium flowers contrasts with what has been described in Arabidopsis. The transcripts of HcPTL were still detected in inter-sepal boundaries. In addition, the expression of HcPTL was detected in other structures, such as in developing meristems, including cincinnus primordia, floral primordia, common primordia and almost all new initiating floral organ primordia. This indicates HcPTL may function not only in modifying calyx architecture but throughout early floral organ development in other whorls. Furthermore, while PTL in Arabidopsis typically represses cellular growth within its domain of expression (boundaries), HcPTL could have acquired a novel role, as it seems to be promoting growth in developing meristems and floral organs in H. coronarium, as well as organ differentiation in common primordia. Notably, the floral organ-specific qRT-PCRs for HcPTL showed a pattern of expression more consistent with data from A. thaliana, namely, a higher expression in sepals and petals at late developmental stages (See Fig. S3). This divergence in function of a boundary gene with respect to its documented role in A. thaliana floral development, where HcPTL seems to promote growth in floral organs rather than limiting it, was documented in another boundary gene related with the CUC1-3 genes, CUPULIFORMIS (CUP) in the Antirrhinum flower (Rebocho et al. 2017). In this context, it might be possible that the slight growth rate difference observed between petal and stamen/inner androecium members in H. coronarium (Fig. 6e-h) might be explained, at least in part, by the differential expression of HcPTL in the dorsal/ventral part of the common primordia, so as to ensure that the petals envelop the stamen/inner androecium. Thus, a distinct mechanism of HcPTL boundary delimitation could have evolved in the morphogenesis of the Hedychium flower.

Conclusions
While patterns of gene expression can provide insights related with the formation of particular floral structures in species like H. coronarium, a detailed study of orthologs and their patterns of molecular evolution across a broad taxonomic sampling in monocots and other lineages, could yield insights into its functional diversification and its role in organ fusion and/or differentiation in various taxa. Furthermore, functional annotation of transcripts sequenced from Zingiberaceae species would greatly benefit from the existence of a reference genome from a species within this family, as Musa acuminata, the closest reference genome available, is in the basal clade of the Zingiberales and also presents no androecial petaloidy or conspicuous floral organ fusions.
The genetic transformation of select experimental taxa across the Zingiberales to generate loss of function mutants for particular genes, could help unravel the molecular modules that underlie floral organ diversification in this group. Due to the suite of developmental novelties present in Zingiberaceae flowers, it will be exciting to establish a model system, such as H. coronarium, where assays based on different genetic modification techniques, such as virus-induced gene silencing in concert with tissue culture or transformation using CRISPR-Cas technology could greatly increase our understanding of the molecular bases of floral organ determination, growth and fusion.