Introduction
Carnation (Dianthus caryophyllus L.) is a major ornamental cut flower, potted plant, and bedding plant. Malformed flowers were detected in commercially cultivated potted ‘Cherie’ carnation plants in 2011 and 2012 in Japan (Yamane et al., 2018). There are several types of floral malformations (e.g., flat receptacles; phyllody-like proliferated sepaloids; undeveloped petals; proliferation of petals, stamens, and pistils; and proliferated flower formation). These abnormal floral phenotypes are often observed more often in spring, particularly during the Mother’s Day season, than during summer (Fig. 1). In a previous study (Yamane et al., 2018), we analyzed the effects of several temperature treatments on the growth of wild-type (WT) and malformed flower (mlf) ‘Cherie’ plants. Normal phenotype flowers of ‘Cherie’ have 5 vice sepals, a calyx tube, 40 petals, 40 stamens, a pistil and 4 stigmas (Fig. 1A). In contrast to the normal flowers of the WT plants, there were various malformed floral phenotypes among the mlf plants; the level of abnormality varied and some flowers showed intermediate phenotypes. Therefore, double flowers with the normal number of organs and shape were defined as normal, and other flowers were defined as malformed. The percentage of mlf plants with malformed flowers was as follows: 0% at 26°C, 1.3% at 25°C/20°C, 3.1% at 20°C, and 92.2% at 15°C (Yamane et al., 2018). This suggests the external temperature influences the development of malformed carnation flowers, especially between 15 and 20°C.
An earlier examination of shoot apices using optical and scanning electron microscopy revealed the morphological differences between WT and mlf plants following sepal formation (Yamane et al., 2018). After the sepal formation stage, various flower malformations were detected among the mlf plants, including undeveloped petals, undeveloped or irregularly developed stamens, formation of secondary flower primordia, and completely irregular arrangement of undeveloped flower organs. According to PCR analysis, phytoplasmas were undetectable in the mlf plants, suggesting they were not responsible for the observed malformations (Yamane et al., 2018).
Blake (1962) analyzed the carnation cultivars ‘Maytime’, ‘Puritan’, and ‘Miller’s Yellow’ and detected diverse floral abnormalities, including adventitious buds within the flower, carpel-like structures in the stamen whorl, excessive development of ovules, and axillary secondary flower formation. Low-temperature stress (5°C) reportedly induces the formation of secondary carnation flowers in the center of the original flower, resulting in a substantial increase in the total number of petals (Garrod and Harris, 1974). The flower malformation rate of mlf plants increases at approximately 15°C, which is within the temperature range suitable for carnation growth. As a result, limiting the occurrence of flower malformations at carnation production sites is difficult.
In carnation, super-double-flower phenotypes can be generally produced by crossing between normal lines. The reported bullhead or super-double-flower phenotypes (Imai, 1938; Saunders, 1917; Yagi et al., 2014a) are similar to the ‘Cherie’ mlf floral phenotypes (Yamane et al., 2018). The carnation floral phenotype has been suggested to be a monogenic trait determined by the D (double) locus. More specifically, the homozygous recessive genotype [dd] results in the single-flower phenotype, whereas the heterozygous genotype [Dd] and the homozygous dominant genotype [DD] result in the double-flower and super-double-flower phenotypes, respectively (Imai, 1938; Saunders, 1917). The self- and cross-pollination of the ‘Spectrum’ carnation cultivar resulted in the monogenic segregation of the double-flower and bullhead flower phenotypes in a 3:1 ratio (Imai, 1938). Thus, Imai (1938) proposed that two kinds of genes are related to the double-flower phenotype. The dominant allele Da (double a) in heterozygous and homozygous genotypes results in the standard double-flower phenotype and the bullhead-double-flower (i.e., super-double-flower) phenotype, respectively (Imai, 1938). Another dominant allele (Db) is associated with flower doubling, but it only results in the standard double-flower phenotype, even in a homozygous genotype (Imai, 1938). Although several markers associated with the double-flower phenotype were developed (Onozaki et al., 2006; Scovel et al., 1998; Yagi et al., 2011, 2014a), the genes responsible for the double-flower phenotype in carnation plants and their contribution to the super-double-flower phenotype in ‘Cherie’ mlf plants remain unknown.
In ‘Baccara’ rose, malformed bullhead flowers contain more petals than normal flowers (Moe, 1971). Cytokinin activities are higher in the flowers of ‘Baccara’ rose plants grown under relatively cool conditions and in malformed bullhead flowers, which develop in such conditions, than in the normal flowers of ‘Baccara’ plants grown at higher temperatures (Zieslin et al., 1979). However, there are few reports on the role of cytokinins in carnation malformed flowers.
Floral homeotic protein DEFICIENS (DEFA), which is a flower-specific homeotic gene in Antirrhinum majus, controls petal organogenesis in the second whorl, as well as stamen organogenesis in the third whorl (Schwarz-Sommer et al., 1992; Sommer et al., 1990; Zachgo et al., 1995). In Antirrhinum, flowers of the temperature-sensitive DEFA mutant have sepaloid petals and carpelloid stamens at 25–26°C, whereas the floral morphology resembles that of the WT control at 15°C (Schwarz-Sommer et al., 1992; Zachgo et al., 1995). The floral homeotic genes AGAMOUS (AG) from Arabidopsis thaliana and DEFA from Antirrhinum have sequences that are very similar to those of transcription factors (TFs) (Ma et al., 1991).
There is no report describing mutant carnation flower phenotypes that are induced by a temperature change of 5°C from the optimal growing temperature. The genetic basis of the physiological and molecular changes in these malformed carnation flowers remains uncharacterized. Moreover, the genetic determinants of the double-flower phenotype have not been identified in carnation (Yagi et al., 2014a). We hypothesized that temperature-dependent floral malformations are due to mutations to genes related to morphological responses to temperature changes. Yagi et al. (2014b) published the carnation genome sequence, which is a valuable resource for further research. In the present study, transcriptome sequencing (RNA-seq) and whole-genome re-sequencing analyses of the young flower buds and leaves of normal plants and plants with malformed flowers were performed based on the carnation genome data to test our hypothesis. Differentially expressed genes (DEGs) were functionally annotated using the Gene Ontology (GO) (Harris et al., 2004) and Kyoto Encyclopedia of Genes and Genomes (KEGG) (Kanehisa and Goto, 2000) databases to characterize the genes associated with flower malformations. To further clarify the super-double-flower phenomenon, the phenotypes of the F1 lines derived from WT and mlf plants were examined. We herein discuss the candidate genes responsible for floral malformations and their contribution to the super-double-flower phenotype.
Materials and Methods
Plant materials
Potted WT and mlf ‘Cherie’ carnation (Dianthus caryophyllus L.) plants supplied by Japan Agribio Co., Ltd. (Shizuoka, Japan) were used in all experiments. The mlf plants, which had malformed flowers, were selected from farmers’ fields and propagated from cuttings. The WT and mlf plants were transplanted into 12 cm pots and grown under cool white fluorescent artificial light (LH-350S; Nippon Medical & Chemical Instruments Co., Ltd., Osaka, Japan), with a 12-h day/12-h night cycle. Ten WT and mlf plants were grown at 15 or 20°C. The growth medium was standard peat moss (pH 5.0–6.0, 180 mg·L−1 N, 120 mg·L−1 P, and 220 mg·L−1 K; Sakata Seed Corp., Ltd., Yokohama, Japan). A commercial tablet fertilizer (Promic, N:P:K 8:12:10; Hyponex Japan Corp., Ltd., Osaka, Japan) and nutrient solution (Hyponex, N:P:K 6:10:5; Hyponex Japan Corp., Ltd.) were added to the growth medium.
DNA isolation and whole-genome re-sequencing
Genomic DNA was extracted from the young leaves of WT and mlf plants according to the CTAB method (Chang et al., 1993). For each sample, 1.0 μg DNA was used as the input material. Sequencing libraries were generated using a TruSeq Nano DNA HT Sample Preparation Kit (Illumina, San Diego, CA, USA), with index codes added to attribute sequences to each sample. The genomic DNA was randomly fragmented (350 bp) using a Covaris cracker (Covaris M220, M&S Instruments Inc., Osaka, Japan). After repairing the ends and adding an A-tail, the DNA fragments were ligated to a full-length adapter for Illumina sequencing and then amplified by PCR. The PCR products were purified using the AMPure XP system (Beckman Coulter, Brea, CA, USA) and then the libraries were analyzed in terms of the size distribution using a 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA) and quantified by real-time PCR. The high-quality libraries were pooled according to their effective concentrations and expected data volumes and then analyzed using a HiSeq 2000 system (Illumina). The FASTQ data were cleaned using FASTX-Toolkit (minimum quality > 28, 80%). Additionally, BWA was used to map the data to the carnation reference genome (ftp://ftp.kazusa.or.jp/pub/carnation/), while Samtools was used for variant calling. The re-sequencing data sets supporting the study results are available in the DDBJ Sequence Read Archive (accession number DRA008509).
RNA isolation
To screen for DEGs between the WT and mlf flower buds, five replicates of 5–7 young flower buds (2–3 mm) at the early stages of initiation were randomly selected from WT plants grown at 15 or 20°C (normal phenotype) and mlf plants grown at 15°C (malformed phenotype) and immediately frozen in liquid nitrogen before being stored at −80°C. The frozen samples were ground to a fine powder in liquid nitrogen. Total RNA was extracted from 100 mg ground material as previously described (Chang et al., 1993), after which its integrity was evaluated in a 1.0% agarose gel stained with ethidium bromide. The total RNA was quantified and examined for protein contamination (A260/A280) and reagent contamination (A260/A230) using a NanoDropTM ND-1000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA). The RNA integrity number was determined using a 2100 Bioanalyzer (Agilent Technologies) (Copois et al., 2007).
Library preparation and Illumina sequencing
After verifying the high quality of the extracted RNA, total mRNA was enriched using oligo (dT) beads. The mRNA was randomly fragmented in a fragmentation buffer and then the first cDNA strand was synthesized using mRNA as the template and random hexamer primers, after which a custom second-strand synthesis buffer (Illumina), dNTPs, RNase H, and DNA polymerase I were added to synthesize the second cDNA strand. After repairing the cDNA ends and adding an A-tail, as well as a sequencing adapter, the size selection and PCR enrichment steps were performed to complete the construction of the cDNA library. The mRNA-seq libraries were constructed and then sequenced on the HiSeq 2000 system (150 bp paired-end reads) (Illumina) according to the manufacturer’s instructions.
Bioinformatics analysis
Transcript expression values were obtained by mapping the sequencing reads to the carnation reference genome (ftp://ftp.kazusa.or.jp/pub/carnation/) using the HISAT2 (version 2.0.5) and StringTie (version 1.3.1c) programs (Pertea et al., 2016). The ballgown and edgeR R packages were used to identify DEGs. More specifically, the fragments per kilobase per million reads (FPKM) values and a false discovery rate threshold of < 0.05 were used to identify any significant DEGs in the WT and mlf samples. The non-redundant ORFs predicted for each library were merged and produced as a gtf file using the StringTie-merge program. A fasta file was also created by the gtf file and the carnation genome (DCA_r1.0_scaffolds.fa) using the gffread program. The sequences in the fasta file were used for BLASTn searches and annotation against an NCBI nr/nt nucleotide database (NCBI non-redundant sequence database). Further, amino acid sequences extracted using the fasta file including code regions by TransDecoder were annotated by BLASTp to protein databases such as Swiss-Prot (E-value cut-off was set at 1e−5). The data sets for RNA-seq are available in the DDBJ Sequence Read Archive (accession number DRA008510).
To predict the biological functions of unigenes, the coding regions were characterized by GO annotations involving the Pfam database prior to the GO enrichment analysis. The results of the GO enrichment analysis were visualized by constructing Z-score bar plots for specific GO terms. High and low Z-scores reflected increases and decreases in GO term enrichment, respectively. The enriched KEGG pathways among the unigenes were determined using KofamKOALA (KEGG Orthology And Links Annotation) (Aramaki et al., 2020).
Quantitative gene expression analysis
A Maxwell® RSC Plant RNA Kit (Promega, Madison, WI, USA) was used to extract total RNA from 100 mg buds that had been stored at −80°C and lysed in CTAB buffer. The extracted RNA was reverse transcribed to cDNA using a ReverTra Ace® qPCR RT Kit (Toyobo Co., Ltd., Osaka, Japan). The cDNA was stored at −20°C until use. For the quantitative real-time PCR (qPCR) analysis, the cDNA samples were diluted to 20 ng·μL−1 in 100 μL Tris-EDTA buffer. The cDNA concentrations were determined using a NanoDropTM One spectrophotometer (Thermo Fisher Scientific). The qPCR analysis was conducted using a KAPA SYBR Fast qPCR Kit (Roche Diagnostics, Basel, Switzerland) and a LightCycler 480 system (Roche Diagnostics). An actin gene was selected as the reference control for analyzing WUSCHEL (WUS) and Heat Shock Cognate 70 kDa (HSC70) expression levels (Table S1). The qPCR program (95°C for 3 min; 45 cycles of 95°C for 3 s and 60°C for 20 s) was followed by melting curve analysis.
Quantification of trans-zeatin contents by LC-MS/MS
The endogenous cytokinin contents in flower primordia were determined by analyzing trans-zeatin. A stereomicroscope was used to obtain six replicates of flower primordia (0.5–6 mg) from the flower buds of WT and mlf plants grown at 15 or 20°C. The samples were immersed in 50 μL methanol containing 50 pg·μL−1 D5-trans-zeatin (internal standard). The extract was purified using an HLB column and eluted using 1 mL 50% acetonitrile. The eluate was concentrated by vacuum filtration. The residue was dissolved in 20 μL 50% acetonitrile and transferred to a vial. The resulting solution was injected into an ultra-performance liquid chromatography (UPLC) system connected to a tandem mass spectrometer. The UPLC separation was performed using a Nexera X2 instrument (Shimadzu, Kyoto, Japan). The chromatographic separation was completed using a C18 column (Kinetex F5, 2.1 × 150 mm, 2.6 μm; Phenomenex, Torrance, CA, USA). The gradient elution was performed as follows using solvent A (acetonitrile) and solvent B (water), both of which contained 0.1% (v/v) acetic acid: 10% A and then an increase to 50% A at 10 min, followed by a 5 min increase to 100% A. The column was operated at 30°C with a flow rate of 0.2 mL·min−1. The mass spectrometry analysis was conducted using a triple quadrupole/linear ion trap instrument (QTRAP 5500; AB Sciex, Framingham, MA, USA) with an electrospray ionization source. For the analysis of trans-zeatin, the MS/MS spectra were recorded in the product ion scan mode. The ion source was maintained at 400°C. trans-Zeatin was quantitatively analyzed on the basis of multiple reaction monitoring, with monitoring transitions of m/z 220 to 136 (trans-zeatin) and m/z 225 to 137 (D5-trans-zeatin).
Phenotypic examination of the F1 progeny of the WT and mlf plants
To confirm the inheritance of the super-double-flower phenotype, ‘Cherie’ was crossed with ‘Fossett rose’ (commercial cultivar; pink and double flowers for potted plants; Japan Agribio Co., Ltd.) and ‘Buryule’ (commercial cultivar; orange with red margin and double flowers for potted plants Japan Agribio Co., Ltd.), which produce substantial amounts of fertile pollen, because ‘Cherie’ flowers lack fertile pollen. The seedlings derived from the F1 seeds were grown in a growth chamber (20°C) or in a greenhouse (ambient temperature). The phenotypes of the F1 progeny were compared on the basis of the floret appearance.
Genotypes of markers closely linked to the flower-type locus
To identify the flower-type locus in the WT and mlf plants, the size of SSR marker loci (Yagi et al., 2014a) were measured using whole-genome re-sequencing data.
Results
The RNA-seq analysis generated 41.28 Gb raw data. After filtering, 39.97 Gb clean data were retained, with 65,744,861 and 67,487,318 clean reads obtained for the WT and mlf flower buds, respectively (Tables 1 and S2).
The DEGs were subjected to a GO enrichment analysis and assigned to the biological process (BP), cellular component (CC), or molecular function (MF) categories (Table S3). In the BP category, transmembrane transport and protein phosphorylation were the enriched GO terms in mlf, whereas translation and protein folding were the enriched GO terms in WT (Fig. 2A). In the CC category, the anaphase-promoting complex and membrane were the enriched GO terms in mlf, whereas the ribosome and intercellular structure were the enriched GO terms in WT (Fig. 2B). In the MF category, metal ion binding, zinc ion binding, and protein kinase activity were the enriched GO terms in mlf, whereas the ribosome structural constituent was the enriched GO term in WT (Fig. 2C).
The KEGG analysis of the upregulated genes in mlf indicated 9-cis-epoxycarotenoid dioxygenase (NCED) was associated with carotenoid biosynthesis (Fig. S1). Other genes with significantly upregulated expression levels in mlf, including Pyrabactin Resistance 1-Like (PYL), protein phosphatase 2C (PP2C), and Calmodulin (CAM) 4, were revealed to be related to MAPK signaling-plant (Fig. S2). In contrast, the expression of Histone H4 was downregulated in mlf (data not shown).
The FPKM values were used to select candidate DEGs related to floral malformations (Table S4). In mlf, the expression levels of TF genes, including WUS and STERILE APETALA (SAP), were upregulated, which was in contrast to the downregulated expression levels of PISTILLATA-like protein (PI), MADS-box protein CMB2 (CMB2), the F-box family member UNUSUAL FLORAL ORGANS (UFO), and the transcriptional activator REVEILLE 8 (RVE8) (Fig. 3; Table S4). In mlf, the HSC70 and Temperature-induced lipocalin-1 (TIL) expression levels were upregulated, whereas the expression levels of genes encoding histones (H2A, H2AX, H2B, H3.2, and H4) and ribosomal proteins were downregulated (Fig. 4; Table S4). The expression levels of NCED1, Abscisic acid receptor PYL8 and 9, Protein DEHYDRATION-INDUCED 19 homolog 3 (DI19-3), NAC domain-containing protein 2 (NAC002), and Probable WRKY transcription factor 57 (WRKY57) were upregulated in mlf (Fig. 5A–F; Table S4). Similarly, the expression levels of Cytokinin riboside 5´-monophosphate phosphoribohydrolase LOG1 (LOG1), Two-component response regulator ARR12-like (ARR12), and AP2/ERF and B3 domain-containing transcription factor RAV1 (RAV1) were also upregulated in mlf (Fig. 5G–I).
The whole-genome re-sequencing analysis generated 24.40 and 23.98 Gb raw data for the WT and mlf leaf samples, respectively (Table 2). Considering the reported carnation genome size (622 Mb) (Yagi et al., 2014b), the coverage was 39.2× and 38.6× for WT and mlf, respectively (Table 2). Sequence variations in TF including WUS and Calmodulin-binding transcription activator 5 and NAC002, heat shock-related protein including HSC70, Histone H2A and H2B, plant hormone, ribosome, metal-binding protein, transporter, and CAM are listed in Table S5.
Single nucleotide polymorphisms (SNPs) and insertions/deletions (InDels) were detected in the upstream regions and exons of WUS (Table S6; Fig. S3), HSC70-1 (Fig. S4), and HSC70-2 (Fig. S5). The relative expression levels of both WUS and HSC70 were higher at 20°C than at 15°C, with no significant differences between WT and mlf (Fig. S6).
According to the LC-MS/MS analysis, the trans-zeatin contents did not differ significantly between the mlf and WT flower primordia at 15 or 20°C (Fig. S7). Additionally, examination of the F1 progeny derived from the crosses between ‘Cherie’ (WT and mlf) and ‘Fossett rose’, as well as ‘Buryule’, which produce fertile pollen grains, showed normal double-flower phenotypes and several malformed double-flower phenotypes, but no single-flower phenotype (Fig. S8). There was a significant difference between WT and mlf in terms of the ratio of malformed phenotypes among the F1 progeny (Table S7). More specifically, the estimated normal-to-malformed flowers ratios for the F1 progeny were 7:1 and 3:1, while the did not fit ratio was 3:1 and 1:1 for WT and mlf, respectively (Table S7). In both the WT and mlf genomes, CES1982 consisted of 151 bp and CES0212 consisted of 265 bp. In other words, these loci markers did not differ between the WT and mlf genomes (Table S8).
Discussion
Roles of WUS and floral homeotic genes
The RNA-seq data indicated the WUS and SAP expression levels were upregulated in mlf, whereas the PI, CMB2, and UFO expression levels were downregulated (Fig. 3A–E). Earlier research showed WUS encodes a homeodomain protein that is required for maintenance of the cell population in Arabidopsis shoots and floral meristems (Battey and Tooke, 2002; Laux et al., 1996; Mayer et al., 1998). Moreover, WUS positively regulates the size of each shoot meristem by maintaining the appropriate number of pluripotent stem cells (Ikeda et al., 2009). In the flower development stage, the primordia for the final floral organs (carpels) develop, after which WUS is no longer expressed, resulting in the termination of floral stem cell maintenance (Smyth et al., 1990). The repression of WUS expression by AG is essential to terminate the floral meristem. Additionally, floral determinacy depends on a negative autoregulatory mechanism involving WUS and AG that terminates stem cell maintenance (Lenhard et al., 2001). In the current study, however, AG was similarly expressed in the WT and mlf flower buds at 15°C (data not shown). In Arabidopsis, SAP belongs to a class of transcriptional regulators essential for a number of floral development-related processes (Byzova et al., 1999). Together with the organ identity gene AG, SAP is required for the maintenance of floral identity, functioning in a manner similar to that of APETALA1 (AP1; i.e., an A class gene) (Byzova et al., 1999; Krizek and Meyerowitz, 1996). It is possible that the upregulated expression of WUS (Fig. 3A) and SAP (Fig. 3B) contributes to the proliferation of sepaloids or petaloids in mlf plants.
In the present study, sequence variants were detected in WUS (Tables S5 and 6; Fig. S3). The ectopic expression of WUS increases the size of shoot meristems and induces the accumulation of stem cells, which leads to the formation of adventitious shoots and somatic embryos in root tissues (Gallois et al., 2004; Zuo et al., 2002). In some mlf plants, shoot proliferation was also observed as the number of floral organs increased (data not shown). The receptacles tended to be larger in the mlf flowers than in the WT flowers, which may explain why the mlf florets were flat. Earlier research suggested the epigenetic network may also regulate WUS expression, thereby influencing floral meristem determinacy (Cao et al., 2015). However, in both WT and mlf, WUS was more highly expressed at 20°C than at 15°C (Fig. S6). The expression of WUS could not be validated, possibly because the timing of expression differed slightly among samples. Nevertheless, WUS was not consistently highly expressed in the mlf flower buds. Moreover, SNPs and InDels were detected in the upstream regions and the exons of WUS in mlf. Notably, there were ‘CAG’ repeats between intron 1 and exon 2 (Table S6). These sequence variants may influence the expression and/or function of WUS.
The B class gene PI encodes a MADS box TF that specifies the petal and stamen identities in specific flower primordia (Fernandez et al., 2013; Goto and Meyerowitz, 1994; Krizek and Meyerowitz, 1996). Plants lacking PI have sepals instead of petals, as well as carpels, filaments, or carpelloid organs instead of stamens (Bowman et al., 1989, 1991; Hill and Lord, 1989). The carnation flower-specific MADS box gene CMB2 is expressed in petals from the initial stages of development until flowers bloom (Baudinette et al., 2000). A previous study indicated CMB2 is most homologous to DEFA, which is a petal and stamen identity-related gene (Baudinette et al., 2000). In Antirrhinum, DEFA (Schwarz-Sommer et al., 1992; Sommer et al., 1990) and GLOBOSA (GLO) are flower-specific homeotic genes that control petal organogenesis in the second whorl, as well as stamen organogenesis in the third whorl (Tröbner et al., 1992; Zachgo et al., 1995). The absence of a functional DEFA results in sepaloids instead of petals and carpelloids instead of stamens (Zachgo et al., 1995). In this context, the downregulated expression of PI (Fig. 3C) and CMB2 (Fig. 3D) may be associated with flower malformations. Moreover, the F-box protein UFO is required for activation of APETALA3 (AP3; encoded by a B class gene), which is responsible for specifying petal and stamen identities (Lee et al., 1997; Sharma et al., 2019). Specifically, UFO functions as a co-factor and mediates the interaction between LEAFY (LFY) and AP3 (Lee et al., 1997), suggesting the observed decrease in UFO expression (Fig. 3E) destabilizes floral development-related activities.
In an earlier study involving mlf plants, most of the flowers were malformed at 15°C, whereas malformed flowers were undetectable at 26°C (Yamane et al., 2018). Interestingly, the flowers of the Antirrhinum DEFA mutant are normal at 15°C, but abnormal at 25–26°C (Zachgo et al., 1995), which is inconsistent with the flower types of the mlf plants examined in this study. In the temperature-sensitive Antirrhinum DEFA mutant, the ability of DEFA to bind to DNA and the formation of the DEFA-GLO protein complex decrease at high temperatures in vitro, which may explain the temperature-dependent lack of a functional DEFA under warm conditions (Zachgo et al., 1995). In the current study, GLO was similarly expressed in WT and mlf (data not shown). In mlf, the expression of CMB2, which is a DEFA homolog, decreased significantly at 15°C (Fig. 3D), suggesting that the CMB2-GLO complex was insufficient for normal petal and stamen development under low-temperature conditions.
Involvement of histone and heat stress-related genes
The expression levels of genes encoding different histones (H2A, H2AX, H2B, H3.2, and H4) were downregulated in mlf at 15°C (Figs. 3F and 4; Table S4). Histone variants have specialized regulatory effects on chromatin dynamics (Kamakaka and Biggins, 2005). Two histone H3 variants (H3.1 and H3.2) differ in terms of the abundance of silencing and activating marks, which is in accordance with the findings of a previous study that indicated that replication-independent histone H3 is enriched in active chromatin (Johnson et al., 2004). A decrease in Histone H3.2 expression may affect the expression of other genes because of the associated changes to methylation sites. Histone acetylation/deacetylation is a dynamic process regulating gene expression (Servet et al., 2010). In Arabidopsis, AtGCN5 (general control non-depressible 5), which is related to histone acetyltransferases, is required for the regulation of floral meristem activities through the WUS/AG pathway (Bertrand et al., 2003). Additionally, RVE8 along with LCL5 (LHY-CCA1-LIKE5) maintain circadian changes in the histone 3 acetylation of the TOC1 (Timing of CAB expression 1) promoter (Farinas and Mas, 2011a, b). Therefore, it is possible that RVE8 expression increases in response to the limited expression of Histone H3. Decreases in the expression of genes encoding diverse histones may affect flower bud morphology.
The HSC70 and TIL expression levels were upregulated in mlf (Fig. 4A, B). The TIL proteins, which belong to the lipocalin family, stabilize the plasma membrane under cold and oxidative stress conditions (Charron et al., 2005; He et al., 2015). The transcription of TIL is induced more in cold-tolerant varieties than in cold-sensitive varieties (Charron et al., 2005). In the present study, TIL expression was upregulated in mlf (Fig. 4B), suggesting the effects of temperature or oxidative stress in the lower buds of mlf.
In carnation leaves, the expression levels of 18 genes encoding 70-kDa heat shock protein (HSP70) forms increase in response to heat stress (42 and 46°C) (Wan et al., 2015). These genes are reportedly involved in protein folding and unfolding, thereby mediating cellular thermotolerance (Hartl, 1996; Su and Li, 2008; Sung and Guy, 2003). Among the HSP family members, HSP70 and its constitutively expressed cognate (HSC70) are abundantly expressed, with the encoded proteins involved in a number of cellular activities (Shi and Thomas, 1992). For example, HSC70 is a constitutively expressed chaperone that contributes to diverse cellular processes (e.g., protein folding and protein degradation). In spinach plants acclimating to cold stress, HSC70s may be important for maintaining cellular homeostasis and suitable protein biogenesis (Anderson et al., 1994). In potato, the HSC70 expression level considerably influences the yield following a moderate increase in the temperature (Trapero-Mozos et al., 2018). Upregulation of HSC70 may have important functions in responses to moderate changes in temperature in mlf. However, the high relative expression of HSC70 was not validated by qPCR (Fig. S6). In addition, SNPs and InDels were detected in the upstream regions and exons of HSC70-1 (Fig. S4) and HSC70-2 (Fig. S5). Because many genes, including 22 HSC70-1 and eight HSC70-2 genes, were detected in the carnation transcriptome, additional research is necessary to determine whether HSC70-1 and HSC70-2 have specific roles related to the development of floral malformations.
Roles of cytokinin-related genes
In this study, the expression levels of cytokinin-activating genes, including LOG1, 4, and 5, and cytokinin-signal transduction genes, including ARR12 and RAV1 were upregulated in mlf (Fig. 5G–I; Table S4). In ‘Baccara’ rose plants, the cytokinin activity was higher in bullhead flowers at relatively low temperatures than in ‘normal’ flowers at higher temperatures (Zieslin et al., 1979). Treating plants grown at high temperatures with cytokinin leads to proliferation of the nectary and promotes the development of adventitious florets, which are characteristic features of malformed flowers (Zieslin et al., 1979). The treatment of Petunia hybrida flower buds with forchlorfenuron (CPPU), which is an inhibitor of cytokinin oxidase/dehydrogenase, results in an enlarged corolla due to an increase in the number of cells (Nishijima et al., 2006). In torenia, a paracorolla and a serrated petal margin are induced by localized extensive cytokinin signaling in flower buds (Niki et al., 2013). Exogenous CPPU treatment caused WT ‘Cherie’ plants to produce flat flowers (data not shown). In an earlier study, benzylaminopurine treatment during the flower development stage resulted in a decrease and subsequent recovery of CLAVATA (CLV) 1 expression and a concurrent increase in WUS expression (Lindsay et al., 2006). Furthermore, WUS and cytokinin signals interact through multiple positive feedback loops that ultimately control the number of stem cells in the SAM (Bartrina et al., 2011; Chickarmane et al., 2012; Gordon et al., 2009). The interplay among cytokinin signals, the WUS/CLV feedback loop, and boundary signals may account for the localization of WUS expression (Adibi et al., 2016; Sablowski, 2009). Cytokinin regulates downstream responses through the action of a multistep phosphorelay that culminates in transcriptional regulation by ARR1, ARR10, and ARR12 (Argyros et al., 2008). RAV1 enhances cytokinin action during primary root development in Arabidopsis (Mandal et al., 2023). Thus, upregulation of ARR12 and RAV1 may increase cytokinin sensitivity and contribute to malformation in mlf (Fig. 5H, I; Table S4). However, there was no significant difference in the trans-zeatin contents of the mlf and WT flower primordia at 15 and 20°C (Fig. S7). Because we only measured the total trans-zeatin content in flower primordia, the localization of trans-zeatin and other cytokinins in the floral meristem remains to be determined.
Roles of abscisic acid-related genes
The NCED1, PP2C, PYL8 and 9, DI19-3, NAC002, and WRKY57 expression levels were upregulated in mlf (Figs. 5A–F, S1, and S2). Notably, NCED is an enzyme that mediates a key step in the ABA biosynthetic pathway in plastids, thereby controlling ABA accumulation in the cell and inducing the mechanism underlying the cellular response to water deficiency (Hwang et al., 2010; Iuchi et al., 2001; Qin and Zeevaart, 1999). Previous research characterized PYR/PYLs (also known as RCARs) as ABA receptors that control ABA signaling by inhibiting PP2C phosphatase activity (Fujii et al., 2009; Ma et al., 2009). Abscisic acid is involved in several stress responses (Shinozaki et al., 2003). In addition, DI19-3, NAC002, and WRKY57 are stress-responsive genes (Jensen et al., 2008; Qin et al., 2014). Previous studies showed AtDI19-3 may participate in plant responses to drought and salt stresses in an ABA-dependent manner (Qin et al., 2014), whereas NAC002 negatively regulates the ABA-activated signaling pathway (Jensen et al., 2008) and WRKY57 modulates drought stress responses (Jiang et al., 2012). Although ABA confers stress tolerance (Iuchi et al., 2001; Shinozaki et al., 2003) and promotes senescence (Trivellini et al., 2010), there are few reports on its effects on floral organ differentiation. In Arabidopsis, drought stress results in the formation of underdeveloped filaments that cannot reach the top of the pistil (Trivellini et al., 2010). Based on these findings, it is possible that the expression of ABA- and stress-related genes is induced by exposure to stress conditions, including water stress, in mlf flower buds.
Roles of calmodulin-related genes
The CAM expression level was upregulated in mlf (Fig. S2; Table S4). Sequence variants were detected in CAM and Calmodulin-binding transcription activator 5 (Table S5). Calmodulin, which is a Ca2+-binding protein, plays important roles in the Ca2+ signaling pathway (Ranty et al., 2006; Snedden and Fromm, 1998). The importance of Ca2+ for mediating plant responses to external abiotic stressors, including light, cold, heat, and drought, has been reported (Snedden and Fromm, 1998). In this context, the expression of CAM, as well as ABA- and stress-responsive genes, may be upregulated by heat or drought stress in mlf. In wheat, the expression of Wheat Calmodulin-Binding Protein 1 (WCBP1) is significantly upregulated in the young spikes of plants exhibiting pistillody (Saraike et al., 2007; Yamamoto et al., 2013). Thus, the possible effects of CAM on flower malformations will need to be examined more thoroughly.
Roles of ribosomal proteins
The GO terms ribosome (CC category) and structural constituent of ribosome (MF category) were enriched in WT (Fig. 2B, C). A comparison of FPKM values indicated several ribosomal protein-encoding genes had downregulated expression levels in mlf (Table S4). Ribosome sequence variations were detected (Table S5). The ribosome is an essential ribonucleoprotein complex that facilitates translation and is indispensable for plant growth (Horiguchi et al., 2012). In Arabidopsis, a pointed first leaf (pfl)1-like mutation in Small ribosomal subunit protein bS18c (RPS18) results in multiple phenotypic changes (e.g., aberrant leaf and trichome morphological features, retarded root growth, and delayed flowering) (Ito et al., 2000). A mutation in Large ribosomal subunit protein uL13x (RPL13aC) leads to a decrease in potassium accumulation in the shoots and roots, resulting in morphologically abnormal roots (Ma et al., 2023). Mutations in the genes encoding ribosomal proteins or decreases in the expression of these genes in mlf flower buds may result in ribosome deficiencies, heterogeneity, or abnormalities, with adverse implications for morphogenesis.
Relationship between the super-double-flower phenotype and floral malformations
The appearance of the malformed mlf flowers at 15°C was consistent with previously reported super-double-flower phenotypes (Imai, 1938; Saunders, 1917; Yagi et al., 2014a). The F1 progeny derived from the hybridizations of ‘Cherie’ (WT and mlf) with ‘Fossett rose’ and ‘Buryule’ produced several malformed flowers (Fig. S8), similar to the F2 progeny (DD genotype) examined in an earlier investigation (Yagi et al., 2014a). According to the model of Saunders (1917), the estimated normal-to-malformed flowers ratio for the F1 progeny is 3:1 for WT, and that for mlf is 1:1; however, neither separation ratio was applicable (Table S7). Imai (1938) proposed that the dominant allele Da in a heterozygous genotype (Dada) results in the standard double-flower phenotype, but the homozygous genotype (DaDa) results in the bullhead-double-flower (super-double-flower) phenotype. The other dominant allele (Db) associated with the double-flower phenotype has a weaker effect than Da (i.e., the homozygous genotype only results in the standard double-flower phenotype) (Imai, 1938). These findings are supported by the relatively high proportion of super-double-flowers among the F1 progeny following the crosses with normal double-flower cultivars similar to ‘Spectrum’ (Imai, 1938). If the mlf and WT ‘Cherie’ plants have the DaDaDbdb and DadaDbdb genotypes, the proportions of normal and malformed flowers among the F1 progeny of mlf and WT fit the 3:1 and 7:1 ratios, respectively (Table S7). To date, temperature-dependent phenotypic variations have not been detected in the F1 progeny, suggesting that super-double-flowers and temperature dependence are phenotypes that are inherited independently. However, analyzing malformed super-double-flower phenotypes only based on appearance is not ideal because the super-double-flower phenotype tends to vary depending on the season and flower position (Yagi et al., 2014a).
An RFLP marker was previously developed to discriminate between the semi-double-flower and double-flower phenotypes in carnations (Mediterranean and American groups) (Scovel et al., 1998). Moreover, RAPD and SSR markers applicable for carnation were developed (Onozaki et al., 2006; Yagi et al., 2011), among which CES1982 and CES0212 were closely linked to the D85 locus, which controls the flower type (single or double) (Yagi et al., 2014a). Loci markers linked to double-flower and single-flower phenotypes were designated as CES1982 (176 and 180 bp) and CES0212 (269 and 264 bp), respectively (Yagi et al., 2014b). Although it is possible that there are substantial differences in the marker loci between the cut-flower cultivars (standard and spray types) (Yagi et al., 2014a, b) and the potted carnation cultivar ‘Cherie’, both WT and mlf genomes had the same CES1982 (151 bp) and CES0212 (265 bp) loci (Table S8), suggesting that WT and mlf plants contain the DD genotype. Yagi et al. (2014b) suggested that carnation plants in the Mediterranean group (DD) have a mechanism that suppresses the development of malformed flower organs. This mechanism is likely inactivated in mlf plants grown at approximately 15°C, thereby causing super-double-flowers to form. Somaclonal variations associated with the super-double-flower phenotype probably occurred during the propagation of mlf lines, but the variants were masked under warm conditions, with the malformed flower phenotypes occurring at temperatures below 15°C in production fields. Because histone acetylation/deacetylation is a dynamic gene regulatory process (Servet et al., 2010), its potential contribution to floral malformations should be investigated in future studies.
Based on previous studies, the miR172-mediated repression of AP2 expression is critical for the timely termination of floral stem cells and the maintenance of floral meristem size (Chen, 2004; Jung et al., 2014; Zhao et al., 2007). The overexpression of miR172-resistant AP2 in Arabidopsis increases AP2 protein levels and results in an indeterminate floral meristem that produces an excessive number of stamens or petals (Chen, 2004; Wollmann et al., 2010). Further research is necessary to clarify the relevant roles of microRNAs.
Conclusion
The upregulated expression of WUS, SAP, and cytokinin-activating genes and the downregulated expression of CMB2, PI, and UFO may destabilize the floral meristem and induce the formation of numerous undifferentiated floral organs. The downregulation of histone-encoding gene expression levels and the upregulation of HSC70 expression may lead to decreased transcription and contribute to the temperature-dependent conversion of phenotypes. The lack of uniformity in the floral malformations in mlf plants is a limitation of this study; therefore, elucidating the involvement of specific genes was difficult. Furthermore, the final phenotypes of the flower buds from which RNA was extracted are unclear. Although it is necessary to consider the influence of the genotypes of the two varieties used as pollen parents, the malformed flower phenotypes are inherited and may be related to the super-double-flower phenomenon. However, further research is needed to identify the trigger for the malformed flower phenotypes and to characterize the mechanism underlying the double-flower and super-double-flower phenotypes in carnation.
