{'fraction_non_alphanumeric': 0.03737255442270598, 'fraction_numerical': 0.014776798015982364, 'mean_word_length': 4.800799200799201, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 1, 'https://': 2, 'lorem ipsum': 0, 'www.': 8, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Validation of suitable reference genes by various algorithms for gene expression analysis in Isodon rubescens under different abiotic stresses

Isodon rubescens (Hemsley) H. Hara (Lamiaceae) is a traditional Chinese medicine plant that has been used to treat various human diseases. Oridonin is one of the main active ingredients, and the route of its molecular biosynthesis remains to be determined. The study of gene expression patterns can provide clues toward the understanding of its biological functions. The selection of suitable reference genes for normalizing target gene expression is the first steps in any quantitative real-time PCR (RT-qPCR) gene expression study. Therefore, validation of suitable reference genes is necessary for obtaining reliable results in RT-qPCR analyses of I. rubescens. Here, 12 candidate reference genes were chosen, and their expression stability in different tissues of I. rubescens and in leaves under different abiotic stresses (NaCl, dehydration, SA, MeJA, and ABA) was evaluated using the ∆Ct, NormFinder, GeNorm, BestKeeper, and RankAggreg statistical tools. Analysis using the comprehensive tools of RankAggreg algorithm showed that GADPH, 18S and eIF were stably expressed in different tissues; UBQ, Apt, and HIS; Cycl, UBQ, and PP2A; GADPH, 18S, and eIF; eIF, UBQ, and PP2A; TUB, Cycl, and UBQ; were the best three candidate reference genes for the samples of Dehydration, NaCl, SA, MeJA, and ABA treatment, respectively. While for the concatenated sets of ND (NaCl and dehydration) and SMA (SA, MeJA, and ABA), UBQ, HIS, and TUA; UBQ, eIF and Apt were the three appropriate candidate reference genes, respectively. In addition, the expression patterns of HMGR in different tissues and under different treatments were used to confirm the reliability of the selected reference genes, indicating that the use of an inappropriate reference gene as the internal control will cause results with a large deviation. This work is the first study on the expression stability of reference genes in I. rubescens and will be particularly useful for gene functional research in this species.

research of Chinese herbal medicines, the genes function of genes related to specific medicines are becoming more scrutinized. The screening of internal reference genes of Chinese herbal medicines has also made some progress, such as with Panax ginseng 8 , Salvia miltiorrhiza 9 , and Ganoderma lucidum 10 . Therefore, with the increase of studies on I. rubescens genes, it is of great significance to study the screening of its internal reference genes.
In the study of gene expression, RT-qPCR has been widely used as a powerful technique to determine the gene expression levels in different samples, due to its high sensitivity, efficiently, specificity, and reproducibility 11 . However, the accuracy of RT-qPCR results strongly depends on the precise normalization of the transcripts using stably expressed reference genes 12 . Reference genes are considered to be those with constant expression in all samples of the organism, including different tissues, or under different treatments, or at different developmental stages. Previous studies have shown that the common reference genes in many plants are those coding for elongation factor 1-α (EF-1α) 13 , actin (ACT ) 14 , glyceraldehyde-3-phosphate dehydrogenase (GADPH) 15 , and 18S ribosomal RNA (18S) 16 , to name a few. However, it can be seen from previous studies that there is no universal reference gene across different species, and the most stable reference genes are different in different plants or under different treatments of the same plant 8 . And the use of inappropriate reference genes can result in target gene expression profile bias. Therefore, to get accurate and reliable gene expression profiles, it is important to identify the optimum reference genes among different tissues and treatments when considering evaluation by RT-qPCR.
In order to better screen out the most suitable internal reference genes, various statistical algorithms are used to calculate or evaluate the stability of reference genes in a set of samples. Based on previous reports, the main calculation tools for qRT-PCR data normalization included the ∆Ct method 17 , geNorm 18 , NormFinder 19 , and BestKeeper 20 . These methods select the most suitable gene or gene combination based on the calculation of expression stability for normalization. In addition, the calculation tools listed above use different algorithms and therefore might generate different stability ranking results 21,22 . RankAggreg, an R package for weighted rank aggregation, was used in a comprehensive sequencing study on genetic stability 23 .
In addition, a large number of previous studies have shown that the main components of plants consist of secondary metabolites, and the secondary metabolites in different tissues of plants are significantly different 24 . And the accumulation of plant secondary metabolites produces a stress effect that can change the content of plant secondary metabolism 25 . Therefore, for the expression of genes related to the synthesis of secondary metabolites, the expression profiles may be difference in different tissues or under different stress treatments, and revealing their expression profiles is of great significance for the analysis of the genes. Previous studies in our laboratory showed that the content of oridonin, the main active component of I. rubescens, was significantly different in different tissues, and its content was also different in different environments. Therefore, the selection of reference genes in different tissues and under different stress treatments is of great significance for the study of genes in I. rubescens. Here, we have for the first time sought to evaluate the stability of 12 candidate reference genes for screening the most suitable reference genes for RT-qPCR data normalization in different sets of I. rubescens samples, including different tissues and different stress treatments. The selection of reference genes will provide a foundation for the benefit of future studies of gene expression in I. rubescens.

Results
Candidate reference genes selection and amplification specificity. Combined with the transcriptome data of I. rubescens, a total of 12 candidate reference genes (18S, TUB, ACT , Apt, HIS, UBQ, TUA , Cycl, GAPDH, eIF, PP2A, and EF-1α) that have been reported to be good potential reference genes in previously published papers were selected for validation as potentially suitable reference genes in I. rubescens from different tissues and under different abiotic stresses. The sequences of 12 candidate reference genes were obtained from transcriptome data of I. rubescens (Supplementary Table S1), and the primer sequences and descriptions of the candidate reference genes are shown in Table 1. To test the accuracy and specificity of the primers, agarose gel electrophoresis showed a single DNA band in the range of 110-200 bp, which was consistent with the predicted PCR products according to the various primer sets (Fig. 1A, Figure S1). Furthermore, the melting curves by RT-qPCR for the amplified products of all 12 candidate reference genes exhibited a single peak, corresponding to a specific melting temperature which ranged from 79.70 °C (UBQ) to 86.67 °C (HIS) (Fig. 1B). These two results indicated good specificity of all the primer pairs used in our study. In addition, the PCR efficiency of each candidate reference gene ranged from 1.8900 to 2.0528, among which 18S had the lowest amplification efficiency, PP2A had the highest amplification efficiency, and TUB showed the optimal efficiency (Table 1). In conclusion, these results suggest that the designed primers were suitable for RT-qPCR analysis of gene expression. Expression profiling of candidate reference genes in I. rubescens. The threshold cycle (Ct)-values were used to determine the expression levels of 12 tested candidate reference genes in all 24 samples. The Ct values are listed in Supplementary Table S2 and also exhibited as box-plot graphs ( Fig. 2A). The results showed that 18S were the highest expressed genes, with the lowest mean Ct-values of 10.55 in the different samples of I. rubescens. The mean Ct-values of others candidate reference gene ranged from 23.90 to 28.11, among which TUB showed the highest Ct-values of 28.81, and was the lowest expressed gene in all samples tested ( Fig. 2A). Obviously, although each candidate reference gene showed expression difference, 18S displayed high expression. This suggests that 18S, a commonly used traditional reference gene, can be selected as an internal reference gene for genes with high expression. Furthermore, the 2 −∆Ct method was used to analyze the expression stability of 12 tested candidate reference genes in all 24 samples. In this study, the fold change of not treated leaf sample was selected as a calibrator. The heatmap directly revealed the stability of gene expression in all samples (Fig. 2B). The expression levels of 18S, TUB, HIS, and eIF were stable, while EF-1α, TUA , and Cycl showed large variations in expression levels. www.nature.com/scientificreports/ Expression stability analysis of candidate reference genes in I. rubescens. In this study, four major computational algorithms were used to screen suitable candidate reference genes: the comparative ∆Ct method, BestKeeper, Normfinder, and geNorm. In the method of comparative ∆Ct, the stability of the candidate reference genes was assessed by the index of the mean standard deviation (STDEV). A low STDEV value of a candidate reference gene indicated that the gene was expressed stably in different samples. According to the results of this method, eIF was the most stable reference gene showing the lowest STDEV for all the samples and methyl jasmonate (MeJA) treatment of I. rubescens (Fig. 3A,F). GADPH was the most stable reference gene for the different tissues (DT) and salicylic acid (SA) treatment samples (Fig. 3B,E). Cycl and TUB exhibited the most stable reference gene for the NaCl and abscisic acid (ABA) treatment samples, respectively (Fig. 3D,G). UBQ was the most stable reference gene for dehydration treatment samples, ND set (combine NaCl and dehydration treatment samples) and SMA set (combine SA, MeJA, and ABA treatment samples) (Fig. 3C,H,I). The BestKeeper software tool is an Excel-based tool that uses the values of the standard deviation (STDEV) and the coefficient of variation (CV) of each candidate reference gene to estimate the gene expression stability in different samples. According to the standard values determined by BestKeeper, reference genes with STDEV > 1 are considered unstable and should be avoided 26 . In this study, there were four genes with STDEV < 1, including PP2A, HIS, TUB, and eIF, which showed that these genes possessed stable expression in total samples with SDvalues of 0.68, 0.73, 0.75, and 0.78, respectively. The most stable gene in total samples was PP2A, followed by HIS, TUB, and eIF, while the least stable gene was TUA (Fig. 4A). For the DT and SA samples, TUB was deemed the most stable gene (Fig. 4B,E). EF-1α was determined to be the most stable gene both in the MeJA samples and ABA samples (Fig. 4F,G). UBQ was deemed the most stable gene in the dehydration samples (Fig. 4C). PP2A was determined to be the most stable gene in the NaCl, ND and SMA samples by the BestKeeper software (Fig. 4D,H,I).
Normfinder is another algorithm for identifying the optimal reference gene among a set of candidates, and this algorithm can also calculate the variation between subgroups of the total sample set 27 . According the result of Normfinder, as the stability values (SV) decreased, the stability of the candidate reference genes increased. In the total samples of I. rubescens, eIF ranked first with an SV of 0.0.564, followed by Apt and UBQ with SV of 0.599 and 0.715, while, EF-1α ranked last with an SV of 1.974 (Fig. 5A). GADPH was the most suitable candidate genes for normalization with an SV of 0.245 in the DT samples (Fig. 5B). UBQ, Cycl, GADPH, eIF, TUB, and TUB ranked first and deemed as the most stable gene in the dehydration, NaCl, SA, MeJA, and ABA samples, respectively ( Fig. 5C-G). In addition, in the composite group of ND samples and SMA samples, UBQ ranked first with an SV of 0.292, and 0.217, respectively (Fig. 5H,I). The results of Normfinder were similar when compared with the results of the ∆Ct method.
GeNorm, an applet for Microsoft Excel, was used to calculate the gene expression stability measure M (M-value) for a reference gene as the average pairwise variation V for that gene with all other tested reference genes 28 . According to the algorithm used by GeNorm, the two most stable reference genes were determined based on the M-value. The M-value is inversely related to stability of the candidate reference genes, the lower Table 1. Gene description, primer sequences, amplicon length, and PCR efficiency for I. rubescens candidate reference genes.

Gene
Gene annotation arabidopsis ortholog www.nature.com/scientificreports/ the M-value is, the higher the stability will be, conversely, the higher the M-value is, the lower the stability will be. In our analysis, the most stable genes were ACT and GADPH with an M-value of 0.619 in   www.nature.com/scientificreports/ GeNorm was also used to revealed the optimal number of reference genes for normalization 21 . The pairwise variation between the normalization factors NF n and NF n+1 of two sequences (V n/n+1 ) was calculated. If V n/n+1 was more than 0.15, the number of genes used in normalization was n + 1, until V n/n+1 was less than 0.15, the number of genes used in normalization was n. In this study, the value of n was different in the four sample sets. For the total samples, six reference genes were needed for normalization since the V 6/7 -value (0.1451) was inferior to the 0.15 cut-off level. For ND samples, three reference genes were required for accurate RT-qPCR data normalization since the V 3/4 -value (0.1236) was < 0.15. And in the DT and SMA samples, four reference genes were required since the V 4/5 -values were 0.1423 and 0.1406, respectively (Table 3). In addition, for the remaining samples of dehydration, NaCl, SA, MeJA, and ABA samples, only two reference genes were required since the V 2/3 -values were all inferior to the 0.15 cut-off level.
Comprehensive ranking of the reference genes by RankAggreg. The ∆Ct method, BestKeeper, Normfinder, and geNorm exhibited different ranking results in their stability analysis of candidate reference gene expression in this study. For the total samples, the comparative ∆Ct method (Fig. 2) and NormFinder (Fig. 4) analysis ranked eIF as the most stable, while BestKeeper suggested PP2A was the most appropriate reference gene (Fig. 3), and GeNorm analysis showed ACT and GADPH to be the two most stable genes ( Table 2). The results of the other sets showed similar variable results. Thus, a comprehensive tool RankAggreg, known as a non-weighted unsupervised rank aggregation method, was used to evaluate the comprehensive ranking list of all genes based on the merged results of the ∆Ct, Bestkeeper, NormFinder and geNorm methods 23 . The results analysed by RankAggreg indicate that the most suitable reference genes tested for normalization in Total samples were eIF, Apt, and UBQ. In DT set, the best three candidate reference genes were GADPH, 18S, and eIF. UBQ, Apt, and HIS; Cycl, UBQ, and PP2A; GADPH, 18S, and eIF; eIF, UBQ, and PP2A; TUB, Cycl, and UBQ; were the best three candidate reference genes for the samples of Dehydration, NaCl, SA, MeJA, and ABA treatment, respectively. While in ND set, UBQ, HIS, and TUA were the three appropriate candidate reference genes. In addition, UBQ, eIF and Apt were top-ranked in the SMA set (Fig. 6).
Validation of the selected reference genes. Previously, only one ACT 7 or GADPH 29 gene was used as a reference gene in I. rubscences for data normalization. However, the results of our study show that ACT appears to be not most stably expressed reference gene in I. rubscences, especially in different tissues. GADPH appears stably expressed in different tissues, but not under other conditions (Fig. 3). Therefore, in order to validate the reliability of the three most stable candidate reference genes determined by RankAggreg (Fig. 6), all were either individually or in combination used as internal controls to analyze the expression pattern of 3-hydroxy-3-meth- The results were then compared to ACT or GADPH for data normalization. As shown in Fig. 7, the results showed the same expression patterns when using the combination of the two or three most stable genes or ACT or GADPH reference genes to normalize RT-qPCR data. However, the accumulation of HMGR transcripts was under-or over-estimated when the data were normalized using ACT or GAPDH as internal controls alone (Fig. 7). This result suggests that the determination of gene expression levels was obviously affected when unstable reference genes were used for normalization, and two or more reference genes contributed to the reliability of the results. Taken together, considering the cost and operation process, we recommend using GAPDH or combined with 18S in different I. rubscences tissues to normalize RT-qPCR data of gene expression. UBQ alone, or combined with Apt, and HIS; Cycl alone, or combined with UBQ, and PP2A; GADPH alone, or combined with 18S, and eIF; eIF alone, or combined with UBQ, and PP2A; TUB alone, or combined with Cycl, and UBQ; to normalize RT-qPCR data of gene expression. for the samples under Dehydration, NaCl, SA, MeJA, and ABA treatment, respectively.

Discussion
RT-qPCR is an accurate and sensitive technology used for the detection of gene expression levels. To obtain successful gene expression by RT-qPCR, it is crucial to select an appropriate reference gene. In the past few years, with the increase in research on medicinal plant genes, studies aimed at selecting appropriate reference genes as internal controls for RT-qPCR analysis have been reported in many medicinal plant species [14][15][16] . Furthermore, more evidence shows that the expression stability of reference genes may be altered across different tissues or under different experimental conditions. In Panax ginseng, GAPDH/30S RPS20, CYP/60S RPL13, and CYP/QCR were the optimum pair of reference genes in the roots, stems, and leaves 30 . In Metasequoia, ACT2, HIS, and TATA were stably expressed in different tissues. While EF-1α, HIS, and TATA were optimal for abscisic acid-response 21 .
In the Populus euphratica, the genes selected as optimal reference genes in different treatments were varied: RPL17 in ABA; EF-1α in cold; HIS in dehydration; GIIα in drought and short-duration salt; and TUB in long-duration salt 31 . However, as for I. rubscences, ACT 7 or GADPH 29 was frequently used as the internal controls in RT-qPCR without any form of stable reference gene examination, and this results in a lack of a certain degree of rigor. Thus, it was necessary to investigate the expression stability of candidate reference genes under different experimental conditions prior to using them for RT-qPCR data normalization. In our work, 12 genes commonly used as internal reference genes in previous articles were selected as candidate reference genes, and six statistical algorithms (ΔCt, Bestkeeper, NormFinder, geNorm, RefFinder, and RankAggreg) were used in gene expression studies of I. rubscences from different tissues or under different abiotic stresses ( Figure S2). We analyzed the candidate reference genes of each group under different treatment separately. Furthermore, in order to avoid the complexity of the experiment caused by too much subdivision of the six sets, we also set up three groups, including different tissue samples (DT set), NaCl and dehydration treatment samples (ND set), SA, MeJA, and ABA treatment samples (SMA set), so as to facilitate the precise application of internal reference genes under different environments conditions. In addition, in order to comprehensively compare the expression levels of genes in different samples, we also carried out screening and analysis of internal reference genes for all samples. As shown in Fig. 6, Comprehensive Ranking of the Reference Genes Table 3. Pairwise variation by geNorm determine the optimal number of reference genes.  www.nature.com/scientificreports/ by RankAggreg, which used the output data of ∆Ct, BestKeeper, NormFinder, and GeNorm methods to determine suitable reference genes in different sets, and the results indicate eIF, Apt, and UBQ were optimal for the whole sample analysis. GADPH, 18S, and eIF were stably expressed in different tissues. UBQ, Apt, and HIS; Cycl, UBQ, and PP2A; GADPH, 18S, and eIF; eIF, UBQ, and PP2A; TUB, Cycl, and UBQ; were the best three candidate reference genes for the samples of Dehydration, NaCl, SA, MeJA, and ABA treatment, respectively. In addition, UBQ, HIS, and PP2A were appropriate reference genes in the combination of NaCl or in the dehydration treatment, while UBQ, eIF and PP2A were optimal for the combination of SA, MeJA or ABA treatment. Therefore, the overall analysis results provide a basis for the appropriate selection of internal reference genes to study gene expression levels under different treatments in samples of I. rubscences.

Pairwise variations Total (values) DT (values) Dehydration (values) NaCl (values) SA (values) MeJA (values) ABA (values) ND (values) SMA (values)
The expression values of stable internal reference genes are not constant in different samples. As shown in Fig. 2, each candidate reference gene showed different Ct ranges in all samples, and the stability of each gene in different samples was clearly different. As a single software package may introduce bias, ∆Ct and geNorm, BestKeeper, and NormFinder were used for evaluation of candidate reference genes 32 . It can be seen that the different methods mentioned above yielded different ranking patterns. For the total samples, eIF, Apt, and UBQ were suitable genes by the ∆Ct (Fig. 3) and NormFinder (Fig. 5) methods, PP2A, while HIS and TUB were more stable than others by BestKeeper (Fig. 4), and ACT and GADPH were the most stable genes based on geNorm ( Table 2). For the analysis of other sets, there were similar phenomena; the results of different calculation methods showed different stability ranking lists of genes. The above results show that a single calculation method may introduce bias. Thus, RankAggeg, a comprehensive evaluation tools, were used to determine the final comprehensive ranking list of all genes 14 . As shown in Fig. 6, UBQ was relatively stable, especially in the ND set and SMA set, while it was less stable in the DT set. This inconsistency might be due to the different tissues used. UBQ is one of the most commonly used reference gene in many species, such as in NaCl-treated and ABA-treated Platycladus orientalis 33 , during somatic embryogenesis in longan tree 34 , and across all samples of Gossypium hirsutum 35 . However, our study demonstrated that UBQ was not stable in different tissues when compared with other reference genes, and GAPDH was the most stable in different tissues of I. rubscences (Figs. 3B, 5B). These results show that the stability of different internal reference genes is also affected by species. In Metasequoia, UBQ was also not stable in different tissues 21 . A traditional reference gene GAPDH was also one of the most used reference gene in many species. In Eremosparton songoricum, EsGAPDH was one of the most stable genes across multiple adult tissue samples 36 , as the same results were seen in Galium aparine with different tissues 8 , and at different developmental stages of Jatropha curcas 37 . These results demonstrate that different reference genes exhibit variable gene expression patterns in different species.
In most gene expression studies, only one reference gene was adopted for RT-qPCR data normalization. However, the use of one reference gene was more than likely to result in bias 38 . Therefore, the use of multiple reference genes has been adopted in more and more studies 39 . Additionally, some journals require at least two validated reference genes for normalization, such as Journal of Experiment Botany. In our study, GeNorm was www.nature.com/scientificreports/ used to select the optimal number of reference genes for normalization: for total samples, DT, ND and SMA sets, six, four, three and four reference genes were the optimal number, respectively (Table 3). In addition, to validate the selected reference genes, the three most stable genes were used as internal controls to normalize the RT-qPCR data from HMGR and were compared with ACT and GADPH normalization (Fig. 7). The results showed that the expression levels of the target gene HMGR were significantly changed when different reference genes were used for normalization, leading to unreliable experimental results. Our results support that a combination of multiple reference genes is preferred rather than a single internal reference gene. Similar results were also reported in other plant species using multiple reference gene normalizers 14 . Above all, the results acquired in this study will improve the accuracy of quantification of target gene expression with RT-qPCR analysis in I. rubscences and support the importance of selection of suitable reference genes for future gene expression studies in this species.  Analysis of stability of candidate reference genes. The Ct-values of the 12 reference genes in 24 different samples of I. rubescens were used to analyze the gene expression stability levels generated by the QuantStudio 5 software. Four algorithms were used to evaluate the stability of the eight candidate reference genes: the comparative ∆Ct method 17 , NormFinder 19 , BestKeeper 20 , geNorm 18 , and RankAggreg (https:// cran.r-proje ct. org/ web/ packa ges/ RankA ggreg/ index. html). The other statistical analyses software used in this study included SPSS v 19.0 (https:// www. ibm. com/ produ cts/ spss-stati stics), Origin8.0 (http:// www. origi n8. com/), and Heml 1.0 (http:// hemi. biocu ckoo. org/) 43  www.nature.com/scientificreports/ reference genes, the relative expression profiles of 3-hydroxy-3-methylglutaryl-CoA reductase (HMGR, accession number: MK344323.1), the first rate-limiting enzyme of mavalonic acid pathway, in samples of I. rubescens different tissues, and leaves under different stresses were analyzed and normalized to the most stable reference genes identified by RankAggreg. By comparison, ACT 7 or GADPH 29 reference genes were also used. The relative expression of HMGR was calculated using the 2 -∆∆Ct method 44 . Three technical replicates were performed for each biological sample. Values are reported as the mean ± SD from three independent experiments. Statistical analyses were performed using variance (ANOVA) followed by Duncan's new multiple range tests with SPSS version 19.0 (SPSS, Chicago, IL, USA). Different letters represented significant differences at p < 0.05.

Materials and methods
Ethical approval. We compliance of field studies for collection of plant or seed specimens with relevant guidelines by the IUCN Policy Statement on Research Involving Species at Risk of Extinction and the Convention on the Trade in Endangered Species of Wild Fauna and Flora. And we does not contain any studies with human or animal subjects.

Conclusions
In this study, we conducted the first comprehensive analysis of the selection of stable reference genes for RT-qPCR analysis of target gene expression in different tissues or under different treatments of I. rubscences. The combination of GAPDH, 18S and eIF was suitable for gene quantification in different tissues of I. rubscences.UBQ, Apt, and HIS; Cycl, UBQ, and PP2A; GADPH, 18S, and eIF; eIF, UBQ, and PP2A; TUB, Cycl, and UBQ; were the best three candidate reference genes for the samples of Dehydration, NaCl, SA, MeJA, and ABA treatment, respectively. While for the concatenated sets of ND (NaCl and dehydration) and SMA (SA, MeJA, and ABA), UBQ, HIS, and TUA; UBQ, eIF and Apt were the three appropriate candidate reference genes, respectively. The results acquired in this study will improve the accuracy of quantification of target gene expression with RT-qPCR analysis in I. rubscences and support the importance of selection of suitable reference genes for future gene expression studies in this species.

Data availability
The sequences data presented in this study are available under NCBI accession numbers OP006153-OP006163 and ON993894. And all the sequences data are also included as Supplementary Table s1.

Material availability
The seeds of I. rubescens were collected from Jiyuan, Henan province. And we statement that I. rubescens is not a protecting plant species in China, so we can harvest its seeds or other materials.