Identification of microRNAs and their target genes related to the accumulation of anthocyanin in purple potato tubers (Solanum tuberosum).

MicroRNAs (miRNAs) are types of endogenous non-coding small RNAs found in eukaryotes that are 18-25 nucleotides long. miRNAs are considered to be key regulatory factors of the expression of target mRNA. The roles of miRNAs involved in the regulation of anthocyanin accumulation in pigmented potatoes have not been systematically reported. In this study, the differentially expressed miRNAs and their target genes involved in the accumulation of anthocyanin during different developmental stages in purple potato (Solanum tuberosum L.) were identified using small RNA (sRNA) and degradome sequencing. A total of 275 differentially expressed miRNAs were identified in the sRNA libraries. A total of 69,387,200 raw reads were obtained from three degradome libraries. The anthocyanin responsive miRNA-mRNA modules were analyzed, and 37 miRNAs and 23 target genes were obtained. Different miRNAs regulate the key enzymes of anthocyanin synthesis in purple potato. The structural genes included phenylalanine ammonia lyase, chalcone isomerase, flavanone 3-hydroxylase, and anthocyanidin 3-O-glucosyltransferase. The regulatory genes included WD40, MYB, and SPL9. stu-miR172e-5p_L-1R-1, stu-miR828a, stu-miR29b-4-5p, stu-miR8019-5p_L-4R-3, stu-miR396b-5p, stu-miR5303f_L-7R + 2, stu-miR7997a_L-3, stu-miR7997b_L-3, stu-miR7997c_L + 3R-5_2ss2TA3AG, stu-miR156f-5p_L + 1, stu-miR156a, stu-miR156a_R-1, stu-miR156e, stu-miR858, stu-miR5021, stu-miR828 and their target genes were validated by qRT-PCR. They play important roles in the coloration and accumulation of purple potatoes. These results provide new insights into the biosynthesis of anthocyanins in pigmented potatoes.

INTRODUCTION

Potato (Solanum tuberosum L.) is the fourth most widely cultivated food crop in the world. China produces the largest number of potatoes in the world, and Inner Mongolia is an important area for the breeding of seed potatoes and the production of commercial potatoes in China. However, the production of pigmented potatoes lag far behind the market demand for them. Pigmented potatoes are special types of cultivated potatoes. Its skin and flesh are purple, red, pink, and other colors. They not only contain starch, protein, a variety of trace elements, and a small amount of fat, but they are also more highly nutritious because they are rich in anthocyanins, carotenoids, and other antioxidants (Vaitkevičienė, 2019; Yin et al., 2016). A total of 26 anthocyanins were isolated from pigmented potatoes, and 24 were identified. Five cis isomers were identified, and four of them were reported for the first time. They include cis‐petanin, cis‐peonanin, petunidin 3‐cis caffeoylrutinoside‐5‐glucoside, and petunidin 3‐cis‐feruloylrutinoside‐5‐glucoside (Kim et al., 2018). Potatoes with red flesh are unusual in having a significant content of pelargonidin‐3‐feruloylrutinoside‐5‐glucoside, and those with purple flesh contained a significant amount of cyanidin‐3‐rutinoside (Rytel et al., 2019). Anthocyanins have been found to improve human health because of their free radical scavenging and their antioxidant, anticancer, antimicrobial, and antiviral activities (Stone et al., 2007), anti‐aging (Poulsen et al., 2020), and therapeutic effects on a variety of diseases, including dementia (Khalifa et al., 2020) and diabetes (Markovics et al., 2020). Therefore, the study of regulation of anthocyanin metabolism is highly significant to develop excellent gene resources and improve the varieties of colored potatoes.

The genes that participate in anthocyanin biosynthesis include both structural and regulatory genes. The structural genes directly encode an assemblage of enzymes, including phenylalanine ammonia lyase (PAL), chalcone synthase (CHS), chalcone isomerase (CHI), flavonoid 3‐hydroxylase (F3H), dihydroflavonol reductase (DFR), anthocyanin synthase (ANS), and anthocyanin 3‐O‐glucosyltransferase (UFGT) (Hara et al., 2003; Li et al., 2014). The genes for these enzymes can be regulated by the transcription of regulatory genes, such as WRKY, bHLH (Arlotta et al., 2020), MYB (Khan et al., 2022; Li et al., 2020), WD40 (Wang et al., 2020), and HY5 (Bustamante et al., 2018). In addition, some regulatory genes have also been reported to be related to the biosynthesis of anthocyanin, including squamosa promoter‐binding protein‐like (SPL) (Yang et al., 2021), auxin response factor (ARF) (Wang et al., 2020), and Jasmonate Zim‐domain (JAZ) (Adrian et al., 2020). In addition to the involvement of structural and regulatory genes, recent studies show that microRNAs (miRNAs) can also play important roles in mediating the biosynthesis of anthocyanin in plants.

MicroRNAs (miRNAs) are endogenous non‐coding small RNAs (sRNAs) in eukaryotes that are 18–25 nucleotides long (Marc & Filipowicz, 2010). In plants, the mature microRNA binds the Argonaute (AGO) protein and cleaves it to target mRNA. Moreover, it cuts the nucleotides at position 10–11 at the 5′ end to the 3′ end of the target gene. Once the mRNA is degraded into fragments, the encoded protein is no longer functional (Iwakawa & Tomari, 2015). miRNAs are key regulatory factors that target the expression of mRNA. The MiRBase database (mirbasev22.0) contains 38,589 miRNA precursors and 48,860 mature miRNAs from 271 species. Many of them are widely involved in the plant growth and development, metabolism, hormone regulation, and various biological and abiotic stresses (Li et al., 2020). Presently, miR828, miR858, miR156, miR159, miR157, miR160, miR172, miR393, miR396, miR824, and miR870 in Litchi chinensis (Liu et al., 2017), Punica granatum (Saminathan et al., 2016), Hylocereus monacanthus (Chen et al., 2020), Raphanus sativus (Gao et al., 2019), Lycium ruthenicum (Qi et al., 2020), Diospyros kaki (Luo et al., 2015), Vaccinium ashei (Li et al., 2018), Actinidia arguta (Li et al., 2019), and Lonicera edulis (Li et al., 2018), respectively, have been found to be involved in the regulation of anthocyanin metabolism. Both miR156h‐3p and miR396a‐3p can target UFGT in L. edulis and significantly negatively correlate with anthocyanin biosynthesis (Cui et al., 2020). The upregulation of miR156 leads to the downregulation of LcSPL1, which negatively regulates anthocyanin biosynthesis by interacting with LcMYB1 in litchi (L. chinensis) (Liu et al., 2017). Light‐induced downregulation of BrmiR828 can target BrTAS4, BrPAP1 (Bra039763), and MYB82 (Bra022602) to negatively regulate their transcript levels leading to the accumulation of MYB transcription factors that positively regulate anthocyanin biosynthesis in light‐exposed seedlings of turnip (Brassica rapa) (Zhou et al., 2020). miR828 exists in pineapple (Ananas comosus) and directs post‐transcriptional gene silencing of mRNAs encoding MYB family members with inferred function to regulate the conspicuous red fruit trait in var (Christopher, 2020). In carmine radish (R. sativus), miR165a‐5p, miR172b, miR827a, miR166g, and miR1432–5p participate in the biosynthesis of anthocyanins (Gao et al., 2019). miR828, miR858, miR165/166, and miR156 were identified to be involved in the biosynthesis of anthocyanins in tomato (Solanum lycopersicum) in the Solanaceae family (Jia et al., 2015). The biological function of miRNA is closely related to the function of its target genes. Therefore, it is important to explore and identify miRNAs and its target genes to elucidate the complex regulatory mechanism of anthocyanin synthesis mediated by miRNA. However, the roles of miRNAs involved in the regulation of the biosynthesis of anthocyanin in pigmented potatoes have not been systematically reported. To our knowledge, this is the first report on the generation of the miRNA and respective targets that are potentially involved in the biosynthesis of anthocyanin in purple potato. Small RNA and degradome libraries from the tubers of purple potato S1 (tuber formation stage), S2 (tuber bulking stage), and S3 (tuber maturation stage) were constructed and sequenced. Our goal was to use an integrated analysis of GO/KEGG/qRT‐PCR to comprehensively explore the stu‐miRNA populations in the biosynthesis of anthocyanins in pigmented potatoes.

