Anthocyanin color patterns are often observed in flowers, including spots, stripes, venation, bud-blush, and bicoloration (two colors in a single petal). The color patterns create major variations in flower appearance and often enhance the ornamental value of flowers. Consequently, the exploration of the underlying mechanisms involved in the development of anthocyanin color patterns could facilitate the breeding of floricultural crops (Nakayama, 2014). Mechanisms of development of anthocyanin color patterns have been reported in some plant species, mainly model plants (Davies et al., 2012; Glover et al., 2013; Albert et al., 2014). The primary mechanism is transcriptional regulation of anthocyanin biosynthesis genes. Subgroup 6 members of R2R3-MYB transcription factors in numerous species and subgroup 5 members of R2R3-MYB transcription factors in orchids predominantly regulate anthocyanin biosynthesis in flowers (the grouping of R2R3-MYBs is according to Stracke et al., 2001), and single plant species often have a couple of R2R3-MYB genes grouped into the subgroups. Each of the genes exhibits spatially and temporally distinct expression profiles and causes restricted pigment deposition, resulting in the development of varied color patterns. For example, petunias have four subgroup 6 R2R3-MYB genes: PhAN2 paints the entire petal region, PhAN4 causes pigmentation in the flower tubes, and DEEP PURPLE and PURPLE HAZE are responsible for venation and bud-blush pigmentation patterns, respectively (Albert et al., 2011). In Phalaenopsis, three subgroup 5 R2R3-MYB genes regulate full-red pigmentation, red spots, and venation patterns, respectively, in petals (Hsu et al., 2015).

Lilies (Lilium spp.) are among the most commercially valuable and globally cultivated floriculture crops (Ministry of Agriculture, Forestry and Fisheries of Japan, 2019). Interspecific hybridization is the principal method of lily breeding. The Lilium species are classified into several sections, and Asiatic hybrid lilies (Lilium spp.) are derived from crosses among species in Sinomartagon and Daurolirion section (Yamagishi and Nakatsuka, 2017; Marasek-Ciolakowska et al., 2018). Since lilies are a major floricultural crop, genetic evaluation of their agricultural traits facilitates the improvement of the character traits of cultivars. However, lilies are not model species, and their genetic evaluation is challenging considering Lilium species have a large genome size (33–36 Gb; Bennett and Leitch, 2011) and long life cycles (3–7 years from sowing to anthesis), and its stable transformation remains challenging.

A high variation in flower color hues is the major selling point for Asiatic hybrid lilies, which accumulate anthocyanins (cyanidin 3-rutinoside, Nørbæk and Kondo, 1999) and carotenoids (Wang and Yamagishi, 2019) in flowers. In addition, a variety of anthocyanin color patterns, such as restricted pigmentation at tepal bases (large spots or blotches) or the upper halves (bicolor), and several types of spots, including raised spots, splatter-type spots, and brush mark spots, are often observed in Asiatic hybrid lily cultivars (Yamagishi, 2013). The underlying mechanisms of such color patterns have been evaluated. Several subgroup 6 members of the R2R3-MYB genes are expressed in lily flowers: MYB12 usually paints entire tepal regions (Yamagishi et al., 2010), MYBSPLATTER (formerly LhMYB12-Latvia, Yamagishi, 2020a) causes splatter-type anthocyanin spots (Yamagishi et al., 2014), MYB19Short and MYB19Long are responsible for raised-spot pigmentation (Yamagishi, 2020b), MYB19Long is indispensable for brush mark pigmentation (Yamagishi, 2020c), and MYB18 is involved in the development of large spots (blotches) at tepal bases (Yamagishi, 2018). Therefore, transcriptional regulation, governed by subgroup 6 members of R2R3-MYB genes, is the primary mechanism of the establishment of a variety of color patterns in lilies.

Some Asiatic hybrid lily cultivars, including Lollypop, develop bicolor tepals where anthocyanin accumulates in the upper halves of tepals but scarcely accumulates at the lower halves of tepals (Figure 1). We have previously compared transcriptome profiles between the upper and lower halves of Lollypop tepals and clarified that the entire anthocyanin biosynthesis pathway is suppressed at the lower halves (Suzuki et al., 2016). Therefore, the bicolor pattern in Lollypop is caused by transcriptional regulation of anthocyanin biosynthesis genes. Putative transcription factors that positively or negatively regulate anthocyanin biosynthesis are surveyed using transcriptome data, among which, certain negatively regulating transcription factors are expressed in Lollypop tepals, including R3-MYB suppressors (Sakai et al., 2019). If such negative regulators are involved in bicolor formation, their expression would be expected to be high in white tepal regions. To date, however, there have been no reports of putative suppressors exhibiting higher expression in the unpigmented regions of flowers (Suzuki et al., 2016; Sakai et al., 2019). Three R2R3-MYB sequences belonging to subgroup 6 are included in the Lollypop transcriptome. Among the three R2R3-MYB sequences, MYB12 is expressed mainly in the tepals, and its expression levels are higher in the upper halves than in the lower halves of the tepals, suggesting that MYB12 is the main transcription factor causing varying anthocyanin biosynthesis activities, and mediating bicolor anthocyanin pigmentation patterns (Suzuki et al., 2016). However, since MYB12 is often responsible for the pigmentation of whole tepal regions in many lilies (Yamagishi et al., 2010, 2012; Yamagishi, 2011), it remains unclear why and how MYB12 expression is suppressed in the lower halves of bicolor tepals.