MATERIALS AND METHODS

Plant materials and sample preparation

A tetraploid potato genotype with purple tuber skin and flesh, Huasong 66, which was used in this study was kindly provided by HuaSong Seed Industry Co., Ltd. (Beijing, China). Seeds were sown on the farm of Agricultural College of Inner Mongolia Agricultural University (Hohhot, China) (40°46`N, 110°45`E) in early May 2018. The potato plants were weeded, fertilized, and watered during their developmental phases. The potato tubers of five different developmental phases (Sa, Sb, Sc, Sd, and Se) were collected between July 29 and September 7, 2018 at 10‐day intervals. When all the stems and leaves had withered, the final samples (Sf) were collected on September 26, 2018. According to the classification standards of potato growth stages (Liu et al., 2020), the sample of tubers was mixed at the tuber formation stage (S1, containing Sa and Sb), tuber bulking stage (S2, containing Sc, Sd, and Se), and tuber maturation stage (S3, Sf), respectively.

Three tubers from the same plant were used as three biological replicates. Anthocyanin was extracted from some of the fresh samples, while other samples were immediately frozen in liquid nitrogen and stored at −80°C until use.

Analysis of the anthocyanin content

The anthocyanins were extracted from purple potato tubers using a modification of the method described by Huan et al. (2020). One‐gram samples of fresh purple potato tubers were collected in triplicate at each stage (Sa–Sf). The samples were ground and extracted twice by adding 15 mL of leaching solution (a 1:1 mixture of .1 mol/L HCl: 95% ethyl alcohol [v/v]) at 60 °C for 1.5 h. The two extracts were mixed and collected by centrifugation at 12,000 × g for 8 min. Each solution was extracted with the buffers that were used in the leaching solution at pH 1.0 and pH 4.5 and measured using visible spectroscopy (TU‐1810SPC; Beijing, China) at 530 nm and 700 nm. The content of anthocyanins was calculated using a pH differential method.

Small RNA, degradome sequencing

The total RNA was extracted from the tubers using a TRK‐1002 Total RNA Purification Kit (TPK‐1002; LC Sciences, Hangzhou, China) according to the manufacturer's instructions. Equal volumes of RNA extracts from Sa to Sf were pooled for use as S1 or S3 to construct and sequence the sRNA library. Three libraries were denoted S1 (Sa and Sb), S2 (Sc, Sd, and Se), and S3 (Sf).

Analysis of the small RNA data

The sRNA library was constructed using TruSeq Small RNA Sample Prep Kits (Illumina, San Diego, CA, USA), and the libraries were sequenced on an Illumina HiSeq 2,500 platform to prepare 1 × 50 bp single end reads. The ACGT101‐miR software that was used to analyze the miRNA data was obtained from Lianchuan Bio‐Technology headquartered as LC Sciences (Houston, TX, USA). The analytical process was as follows: clean reads were obtained after quality control treatment of the original data, and the 3′ linker was removed. The reads were screened for length, so that the bases were 18–25 nt long. The remaining sequences were compared with those in various RNA databases (excluding miRNA) that included rRNA, tRNA, small nuclear RNA (snRNA), and small nucleolar RNA (snoRNA) among others. These databases included the mRNA and Rfam databases and the Repbase database (repetitive sequence database). The sequences were then filtered, and the final data obtained were valid and could subsequently be used to analyze the small RNA data. Suitable difference test methods were used to screen significant genes based on the different experimental designs. The differentially expressed genes were screened based on P ≤ 0.05. Computational target prediction algorithms, such as TargetFinder, were used to identify the miRNA‐binding sites to determine the target genes for the most abundant miRNAs. The Gene Ontology (GO) terms and the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway of the most abundant miRNA targets were also annotated.

The degradome data process and target identification

The degradome library was constructed as previously described by Ma and Axtell (2010) with some modifications. Magnetic beads were used to capture the mRNA, which were connected using a 3.5′ adaptor. The degradome libraries were sequenced using an Illumina HiSeq 2,500 platform for 1 X 50 bp single end reads. These libraries were constructed at Lianchuan Bio‐Technology (Hangzhou, China).

The software to analyze the degradome data was provided by the Lianchuan Biology Institute (ACGT101‐DEG; LC Sciences). A series of data processing steps enabled the use of original data obtained by sequencing for the subsequent analyses of comparable sequences. The comparative sequences were compared with those from the cDNA database of sequenced species to generate the degradome density file. CleaveLand version 3.0 (Quaye et al., 2009) was used to predict the target gene mRNA sequence matched with the small RNA sequence of the species sequenced. The target genes that corresponded to the predicted miRNA and the mRNAs in the degradome density file were combined to identify the common mRNAs that were the target genes of miRNA. T‐plots were built to analyze the target genes for miRNA and the patterns of RNA degradation. The GO_ID of anthocyanin biosynthesis was determined by the description of the GO_term in GO database. The pathway_ ID of anthocyanin biosynthesis was screened through the path description information in the KEGG database. The calculating formula of P‐value (p < .05) is as follows: Here N is the number of all genes with GO/KEGG annotation; n is the number of target gene in N; M is the number of all genes that are annotated to the certain GO terms/pathways; and m is the number of DEGs in M. N stands for total background gene (TB gene number); n stands for total significant gene (TS gene number); M stands for background gene (B gene number); and m stands for significant gene (S gene number).

The genes related to anthocyanin biosynthesis were then screened via GO_ID and pathway_ID of anthocyanin biosynthesis from significantly differentially expressed genes (DEGs) at different developmental stages of purple potatoes (Kanehisa, 2016; Michael et al., 2000).

Analysis of the levels of expression of miRNAs and corresponding target genes by qRT‐PCR

Total RNA was extracted from six different periods of the tubers (Sa‐Sf) using TRNzol (Tiangen, Beijing, China). The Mir‐X miRNA First‐Strand Synthesis Kit that was included in the Mir‐X miRNA qRT‐PCR TB Green Kit (TaKaRa Bio, San Jose, CA, USA) was used to convert miRNAs to cDNA to enable the quantification of specific RNAs using real‐time quantitative reverse transcription PCR (qRT‐PCR). The RNA molecules were polyadenylated and reverse transcribed using poly (A) polymerase and SMART® MMLV Reverse Transcriptase that was included in the mRQ Enzyme Mix. Then, TB Green Advantage® qPCR Premix and mRQ 3′ primers were used in qRT‐PCR in concert with miRNA‐specific 5′ primer(s) to quantify specific miRNA sequences in the cDNA. The total RNA was used to synthesize the cDNA of the target genes using the FastQuant cDNA Kit KR106 (Tiangen). The qRT‐PCR analysis was conducted using a miRcute Plus miRNA qPCR Kit FR411 (Tiangen). The gene‐specific primers are listed in Table S1. The U6 gene was used as an internal reference to accurately normalize each reaction for the miRNA and target genes of the potato, respectively. The relative levels of expression were calculated using the 2−ΔΔCt method (Livak & Schmittgen, 2001). Each reaction was conducted in triplicate with different cDNAs synthesized from three biological replicates.

RESULTS

Sequencing analysis of sRNAs

To identify the miRNAs involved in the biosynthesis of anthocyanins in purple potatoes, three sRNA libraries were constructed from the purple potato tubers at different developmental stages and sequenced (Table 1). After the low‐quality tags had been removed, the samples were compared and filtered with mRNA. The use of Rfam, including rRNA, tRNA, snRNA, and snoRNA among others, and the Repbase database provided 1,161,321 (S11), 11,019,800 (S12), 11,462,927 (S13), 9,684,037 (S21), 10,250,487 (S22), 11,103,152 (S23), 9,546,518 (S31), 10,437,838 (S32), and 11,855,603 (S33) valid reads from the purple potato tubers at different developmental stages, respectively.

TABLE 1

Data summary of three sRNA sequencing from purple potato

Sample	S11		S12		S13		S21		S22		S23		S31		S32		S33
	Total	% of Total	Total	% of Total	Total	% of Total	Total	% of Total	Total	% of Total	Total	% of Total	Total	% of Total	Total	% of Total	Total	% of Total
Raw reads	13,311,119	100.00	12,417,359	100.00	13,273,466	100.00	11,333,796	100.00	13,379,878	100.00	13,064,472	100.00	11,744,626	100.00	13,466,778	100.00	13,937,257	100.00
3ADT and length filter	1,303,981	9.80	862,828	6.95	1,138,160	8.57	948,467	8.37	2,105,126	15.73	1,237,310	9.47	1,393,612	11.87	2,175,125	16.15	1,124,424	8.07
Junk reads	44,022	0.33	56,462	0.45	54,648	0.41	40,435	0.36	45,551	0.34	62,869	0.48	46,180	0.39	44,900	0.33	53,889	0.39
Rfam	794,419	5.97	471,332	3.80	607,396	4.58	652,857	5.76	969,121	7.24	653,771	5.00	752,023	6.40	799,737	5.94	896,360	6.43
Repeats	17,330	0.13	14,379	0.12	23,409	0.18	19,499	0.17	27,015	0.20	18,643	0.14	28,007	0.24	23,332	0.17	20,166	0.14
Valid reads	11,161,321	83.85	11,019,800	88.75	11,462,927	86.36	9,684,037	85.44	10,250,487	76.61	11,103,152	84.99	9,546,518	81.28	10,437,838	77.51	11,855,603	85.06
rRNA	500,507	3.76	301,562	2.43	358,676	2.70	371,851	3.28	509,767	3.81	415,638	3.18	326,803	2.78	470,827	3.50	530,974	3.81
tRNA	188,232	1.41	89,437	0.72	149,922	1.13	183,652	1.62	335,785	2.51	121,876	0.93	342,338	2.91	205,787	1.53	261,557	1.88
snoRNA	17,428	0.13	11,387	0.09	14,899	0.11	16,256	0.14	21,354	0.16	26,647	0.20	13,911	0.12	15,507	0.12	17,955	0.13
snRNA	45,650	0.34	35,392	0.29	44,420	0.33	35,122	0.31	42,239	0.32	36,660	0.28	33,549	0.29	39,626	0.29	43,789	0.31
Other Rfam RNA	42,602	0.32	33,554	0.27	39,479	0.30	45,976	0.41	59,976	0.45	52,950	0.41	35,422	0.30	67,990	0.50	42,085	0.30

Analysis of the miRNAs in purple potato tubers

Based on the length of miRNA in each developmental stage (Figure 1), most of the miRNAs were 18–25 nt, with the highest distribution at 24 nt, comprising 48.33%. There were 1,753 known miRNAs, 1,744 conserved miRNAs, and 492 novel miRNAs (Table S2 ). In this study, the differentially expressed miRNAs were detected using pairwise comparisons at different developmental stages. A total of 1,915, 1,869, and 1,813 miRNAs were identified as differentially expressed in purple potato tubers at the S1, S2, and S3 developmental stages, respectively. A total of 1,334 miRNAs were co‐expressed between the three libraries (Figure 2a). A comparison of the differential expression of miRNAs from the three developmental stages identified 239, 195, and 157 miRNAs that were specifically differentially expressed between any two stages of S1, S2, and S3 (Figure 2a).

FIGURE 1

Length and size distribution of unique sequences in three sRNA libraries

FIGURE 2

Differentially expressed miRNAs. (a) Venn diagram of miRNAs detected (S1 vs. S2 vs. S3). (b) the volcanic diagram of S2 versus S1 differentially expressed miRNAs. (c) the volcanic diagram of S3 versus S1 differentially expressed miRNAs. (d) the volcanic diagram of S3 versus S2 differentially expressed miRNAs. S1, stage 1; S2, stage 2; S3, stage 3.

To define the differentially expressed miRNAs, 275 differentially expressed miRNAs were identified in the three libraries based on the level of P ≤ 0.05 (Table S3 ). There were 125 known miRNAs, 109 conserved miRNAs, and 41 novel miRNAs, which could match 58 known miRNA families. There were significant differences in the levels of expression of purple potato miRNAs in S1 (tuber formation stage), S2 (tuber bulking stage), and S3 (tuber maturation stage) (Figure S1). The levels of expression of the genes for 50 miRNAs were upregulated during the whole growth period, including stu‐miR156a and stu‐miR8043. The levels of expression of the genes of 70 miRNAs were downregulated throughout the growth period and included stu‐miR8019‐3P_L‐2 and stu‐miR397‐5P. In addition, the genes of 155 miRNAs, such as stu‐miR398b‐3p_L‐1 and stu‐miR8039_L‐1R + 4_1ss18CT1 among others, were expressed at different levels during the three stages of tuber development. A total of 54 miRNAs were upregulated, and 48 miRNAs were downregulated in S2 versus S1 (Figure 2b). A total of 68 miRNAs were upregulated, and 70 miRNAs were downregulated in S3 versus S1 (Figure 2c). A total of 35 miRNAs were upregulated, and 50 miRNAs were downregulated in S3 versus S2 (Figure 2d ).

Prediction and enrichment analyses of the miRNA target genes