MicroRNAs are non-coding RNA species that play key roles in plant development and stress responses. MicroRNAs usually suppress the activity of target genes post-transcriptionally through cleavage of the transcripts or translational attenuation. Many transcription factor genes, including MYB family genes, are targets of microRNAs. MicroRNA production starts from primary microRNA (pri-miRNA), which is cleaved in the nucleus into smaller precursor microRNA (pre-miRNA) with a hairpin structure. These pre-miRNAs are then exported to the cytoplasm and cleaved by DICER-LIKE 1 (DCL1) to form 20–24 nucleotide-long mature microRNAs. Subsequently, the mature microRNAs are loaded into the ARGONAUTE (AGO) effector protein, a part of RNA-induced silencing complexes (RISCs). Afterward, the mature microRNAs guide the RISCs to target mRNA(s) through sequence complementarity and eventually mediate gene silencing (Zhang, 2015; D’Ario et al., 2017). In addition, some 22 nucleotides-long microRNAs are capable of initiating secondary small interfering RNAs (siRNAs). The RNA fragments downstream (or upstream in rare cases) of the cleaved site are converted into a double-stranded RNA molecule by RNA-dependent RNA polymerase 6, which is processed by the DCL4 enzyme into approximately 21 nucleotide-long siRNAs, which are called phased siRNAs (phasiRNAs). Some phasiRNAs can further regulate the expression of other target genes in trans and are called trans-acting siRNAs (tasiRNAs; Shivaprasad et al., 2012; Fei et al., 2013; Deng et al., 2018).

Among the microRNAs identified in a wide range of plant species, miR828, miR858, and miR159 are capable of targeting R2R3-MYB genes (Li and Lu, 2014; Jia et al., 2015; Sharma et al., 2016). Among them, miR828 often suppresses the expression of R2R3-MYB genes belonging to subgroup 15 and subgroup 4, which negatively regulate anthocyanin accumulation (Xia et al., 2012; Zhu et al., 2012). Therefore, anthocyanin pigment amounts increase in fruits or tubers with an increase in miR828 accumulation (Bonar et al., 2018; Tirumalai et al., 2019). However, in Arabidopsis, miR828 directly targets AtMYB113 and indirectly targets AtMYB75 (AtPAP1) and AtMYB113 through TAS4-siR81(-), one of the tasiRNAs derived from Trans-Acting SiRNA Gene 4 (TAS4) transcripts (Rajagopalan et al., 2006). Since AtMYB75 and AtMYB113 are subgroup 6 members of R2R3-MYB that positively regulate anthocyanin biosynthesis, miR828 and TAS4 suppress anthocyanin biosynthesis and are involved in feedback regulation of anthocyanin accumulation in vegetative tissues under stress conditions (Luo et al., 2012).

An miR828 target sequence was found in the nucleotide sequence of lily MYB12 (shown in the results section); however, sequences potentially recognized by miR159, miR858, or TAS4-siR81(-) were not detected. In addition, partial sequences of pri-MIR828 were included in the Lollypop transcriptome (Suzuki et al., 2016). Therefore, the potential of miR828 to suppress MYB12 expression predominantly in the white tepal region was investigated in the present study. We demonstrated that lily miR828 was highly accumulated in the white tepal region and cleaved MYB12 mRNA. MYB12 transcript-derived secondary siRNA molecules investigated using microRNA-seq were also highly accumulated in the white region. Consequently, we conclude that the miR828/MYB12 module has a key role in the development of bicolor tepals in lilies.