A total of 6,050 miRNA target genes were predicted by 275 differentially expressed miRNAs, including 2,670 target genes annotated by GO and 1,453 annotated by the KEGG (Table S4 ). A GO enrichment analysis was performed to describe the properties of genes and gene products in the purple potato tubers. GO indicated that 2,296 of the small RNAs were divided into three types that included biological processes (1,228), molecular function (757), and cellular components (311) (Figure 3). GO was highly enriched in signal transduction, ADP binding, defense response, and cell differentiation (Figure 4a). GO terms that are involved in the biosynthesis of enzymes that are involved in the production of anthocyanin in potato tubers included phenylalanine, coumarin, chalcone isomerase, chalcone synthase, and malonyl‐CoA biosynthesis. The processes in potato tubers included the L‐phenylalanine catabolic process (GO:0006559), phenylalanine ammonia lyase activity (GO:0045548), response to phenylalanine (GO:0080053), 4‐coumarate‐CoA ligase activity (GO:0016207), flavonoid biosynthetic process (GO:0009813), malonyl‐CoA biosynthetic process (GO:2001295), malonyl‐CoA decarboxylase activity (GO:0050080), naringenin‐chalcone synthase activity (GO:0016210), chalcone biosynthetic process (GO:0009715), chalcone isomerase activity (GO:0045430), and anthocyanin biosynthesis (GO:0033729, GO:0046283).

FIGURE 3

GO annotated cluster map of the sRNA of DEGs classified in three main categories (cellular component, molecular function, and biological process) (S1 vs. S2 vs. S3). DEGs, differentially expressed genes; GO, gene ontology; sRNA, small RNA.

FIGURE 4

Functional enrichment of the miRNA target genes (S1 vs. S2 vs. S3). (a) GO enrichment analysis for sRNA DEGs (p < .05). (b) KEGG pathway enrichment analysis for DEGs in the transcriptome (P < .05). DEGs, differentially expressed genes; GO, gene ontology; KEGG, the Kyoto Encyclopedia of Genes and Genomes; sRNA, small RNA.

Similarly, a KEGG enrichment analysis was performed to describe the properties of genes and gene products in the tubers. The KEGG analysis indicated a high level of enrichment for plant hormone, signal transduction, spliceosome, phenylalanine metabolism, and flavonoid biosynthesis (Figure 4b). Some KEGG terms that are involved in the biosynthesis of anthocyanins in the purple potato tubers included phenylpropanoid biosynthesis (ko00360, ko00400), flavonoid biosynthesis (ko00944, ko00941), and anthocyanin biosynthesis (ko00942).

Thirty‐five miRNAs (Table 2) and target genes in the purple potato tubers (Table S5) were obtained through GO/KEGG annotation related to anthocyanin synthesis. There was a complex regulatory network between miRNA and the target genes (Figure 5).

TABLE 2

Summary of the identified miRNAs from purple potato in sRNA library

Identified miRNAs	Query miRNAs	Len.	Type	CG%	dG	MFEI
stu‐miR8006‐p3_1ss1GA	stu‐MIR8006‐p3_1ss1GA	24	3’	35.1	−26.4	0.8
stu‐miR8006‐p5_1ss11CT_2	stu‐MIR8006‐p5_1ss11CT_2	24	5’	34.7	−8.3	0.5
stu‐miR8006‐p5_1ss16CT	stu‐MIR8006‐p5_1ss16CT	24	5’	37.6	−18.2	0.6
stu‐miR8006‐p3_1ss16CT_2	stu‐MIR8006‐p3_1ss16CT_2	24	3’	39.5	−27.9	0.6
stu‐miR172b‐5p_L‐1R‐2	stu‐miR172b‐5p_L‐1R‐2	18	5’	41.6	−3.3	0.1
stu‐miR5021	csi‐MIR396c‐p5_2ss6CT23TG	24	5’	25.7	−11.7	0.3
stu‐miR8019‐3p_L‐2	stu‐miR8019‐3p_L‐2	22	3’	37.1	−4.4	0.2
stu‐miR11020	mdm‐MIR11020‐p3_2ss1GA17GA	18	3’	23.1	−2.4	0.6
stu‐miR5303a‐3p_R + 2	Nta‐miR5303a‐p3_2ss21GA22TC	24	3’	46.3	−97.1	1
stu‐miR5303f	stu‐miR5303f	24	3’	4.4	−20.3	1.1
stu‐miR6164b	nta‐MIR6164b‐p3_2ss17CG18AG	24	3’	37	−6.6	0.3
stu‐miR6164b_1	nta‐MIR6164b‐p3_2ss15CT17CG_1	21	3’	34.8	−19.2	0.4
stu‐miR6164b_2	nta‐MIR6164b‐p3_2ss15CT17CG_2	24	3’	40	−2.4	0.2
stu‐miR391‐3p_L + 2R‐4	stu‐MIR391‐p3	19	3’	39.8	−36	1
stu‐miR6160	PC‐3p‐42378_149	24	3’	35.3	−123.9	1.7
stu‐miR8043	stu‐miR8043	21	3’	25.6	−73.8	2.3
stu‐miR7711	PC‐3p‐2911_1,262	24	3’	31.3	−8.9	1.7
stu‐miR470	PC‐5p‐6180508_2	25	5’	35.9	−52.2	1.1
stu‐miR9471	PC‐5p‐768978_11	24	5’	35	−48	1.1
stu‐miR6873_1	mes‐MIR172d‐p5_2ss14TA19GC_1	21	5’	48.8	−29.1	0.7
stu‐miR6873	mes‐MIR172d‐p3_2ss14TA19GC	21	3’	33.7	−38.5	0.7
stu‐miR6873_2	mes‐MIR172d‐p5_2ss14TA19GC_2	22	5’	40.2	−18.3	0.4
stu‐miR9409	bol‐MIR9409‐p3_2ss16AC19TA	19	3’	23.4	−1.3	0.1
stu‐miR8006	stu‐MIR8006‐p3_2ss8CT20GA	20	3’	37.2	−34.7	0.6
stu‐miR5253	mtr‐MIR5253‐p3_2ss15AG17CT	19	3’	36.4	−9.9	0.6
stu‐miR8031_L‐7R + 8	stu‐MIR8031‐p3_1ss15TA	24	3’	33.3	−2.2	0.2
stu‐miR156a_R‐1	stu‐miR156a_R‐1	20	3’	35.6	−25.2	0.6
stu‐miR156a	stu‐miR156a	21	5’	38.4	−28.6	0.5
stu‐miR156a_1ss11AC	stu‐miR156a_1ss11AC	21	3’	45.5	−7.1	0.4
stu‐miR156a_R + 1	fve‐miR156f_R + 1	22	5’	43.3	−47.8	1.2
stu‐miR482	sly‐miR482e‐3p_R‐2	20	5’	47.4	−7.4	0.3
stu‐miR5303h_2ss23AG24AT	sly‐MIR5303‐p5_2ss13CT24CT	24	5’	38.6	−45.9	0.9
stu‐miR9605	PC‐3p‐104668_67	24	3’	38.6	−97.4	1.2
stu‐miR5303g_L‐3	sly‐miR5303	21	3’	41.3	−41.6	1
stu‐miR8006‐p5_1ss1GA	stu‐MIR8006‐p5_1ss1GA	24	5’	28.4	−17.1	0.5

FIGURE 5

Map of the interaction between the miRNAs related to the synthesis of anthocyanin and its target gene network. (a) Genes related to anthocyanin annotated by GO. (b) Genes related to anthocyanin annotated by KEGG. Blue, miRNA; yellow, target genes; orange, GO/KEGG. GO, gene ontology; KEGG, the Kyoto Encyclopedia of Genes and Genomes.

Degradome sequencing analysis

A total 69,387,200 raw reads were obtained from three degradome libraries. Of these, 4,655,258 unique raw reads were in S1; 6,032,528 unique raw reads were in S2, and 6,078,819 unique reads were in S3. The mapped reads represented 16,920,374, 26,135,717, and 25,789,524 of purple potato genes in the S1, S2, and S3 stages (Table 3).

TABLE 3

Data summary of three degradome sequencing libraries from purple potato

Sample	S1(number)	S1(ratio)	S2(number)	S2(ratio)	S3(number)	S3(ratio)	Sum (number)	Sum (ratio)
Raw reads	17,049,581	/	26,348,452	/	25,989,167	/	69,387,200	/
Reads < 15 nt after removing 3 adaptor	129,207	0.76%	212,735	0.81%	199,643	0.77%	541,585	0.78%
Mappable reads	16,920,374	99.24%	26,135,717	99.19%	25,789,524	99.23%	68,845,615	99.22%
Unique raw reads	4,655,258	/	6,032,528	/	6,078,819	/	12,777,354	/
Unique reads < 15 nt after removing 3 adaptor	37,336	0.80%	48,377	0.80%	48,457	0.80%	100,797	0.79%
Unique mappable reads	4,617,922	99.20%	5,984,151	99.20%	6,030,362	99.20%	12,676,557	99.21%
Transcript mapped reads	13,369,060	78.41%	21,446,902	81.40%	19,923,925	76.66%	54,739,887	78.89%
Unique transcript mapped reads	3,457,918	74.28%	4,470,681	74.11%	4,305,718	70.83%	8,790,805	68.80%
Number of input transcript	62,784	/	62,784	/	62,784	/	62,784	/
Number of Coverd transcript	47,054	74.95%	46,847	74.62%	46,250	73.67%	50,031	79.69%

A miRNA–mRNA relationship was simultaneously detected by the prediction of target genes, and the results of experimental sequencing identified 11,505 relationships. There were 11,248 genes that differed significantly (P ≤0.05). Cleaved targets that had an alignment score of 4 or less were considered the potential targets. For these genes, 417, 92, 2,516, 246, and 2,118 genes belonged to categories 0, 1, 2, 3, or 4 in the S1 library; 342, 71, 3,256, 666, and 2,389 belonged to categories 0, 1, 2, 3, or 4 in the S2 library, and 238, 44, 2,971, 582, and 2,437 belonged to a category < 4 in the S3 library. A total of 3,238 miRNAs were obtained from the degradation data, and the target genes were predicted by TargetFinder (Table S6). The total number of target genes was 11,505, which belonged to 1,470 miRNAs, after combining the predicted miRNA target genes with the mRNA in the degradation group density file (Table S7). A total of 24,530 transcripts were obtained. The total number of interacting miRNAs that corresponded to all the transcripts was 11,505 after they were predicted by TargetFinder and combined with the mRNA in the generated degradation group density file (Table S8). A T‐plot diagram was used to intuitively display the target genes detected by miRNA. A total of 5,731 target genes were identified from the S1 libraries; 7,157 target genes were identified from the S2 libraries, and 6,751 target genes were identified from the S3 libraries.

The degradome libraries included 3,954 GO terms that were divided into three types, including biological processes (1,199), molecular function (627), and cellular components (314) (Figure 6). GO is highly enriched in the nucleus, protein binding, ATP binding, and defense response (Figure 7a). The GO terms involved in the anthocyanin biosynthetic enzymes in the tubers included the phenylalanine metabolic pathway (GO:0006559 GO:0080130, GO:0045548), coumarate‐CoA (GO:0016207), flavonoid biosynthetic process (GO:0009813), chalcone isomerase activity (GO:0045430), and regulation of the anthocyanin metabolic process (GO:0009718, GO:0046283, GO:0031538, GO:0031539, GO:0031540, and GO:0031537).

FIGURE 6

GO annotated cluster map of the degradome of DEGs that are classified in three main categories (cellular component, molecular function, and biological process) (S1 vs. S2 vs. S3). DEGs, differentially expressed genes; GO, gene ontology.

FIGURE 7

GO and KEGG enrichment analysis for degradome DEGs (S1 vs. S2 vs. S3). (A) GO enrichment analysis for the degradome library genes. (B) KEGG pathway enrichment analysis for the degradome library genes. DEGs, differentially expressed genes; GO, gene ontology; KEGG, the Kyoto Encyclopedia of Genes and Genomes.

A total of 96 KEGG terms were discovered in the degradome libraries. There were five KEGG classifications, including Organismal Systems, Metabolism, Genetic Information Processing, Environmental Information Processing, and Cellular Processes (Figure 8). The KEGG was highly enriched in Ribosome biogenesis in the spliceosome, plant–pathogen interaction, and plant hormone signal transduction (Figure 7b ). Some KEGG terms that are involved in the biosynthesis of anthocyanins in potato include phenylpropanoid biosynthesis (ko00940), flavonoid biosynthesis (ko00941), and anthocyanin biosynthesis (ko00942).

FIGURE 8

KEGG classification of the degradation libraries in purple potato ( ) (S1 vs. S2 vs. S3). KEGG, the Kyoto Encyclopedia of Genes and Genomes.

Multi‐omics analysis: Anthocyanin‐responsive miRNA–mRNA modules

The miRNA–mRNA interactions were identified from two sequencing datasets. The miRNA sequences were identified by the sRNA‐seq results, and the relationship between miRNA and the target gene was identified from the degradome‐seq results. The enrichment of miRNA target genes in the comparison group was analyzed by GO and KEGG, and the annotation information of each function or pathway that corresponded to all the selected miRNA target genes was counted, which enabled the identification of miRNA target genes related to anthocyanin biosynthesis (Table 4). We obtained 37 miRNAs and 23 target genes based on the annotation information (Figure 9). A total of 147 mature miRNAs from tomato, and 37 mature miRNAs from the purple potato tubers were used to construct a phylogenetic tree to help understand their relationship (Figure 10). The mature miRNAs were classified into four different groups, and each group could affect related types of anthocyanin biosynthesis. Stu‐miR156 and sly‐miR156 clustered into one group, and their genetic relationship was relatively close at 96%. This suggests that they have the same function.

TABLE 4

miRNAs and target genes related to anthocyanins in purple potato