Bicolor petal patterns, which consist of pigmented and unpigmented regions in a single petal, are often observed in higher plants. Several mechanisms that create bicolor patterns have been elucidated in some species, mainly petunias. One is a post-transcriptional silencing of CHS genes examined in petunias (Koseki et al., 2005) and dahlias (Ohno et al., 2011). In star and picotee petunia lines, siRNA suppresses CHS mRNA accumulation post-transcriptionally in white petal regions, and tandemly arranged CHS sequences are involved in siRNA generation (Morita et al., 2012). This floral color pattern is further modified by the activation of an endogenous virus. Petunia vein clearing caulimovirus (PVCV) proviruses are activated as host plants age, and CHS gene silencing is suppressed by an RNA silencing suppressor of activated PVCV, resulting in anthocyanin biosynthesis in the white regions of star-type petunia petals (Kuriyama et al., 2020). Another mechanism is s space-specific flavonol synthase (FLS) expression. In morn petunia lines, the central regions of flowers are colorless. FLS levels are high, and flavonols accumulate at much higher levels than anthocyanins in the central region (Saito et al., 2006). In Lollypop tepals, the entire anthocyanin biosynthesis pathway was suppressed in the white region, and the key factor that caused the suppression was the MYB12 transcription factor. In addition, most Asiatic hybrid lily cultivars do not accumulate flavonols and flavones in tepals (Lai et al., 2012), and FLS was scarcely expressed in Lollypop tepals (Suzuki et al., 2016). Therefore, in the present study, we explored additional mechanisms.

Lily MYB12 is a target of miR828 because miR828 cleaves MYB12 mRNA in lily tepals, and MIR828 introduced into tobacco leaves cleaved MYB12 mRNA and suppressed MYB12 activity. In addition, mature miR828 and the MYB12-derived secondary siRNA molecules accumulated more abundantly in the white region than in the pigmented region. Therefore, MYB12 expression was significantly lower in the white tepal regions due to targeting by miR828 and potentially by the secondary siRNAs. Other microRNAs that potentially hybridized into MYB12 mRNA were surveyed using the microRNA-seq data, but such a candidate was not found among the small RNAs of 20–24 nucleotides (data not shown). Therefore, miR828 is potentially a solo microRNA suppressing MYB12 expression, and the miR828/MYB12 module is potentially the key mechanism causing discoloration in lily flowers.

Why does mature miR828 predominantly accumulate in the lower halves of tepals? The amounts of pri-MIR828 transcripts were high in the lower halves and, therefore, spatially variable transcription of pri-MIR828 would be one of the primary mechanisms. However, accumulation profiles of mature miR828 and pri-MIR828 transcripts were different: the amounts of mature miR828 increased during flower bud growth, whereas those of pre-MIR828 decreased with flower bud development. Therefore, processing and stability of microRNA were also potentially involved in the predominant accumulation of mature miR828 in the lower halves, and the post-transcriptional regulation of miR828 accumulation unlikely occurs in the upper halves because the amounts do not increase with bud development. MicroRNA abundance is regulated under multiple levels of control, including processing, RISC assembly, and stability (Wang et al., 2019). Different AGO proteins have different effects on the stability of their bound microRNAs, and spatially and temporally distinct AGO5 expression is indispensable for organ-specific function of miR156 (Roussin-Léveillée et al., 2020).

MYB12 has been found to be the primary transcription factor regulating bicoloration (Suzuki et al., 2016), and we experimentally confirmed the predominant cleavage of MYB12 mRNA in the lower halves of tepals. However, compared with the upper halves, the levels of accumulated MYB12 transcript in tepals were observed to be reduced by 40% in the lower halves, in which the expression of biosynthesis genes was markedly suppressed (Figure 4; Suzuki et al., 2016). We suspect that the observed differences in suppression levels could be attributable to the fact that miR828 not only cleaves the transcripts but also arrests MYB12 translation in lily cells. In addition to the cleavage of target mRNA, translational attenuation of target genes is involved in microRNA-mediated gene regulation (Dalmay, 2013), and it is often reported that suppression of the accumulation of target gene transcripts is weaker than that of downstream genes regulated by the target genes (Gandikota et al., 2007; Wang et al., 2016; Gattolin et al., 2018).