Identified miRNAs	Query miRNAs	Len	Type	CG%	dG	MFEI	Target gene	Describe	GO	KEGG
stu‐miR8045	stu‐miR8045_1ss2TC	21	3’	3.30	−28.80	1.10	PGSC0003DMT400030428	CHI	GO:0009813,GO:0045430	ko00941
stu‐miR8039	stu‐miR8039	21	5’	30.00	0.60	−0.10
stu‐miR8039_1ss21CG	stu‐miR8039_2ss5CA21CG	21	3’	34.20	−17.80	0.50
stu‐miR8039_R + 1	stu‐miR8039_R + 1_1ss4CT	22	5’	32.80	−56.50	1.40
stu‐miR9477‐3p	sly‐MIR9477‐p5_2ss4GA22GA	24	5’	38.20	−8	0.40
stu‐miR9477‐3p_1	sly‐miR9477‐5p_L + 1R‐4_1ss6CA	21	5’	35.70	−102.80	1.30
stu‐miR9477‐5p	sly‐miR9477‐5p_L‐1R + 1_1ss5CA	24	5’	38.20	−63.90	1.10
stu‐miR8039_L + 3	sly‐MIR9477‐p3_2ss3AG24TC	24	3’	34.30	−12.80	0.60
stu‐miR9477‐5p_1	sly‐miR9477‐5p_1ss5CA	24	5’	4.80	−51.60	1.20
stu‐miR9477	PC‐3p‐30742_196	21	3’	33.10	−48.70	1.20
stu‐miR8045_R + 3	PC‐3p‐5439628_3	24	3’	35.70	−9.60	1.20
stu‐miR8039	PC‐3p‐914239_10	24	3’	3.30	−48.90	0.90
stu‐miR8039_L + 3_2	PC‐5p‐250627_30	24	5’	31.80	−69.40	1.10
stu‐miR8014‐3p_L‐3	stu‐miR8014‐3p_L‐2_2ss8CT21AG	22	3’	28.30	−6.10	0.40	PGSC0003DMT400009176	F3H	GO:0009813	ko00941
stu‐miR414	gma‐MIR4348d‐p3_2ss11AG17AT	19	3’	45.20	−7.30	0.40	PGSC0003DMT400060308	PAL	GO:0006559,GO:0045548	ko00360,ko00940
stu‐miR361	mtr‐MI171c‐p3_1ss2TG	19	3’	30.00	−29.60	0.80	PGSC0003DMT400008182	4CL	GO:0016207	ko00360,ko00940
stu‐miR396	hbr‐miR396a_R‐3_1ss17TC	18	5’	37.20	−36	1.10	PGSC0003DMT400062562	UFGT	GO:0009813	ko00942
stu‐miR396b‐5p	stu‐miR396‐5p_R‐3	18	3’	33.60	−22.50	0.50
stu‐miR396b‐5p_1	stu‐miR396‐5p_R‐2	19	5’	39.10	−2.90	.20
stu‐miR828	ptc‐MIR828b‐p5	20	5’	37.40	−16.90	0.40	PGSC0003DMT400078477	MYB	GO:0009718
stu‐miR828a	ptc‐miR828a	22	5’	36.90	−24.70	0.70
stu‐miR172e‐5p_L‐1R‐1	stu‐miR172e‐5p_L‐1R‐1	19	3’	33.70	−9.70	0.30	PGSC0003DMT400060467,PGSC0003DMT400060468,PGSC0003DMT400060469	WD40	GO:0031538
stu‐miR8019‐5p_L‐4R‐3	stu‐miR8019‐3p_L‐4R‐2	18	3’	39.10	−13.10	0.50	PGSC0003DMT400075906,PGSC0003DMT400075907,PGSC0003DMT400075908	MYB	GO:0031537
stu‐miR29b‐4‐5p	gra‐MIR8771d‐p3_2ss14AT18GA_2	19	3’	19.60	−13.80	0.70	PGSC0003DMT400020601,PGSC0003DMT400020605,PGSC0003DMT400020597,PGSC0003DMT400030504	UTGT	GO:0009813	ko00942
stu‐miR5303c_L‐1R‐2_3ss19CA20AG21CT	nta‐MIR5303c‐p3_1ss19GA	21	3’	48.30	−16	0.60	PGSC0003DMT400011627,PGSC0003DMT400011628	UFGT	GO:0009813	ko00942
stu‐miR5303f_L‐7R + 2	stu‐MIR5303g‐p5_1ss19CG	19	5’	43.60	−45.30	1.10
stu‐miR5303f_L + 1R‐3	nta‐MIR5303c‐p3_2ss1CT17GA	22	3’	44.90	−67	1.00
stu‐miR7997a_L‐3	PC‐5p‐831816_11	24	5’	39.40	−47	1.20
stu‐miR7997b_L‐3	PC‐5p‐22694035_1	24	5’	41.60	−57.70	1.10
stu‐miR7997c_L + 3R‐5_2ss2TA3AG	PC‐5p‐315982_25	21	5’	44.90	−58.50	0.90
stu‐miR156f‐5p_L + 1	ta‐miR156a_L + 2R‐1	21	5’	41.00	−55.10	1.30	PGSC0003DMT400029748,PGSC0003DMT400029749,PGSC0003DMT400029750,PGSC0003DMT400029751	SPL
stu‐miR156e_R + 1	Ptc‐miR156a_R + 1	21	5’	42.10	−97.60	1.20
stu‐miR156a	stu‐miR156a	21	5’	38.40	−28.60	0.50
stu‐miR156a_R‐1	stu‐miR156a_R‐1	20	3’	35.60	−25.20	0.60
stu‐miR156e	stu‐miR156e	20	5’	44.40	−55.90	1.30
stu‐miR858	ppe‐R858‐p5_1ss5GA_1	19	5’	36.00	−4.90	0.30	PGSC0003DMT400042007	MYB
stu‐miR5021	ath‐MIR5021‐p5_2ss4AG17AC	18	5’	45.50	−42.20	0.60	PGSC0003DMT400075908,PGSC0003DMT400075907,PGSC0003DMT400075906	MYB	GO:0031537
							PGSC0003DMT400060469,PGSC0003DMT400060468	WD40	GO:0031538

FIGURE 9

Genes related to anthocyanins by GO/KEGG. Blue, miRNA; yellow, target genes; orange, GO/KEGG. GO, gene ontology; KEGG, the Kyoto Encyclopedia of Genes and Genomes.

FIGURE 10

A total of 147 mature miRNAs from tomato ( ) and 37 mature miRNAs from potato ( ) were used to construct a phylogenetic tree

An analysis of the regulatory pattern between miRNA and the target genes related to the biosynthesis of anthocyanin following sRNA‐seq and degradome sequencing indicated that both the miRNAs and target genes had a pattern of negative regulation. There were 27 upregulated miRNAs and 10 downregulated miRNAs in S2 versus S1, 25 upregulated miRNAs and 12 downregulated miRNAs in S3 versus S1, and nine upregulated miRNAs and 28 downregulated miRNAs in S3 versus S2 (Figure 11a). The analysis of the target genes that regulate anthocyanins indicated that they can be divided into two categories. One is the key genes that are involved in anthocyanin biosynthesis, which include the structural genes, PAL, 4CL, CHI, F3H, and UFGT. The other is to regulate the biosynthesis of anthocyanins, which includes the regulatory genes MYB, WD40, and SPL.

FIGURE 11

Analysis of the genes related to anthocyanins. (a) Hierarchical cluster analysis of the miRNA–mRNA related to anthocyanin. Red, upregulated; blue, downregulated. (b) Analysis of the expression of genes related to anthocyanins in different development stages. Orange, high expression; blue, low expression. (c) Anthocyanin contents in the tubers of pigmented potato ( ) at different stages. (Sa and Sb, tuber formation stage; Sc, Sd, and Se, tuber bulking stage; Sf, tuber maturation stage)

Anthocyanin content and related gene expression in purple potato tubers

The total anthocyanin content showed a trend of increasing in pigmented potato tubers during maturation, although there was a slight obvious downward trend in the Se stages (Figure 11c ). Theiler's murine encephalomyelitis virus (TMEV) was used to analyze the differentially expressed miRNA–mRNA (Figure 11b). The levels of expression of PAL, 4CL, CHI, and F3H, which are related to the synthesis of anthocyanin, decreased as the potatoes matured. The three UFGT structural genes, PGSC0003DMT400062562, PGSC0003DMT400011627, and PGSC0003DMT400011628, were expressed in the same manner. Their levels of their expression were the highest in S1 period. However, the PGSC0003DMT400020605 gene was expressed in an opposite manner. The levels of expression of the UFGT genes PGSC0003DMT400020601 and PGSC0003DMT400030504 were the highest in S1 and S3 stages and the lowest in S2. Among the regulatory genes related to anthocyanin biosynthesis, there were two modes of expression of the MYB transcription factor. The levels of expression of PGSC0003DMT400078477 and PGSC0003DMT400042007 were the highest in S1 and S2 stages. However, the levels of expression of PGSC0003DMT400075906, PGSC0003DMT400075907, and PGSC0003DMT400075908 were the highest in S1 and S3 stages. The WD40 gene is expressed in a similar manner and is expressed the most highly in the S2 stage.

The levels of expression of nine miRNAs, including stu‐miR396 and stu‐miR156e, were upregulated during the whole growth period of potato, while the levels of expression of 15 miRNAs, including stu‐miR828a and stu‐miR8045, were downregulated. The levels of expression of miRNAs stu‐miR828 and stu‐miR8045_R + 3 first increased and then decreased as the potatoes matured. Moreover, the level of expression of 11 miRNAs, including stu‐miR172e‐5p_L + 1R‐3 and stu‐miR5303f_L‐7R + 2, did not change significantly during the whole growth period. The regulatory mechanism of 15 miRNAs and their target genes from the anthocyanin biosynthetic pathways at the six growth stages were analyzed using qRT‐PCR to verify the reliability of sequencing data and understand the relationship between the level of expression of miRNA–mRNA genes and the accumulation of anthocyanin in purple potato tubers. The target genes expressed an opposite trend with miRNAs, suggesting that they might be actively cleaved by the miRNAs (Figure 12). The UFGT genes PGSC0003DMT400030504, PGSC0003DMT400020597, PGSC0003DMT400020601, and PGSC0003DMT400020605 were all regulated by stu‐miR29b‐4‐5p and tended to have a decrease in the level of their expression during maturation of the potato, while stu‐miR29b‐4‐5p exhibited a reverse trend in all the developmental stages of potato tubers (Figure 12a, Figure S2A‐C). The PGSC0003DMT400062562 gene was regulated by stu‐miR396b‐5p, and they were significantly negatively correlated (Figure 12b). PGSC0003DMT400011627 and PGSC0003DMT400011628 are all regulated by stu‐miR7997a‐p5_1ss16GA, stu‐miR5303f_L‐7R + 2, stu‐miR7997a_L‐3, and stu‐miR7997b_L‐3, stu‐miR7997c_L + 3R‐5 _2ss2TA3AG (Figure S2D‐M). PGSC0003DMT400078477 was regulated by stu‐miR828a and stu‐miR828 (Figure 12c‐d), while PGSC0003DMT400075906 was regulated by stu‐miR8019‐5p_L‐4R‐3 and stu‐miR5021 (Figure 12e‐f), which had reverse patterns of expression at all the stages of tubers. PGSC0003DMT400060467 was regulated by stu‐miR172e‐5p_L‐1R‐1 and stu‐miR5021 (Figure 12g‐h). The SPL9 genes were regulated by stu‐miR156 (Figure 12i, Figure S2O‐Q), and they were significantly negatively correlated. There was a significant negative correlation between stu‐miR858 and PGSC0003DMT400042007 (Figure 12j).

FIGURE 12

Analyses of the expression of miRNAs and their target genes by qRT‐PCR. (a) Stu‐miR29b‐4‐5p and PGSC0003DMT400020605. (b) Stu‐miR396b‐5p and PGSC0003DMT400062562. (c) Stu‐miR828 and PGSC0003DMT400078477. (d) Stu‐miR828a and PGSC0003DMT400078477. (e) Stu‐miR8019‐5p_L‐4R‐3 and PGSC0003DMT400075906. (f) Stu‐miR5021 and PGSC0003DMT400075906. (g) Stu‐miR172‐5p_L‐1R‐1 and PGSC0003DMT400060467. (h) Stu‐miR5021 and PGSC0003DMT400060467. (i) Stu‐miR156f‐5p_L + 1 and PGSC0003DMT400029750. (j) Stu‐miR858 and PGSC0003DMT400042007. qRT‐PCR, real‐time quantitative reverse transcription PCR.

DISCUSSION

Small RNAs (sRNAs) have recently been recognized as key genetic and epigenetic regulators in various organisms. Their functions range from the modification of DNA and methylation of histone to the modulation of abundance of coding or non‐coding RNAs. Major regulatory sRNAs in plants are classified as miRNA and siRNA, with the miRNAs primarily engaging in posttranscriptional regulation, while the siRNAs are involved in transcription (Chen et al., 2018). miRNAs are involved in a variety of plant regulatory pathways, including cell development, plant secondary metabolism, and an increase in the resistance to plant stress. Anthocyanins are one of the important groups of secondary metabolites in plants. Some miRNAs, such as miR165a‐5p, miR172b, miR827a, miR166g and miR1432–5 (Gao et al., 2019), miR156 (He et al., 2019), miR828 (Bonar et al., 2018), and miR858, are involved in the regulation of anthocyanin synthesis. However, miRNAs that regulate anthocyanins in colored potato have not been systematically reported. Therefore, identification of miRNA and its target genes in pigmented potato will help us to understand the regulatory mechanism of miRNA. In this study, two important high‐throughput methods, small RNA and degradome sequencing, were utilized to investigate the regulatory mechanism of anthocyanin in different developmental stages in pigmented potato. There were 275 differentially expressed miRNAs in the three libraries, which belong to 58 known miRNA families. This is similar to the results of 277 miRNAs identified by Qiao et al. (2017) in four libraries of potato leaves and tubers under light and dark treatments. Compared with Zhang et al. (2013), 16 more miRNAs were identified. Although many miRNA genes are conserved across plant species, the same gene family can vary significantly in size and genomic organization in different species (Li & Mao, 2007). The proportion of miRNA of 24 nt (48.33%) was much higher than that of 21 nt (14.84%). This result is consistent with the length distribution of miRNA in pitaya (H. monacanthus) (Chen et al., 2020) and kiwifruit (A. arguta) (Li et al., 2019) but not in passion fruit (Passiflora edulis) (Paul et al., 2020). The results showed that the miRNAs of different plants differed in the distribution of their lengths.