Sequences of mature miR828 are highly conserved in many plants, including lilies, with an exception in Arabidopsis, in which one nucleotide is different. Major targets of miR828 are R2R3-MYB transcription factor genes and TAS4: MdWEREWOLF (WER, subgroup 15 member of R2R3-MYB), MdC1 (R2R3-MYB), and MdTAS4 in apple (Xia et al., 2012), AtMYB113 (subgroup 6 member of R2R3-MYB), AtMYB82 (ungrouped R2R3-MYB), and AtTAS4 in Arabidopsis (Rajagopalan et al., 2006), VvMYB114 (WER like MYB having EAR suppressor motif) and VvWER in Vitis vinifera (Tirumalai et al., 2019), PpMYB (subgroup 15 member of R2R3-MYB), two PpMYB (ungrouped MYB), and PpTAS4 in peach (Zhu et al., 2012), and R2R3-MYB and StTAS4 in potato (Bonar et al., 2018). R3 repeat sequences in R2R3-MYB genes are highly conserved, and the target site of miR828 exists near the end of the R3 repeat. Therefore, miR828 often targets several MYB genes (Xia et al., 2012; Zhu et al., 2012; Li and Lu, 2014). In Lollypop, 15 sequences were predicted to be targets, although whether they were cleaved by miR828 has not yet been evaluated, and five of the 15 sequences were R2R3-MYB (Supplementary Table S2). MYB15-like exhibited homology with LrMYB15 in Lilium regale, which regulates bud blush anthocyanin pigmentation (Yamagishi, 2016), while, in Asiatic hybrid lilies, its expression levels are not high and its function is unknown (Yamagishi, 2018). MYB16 is also expressed weakly in tepals, and its function is unknown (Yamagishi, 2018). Zm1-like Myb-related protein is a subgroup 2 member of R2R3-MYB (Suzuki et al., 2016). Although its function in lily is unknown, subgroup 2 members are not involved in pigmentation in other species. MYB19S is predominantly expressed in stage 2 flower buds and is involved in the pigmentation at raised spots developed at the basal parts of tepals (Yamagishi, 2020b). Therefore, the R2R3-MYB genes, other than MYB12, are unlikely to be involved in bicoloration.

The involvement of miR828 in pigment accumulation has been reported. MiR828 highly accumulates in pigmented grape fruits and potato tubers and inhibits the expression of grape MYB114 and potato R2R3-MYB, which act as suppressors of anthocyanin biosynthesis (Bonar et al., 2018; Tirumalai et al., 2019). In lilies, miR828 accumulates in the unpigmented regions of tepals and suppresses the positive regulation of anthocyanin biosynthesis. Therefore, the phenomenon in lilies is distinct from those observed in grapes and potatoes.

Subgroup 6 members of R2R3-MYB in numerous species and subgroup 5 members of that in Poaceae and orchids are major positive regulators of anthocyanin biosynthesis. To determine whether miR828 targets these positive regulator genes, putative miR828 target sites found in the genes were compared (Supplementary Figure S4). Base-pairing between miR828 and the putative target sites showed that the number of mismatched nucleotides was three or fewer in R2R3-MYB genes of oil palm, onion, Arabidopsis, snapdragon, and Salvia. In Arabidopsis, AtMYB113 is the direct target of miR828 (Luo et al., 2012), and AtPAP1 and AtMYB113 expression is suppressed in miR828-overexpressed Arabidopsis plants (Yang et al., 2013), indicating that miR828 regulates anthocyanin accumulation in Arabidopsis through the suppression of the positive regulators. However, as miR828 predominantly accumulates in siliques (Rajagopalan et al., 2006), miR828 is unlikely to be involved in flower coloration in Arabidopsis. Deep sequencing analysis revealed that SmMYB36, a subgroup 6 R2R3-MYB in Salvia, is a target of miR828 (Li and Lu, 2014), although the function of SmMYB36 has not been clarified yet. Numerous mismatched nucleotides have been observed in the sequences of other eudicot species, including petunia, grape, and apple, suggesting that miR828 does not target such positive regulators directly in these species. Consequently, the bicolor formation mediated by the modules of miR828 and subgroup 6 members of R2R3-MYB, observed in lily flowers, should be rare in eudicots. Conversely, in addition to the lily MYBs, two of the three subgroup 6 members of R2R3-MYB (EgVIR in oil palm and AcMYB1 in onion) in monocots are putative targets of miR828, although the number of subgroup 6 members of R2R3-MYB genes isolated and characterized in monocots is limited. Therefore, further investigation of miR828-mediated pattern formation should provide novel insights in monocots. Subgroup 5 members of the R2R3-MYBs are unlikely targets of miR828.

Several mechanisms of the development of anthocyanin color patterns have been clarified in flowers. The major mechanisms are spatially distinct expression of R2R3-MYB positive regulators in single petals (Fattorini and Glover, 2020), post-transcriptional gene silencing of CHS (Morita et al., 2012), and the insertion of transposable elements into anthocyanin biosynthesis genes (Momose et al., 2013; Uchiyama et al., 2013). MicroRNA-mediated pattern formation, presented in this study, is likely the novel mechanism of development of the attractive color patterns.

These deep sequencing files have been deposited in the DDBJ Sequence Read Archive under the accession number DRA010576.

MY designed the research, conducted the experiments, and wrote the manuscript. MS and MY analyzed the data. Both the authors have read and approved the submission of the manuscript.