The biosynthesis of anthocyanin is a complex process, which is regulated by multiple genes. In this study, we obtained 35 miRNAs related to anthocyanins in the sRNA libraries, 37 miRNA genes in multi‐omics analyses, and two categories of target genes that are related to anthocyanins. Similarly, in sweet potatoes (He et al., 2019), 26 differentially expressed miRNAs and 36 corresponding targets were potentially involved in the anthocyanin biosynthesis. In this study, structural genes related to anthocyanin include stu‐miR396b‐5p target UFGT (PGSC0003DMT400062562) (Figure 13a), stu‐miR8039 target CHI (PGSC0003DMT400030428) (Figure 13b ), stu‐miR8014‐3p_L‐3 target F3H (PGSC0003DMT400009176) (Figure 13c ), stu‐miR414 target PAL (PGSC0003DMT400060308) (Figure 13d), and stu‐miR361 target 4CL (PGSC0003DMT400008182) (Figure 13e). Previous studies have found that PAL, CHI, 4CL, F3H, and UFGT can promote or inhibit the synthesis of anthocyanins (He et al., 2020; Leng et al., 2020; Wang et al., 2019). The patterns of expression of PAL, CHI, 4CL, F3H, and UFGT were the same in the solanaceous black nightshade (Solanum nigrum Linn.) (Saophea et al., 2020), grape (Vitis vinifera L. × Vitis labrusca L.) berry skin (Wang et al., 2013), Taxus chinensis (Zhang et al., 2019), and the purple‐head trait of Chinese cabbage (B. rapa L.) (He et al., 2020), which promoted the accumulation of anthocyanins. However, in black and white mulberry (Morus nigra and Morus alba, respectively), the expression of F3H was upregulated in white and downregulated in black mulberry, which indicated that the expression of F3H inhibited the synthesis of anthocyanin (Huang et al., 2020). There were two patterns of expression of UFGT (3GT) in Freesia hybrida. FH3GT2 decreased with flower development and inhibited the accumulation of anthocyanin, while FH3GT1 increased with flower development and promoted the accumulation of anthocyanin (Meng et al., 2019). In this study, the levels of expression of PAL, CHI, 4CL, F3H, and UFGT negatively correlated with those of the miRNA, which could be involved with the biosynthesis of anthocyanin in purple potato (Figure 14).

FIGURE 13

Target plots that are cleaved by miRNAs. (a) Stu‐miR396b‐5p target UFGT (PGSC0003DMT400062562). (b) Stu‐miR8039 target CHI (PGSC0003DMT400030428). (c) Stu‐miR8014‐3p_L‐3 target F3H (PGSC0003DMT400009176). (d) Stu‐miR414 target PAL (PGSC0003DMT400060308). (e) Stu‐miR361 target 4CL (PGSC0003DMT400008182). (f) Stu‐miR828 target MYB (PGSC0003DMT400078477). (g) Stu‐miR172e‐5p_L‐1R‐1 target WD40 (PGSC0003DMT400060467). (h) Stu‐miR156e target SPL9 (PGSC0003DMT400029478). 4CL, 4‐coumarate CoA ligase; CHI, chalcone isomerase; F3H, flavanone 3‐hydroxyalse; PAL, phenylalanine ammonia lyase; UFGT, anthocyanidin 3‐O‐glucsoyltransferase.

FIGURE 14

Biosynthetic pathway of anthocyanin and related genes in purple potato ( )

Simultaneously, we found two types of transcription factors, MYB and WD40, and regulatory genes, such as SPL. This is similar to the results in sweet potato (He et al., 2019). In this study, stu‐miR828 cleaved MYB (PGSC0003DMT400078477) (Figure 13f), while stu‐miR858 cleaved MYB (PGSC0003DMT400042007). Typically, conserved miRNAs have the same or homologous target genes with other plants, and most miRNAs have similar functions. BrmiR828 was found to negatively regulate the transcription of MYB82 (bra022602) and promoted the accumulation of anthocyanin following the treatment of turnip (B. rapa) seedlings by light (Zhou et al., 2020). In the red flesh potato, the MYB transcription factor was predicted to be a potential target gene of miR828, and its level of expression decreased (Bonar et al., 2018). miR828 targets MYB and may regulate the biosynthesis of anthocyanin in the purple potato tubers.

The WD40 transcription factor was identified in apple (Malus domestica) (An et al., 2012), the flowers of pagoda tree (Sophora japonica L.) (Guo et al., 2021), sweet potato (Ipomoea batatas) cultivars (Dong et al., 2014), and pigmented potato tissues (Liu et al., 2020), which promoted the biosynthesis of anthocyanin. In this study, the target gene of stu‐miR172e‐5p_L‐1R‐1 encodes the WD40 transcription factor (Figure 13g). stu‐miR172e‐5p_L‐1R‐1 was downregulated at the Sc stage, while the pattern of expression of the WD40 gene was significantly upregulated. The increased level of expression of the WD40 gene might promote the accumulation of anthocyanin. In addition to MYB and WD40 transcription factors and structural genes, SPL is also involved in the synthesis of anthocyanin. VvmiR156b/c/d‐mediated VvSPL9 participated in the formation of grape (V. vinifera) color in response to multi‐hormone signals (Su et al., 2021). The expression of VcSPL12 significantly enhanced the accumulation of chlorophyll and altered the level of expression of several genes associated with chlorophyll that are involved with the coloration of blueberry (Vaccinium spp.) fruit (Li, Wu, et al., 2020). These findings provide novel insight into the functional roles of miR156‐SPLs in plants, particularly in purple potato. In this study, SPL9 that was targeted by stu‐miR156 was found to be involved in the biosynthesis of anthocyanin (Figure 13H). Similarly, miRNA–mRNA related to anthocyanin biosynthesis was also identified in sweet potato (He et al., 2019). miR858, miR156, miR172, and miR396 related to anthocyanins were identified. The mRNAs of miR858 and miR156 were identical, which were miR858‐MYB and miR156‐SPL. The mRNAs of miR172 and miR396 were different. miRNA–mRNA has some differences and conservation in species.

CONFLICT OF INTEREST

The authors Jie Liu and Wei Wei are employed by HuaSong Seed Industry Co., Ltd. (Beijing). The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

AUTHOR CONTRIBUTIONS

YM conceived and designed the research. XW and YM performed the experiments. JL and WW provided the potato seed. XW, YM, JW, PW, ZZ, RX, BF, LN, and XL analyzed the data and wrote the manuscript. All the authors read and approved the final manuscript.

Table S1. qRT‐PCR primers for the miRNAs and target genes related to anthocyanin synthesis in purple potato.

Click here for additional data file.

Table S2. All the expressed miRNAs – a summary of known and predicted miRNA in this study.

Click here for additional data file.

Table S3. A total of 275 differentially expressed miRNAs. P ≤ .05.

Click here for additional data file.

Table S4. Target genes of 275 differentially expressed miRNAs obtained by GO and KEGG annotation. GO, gene ontology; KEGG, the Kyoto Encyclopedia of Genes and Genomes.

Click here for additional data file.

Table S5. Thirty‐five miRNAs related to anthocyanins were obtained by GO and KEGG in the miRNA library. GO, gene ontology; KEGG, the Kyoto Encyclopedia of Genes and Genomes.

Click here for additional data file.

Table S6. The degradation data provided 3,238 miRNAs and target genes.

Click here for additional data file.

Table S7. A combination of the predicted miRNA target genes with the mRNA in the degradation group density file. The total number of interacting miRNAs corresponding to all the transcripts.

Click here for additional data file.

Table S8. The data after prediction by TargetFinder software combined with the mRNA in generated degradation group density file. The total number of interacting miRNAs corresponds to all the transcripts.

Click here for additional data file.

Figure S1. A total of 275 differentially expressed miRNAs. P ≤ .05.

Click here for additional data file.

Figure S2. Analyses of the expression of miRNAs and their target genes by qRT‐PCR.

Click here for additional data file.