Systematic identification of reference genes for qRT-PCR of Ardisia kteniophylla A. DC under different experimental conditions and for anthocyanin-related genes studies

Ardisia kteniophylla A. DC, widely known as folk medicinal herb and ornamental plant, has been extensively investigated due to its unique leaf color, anti-cancer and other pharmacological activities. The quantitative real-time PCR (qRT-PCR) was an excellent tool for the analysis of gene expression with its high sensitivity and quantitative properties. Normalizing gene expression with stable reference genes was essential for qRT-PCR accuracy. In addition, no studies have yet been performed on the selection, verification and stability of internal reference genes suitable for A. kteniophylla, which has greatly hindered the biomolecular researches of this species. In this study, 29 candidate genes were successfully screened according to stable expression patterns of large-scale RNA seq data that from a variety of tissues and the roots of different growth stages in A. kteniophylla. The candidates were then further determined via qRT-PCR in various experimental samples, including MeJA, ABA, SA, NaCl, CuSO4, AgNO3, MnSO4, CoCl2, drought, low temperature, heat, waterlogging, wounding and oxidative stress. To assess the stability of the candidates, five commonly used strategies were employed: delta-CT, geNorm, BestKeeper, NormFinder, and the comprehensive tool RefFinder. In summary, UBC2 and UBA1 were found to be effective in accurately normalizing target gene expression in A. kteniophella regardless of experimental conditions, while PP2A-2 had the lowest stability. Additionally, to verify the reliability of the recommended reference genes under different colored leaf samples, we examined the expression patterns of six genes associated with anthocyanin synthesis and regulation. Our findings suggested that PAP1 and ANS3 may be involved in leaf color change in A. kteniphella. This study successfully identified the ideal reference gene for qRT-PCR analysis in A. kteniphella, providing a foundation for future research on gene function, particularly in the biosynthesis of anthocyanins.


Introduction
Ardisia kteniophylla A. DC, an evergreen shrub from the genus Ardisia (Primulaceae family), was mainly found in South China, North Vietnam, Laos and Thailand, and commonly known as a traditional folk medicine to treat rheumatism and injuries from muscles and fractures in China (Guan et al., 2016;Huang et al., 2017).It was first documented in the classical medical book "Sheng Cao Yao Xing Bei Yao" in the Qing dynasty (written in 1711), generally utilized by ethnic groups such as Yao, Zhuang, Dong, Jing, Li, Miao (Tian et al., 2023).Recently, phytochemical studies have found various chemical components that mainly contained triterpenoids, coumarins, flavonoids, phenols, quinones, sterols, and volatile oils, which active in a wide variety of biological processes such as anti-cancer, anti-HIV, anti-thrombotic antiinfammation, anti-oxidant, and cough expectorant activities (Mu et al., 2015;Mu et al., 2017;Mu et al., 2012;Tian et al., 2023).In addition, in the seedling stage, this plant is also served as an ornamental plant for its leaves are brightly colored, including red, semi red, and green (Tian et al., 2023;Zhao et al., 2021).Current research indicates that the biosynthetic pathways of anthocyanins have been identified in other plants, and found that anthocyanidin synthase (ANS) catalyzes leucocyanidin and catechin for conversion to cyanidin (Jun et al., 2018).Additionally, transcription factors PAP1 often regulate structural genes in the anthocyanin pathway in other plants, including strawberry (Zhang et al., 2020), carrot (Xu et al., 2020), tomato (Liu et al., 2016), and apple (An et al., 2020).Thus, ANS and PAP1 can regulate the expression of anthocyanin biosynthesis genes to influence leaf color.However, the molecular mechanisms involved in the color formation and transformation of leaves in A. kteniophylla are yet to be fully understood.Gene expression analysis has so far been a foundation for understanding gene function and molecular mechanisms in the development of organisms.
Gene expression analysis was critical approach for revealing gene function and enhancing understanding understanding the molecular mechanisms underlying diverse plant biological processes.At present, a variety of methods were available for assessing gene expression, including northern blotting, in situ hybridization, gene chips, sitedirected mutagenesis, and qRT-PCR (Gao et al., 2017;Wang et al., 2021;Yang et al., 2020).Among them, qRT-PCR, a particularly popular and effective method, has gained wide acceptance for its high sensitivity, accuracy, reliability, and specificity for detecting gene expression (Pfaffl et al., 2004).The reliability of gene expression measurements via qRT-PCR analysis, however, was strongly influenced by several factors, including RNA recovery and integrity, PCR amplification, PCR efficiency, primer design and the sample transcriptional activity (Pfaffl et al., 2004;Udvardi et al., 2008;Wang et al., 2021;Yang et al., 2020).Thus, to dispel these negative influences, technical differences from different samples and experimental environments were minimized by normalizing transcript levels of stable reference genes.The search for ideal reference genes that can be expressed under a variety of experimental conditions has been undertaken by several researchers.Unfortunately, most of these reference genes are not universally applicable to all plant species (Luo M. et al., 2018;Nolan et al., 2006;Udvardi et al., 2008;Wang et al., 2016).Therefore, a systematized approach to identify and evaluate multiple stable reference genes in A. kteniophylla was necessary to perform accurate transcript normalization experiments.
In recent years, Next Generation Sequencing (NGS) has been widely used in fields such as genome sequencing, RNA seq, ChIP sequencing, and has become a commonly used biotechnology for its high resolution, low cost, and high sensitivity (de Magalhães et al., 2010;Stone and Storchova, 2015).Most recently, RNA-seq data has been extensively utilized to identify new reference genes in plants, thereby serving as a potential source of candidate genes for various plant species, such as soybean (Yim et al., 2015), grape (Gonzaĺez-Agüero et al., 2013), Lycoris aurea (Ma et al., 2016), Alternanthera philoxeroides (Yin et al., 2021), Zingiber officinale (Li et al., 2022) and apple (Zhou et al., 2017).To screen the stable reference genes, four frequently statistical processing tools: geNorm (Vandesompele et al., 2002), NormFinder (Andersen et al., 2004), BestKeeper (Pfaffl et al., 2004), andDCt methods (Silver et al., 2006), were commonly used to evaluate the best suitable candidate, which were then analyzed from different perspectives.In addition, RefFinder, as a web-based comprehensive analytics tool, which can combine all four of these methods to compare and rank the optimal candidate reference genes (Xie et al., 2012).Consequently, these statistical software have been successfully applied to selecting the reference genes from a range of candidate genes in various plants, such as Codonopsis pilosula (Yang et al., 2020), Scutellaria baicalensis (Wang et al., 2021), Allium sativum (Wang et al., 2023), and so on.However, no systematic studies have been done on A. kteniophylla reference genes for qRT-PCR standardization, which greatly limits the research of functional genes in A. kteniophylla.
Herein, we used a comprehensive expression profile of largescale RNAseq data that from a variety of different tissues and the roots of different growth stages to select for the most stable genes expressed in A. kteniophylla.Subsequently, 29 candidate genes were successfully screened and their expression levels were further determined via qRT-PCR in 14 experimental samples.Finally, five commonly used statistical algorithms-delta-CT, geNorm, BestKeeper, NormFinder, and the comprehensive tool RefFinder were used to determine the stability of those candidates and rank them.Furthermore, to confirm the reliability and applicability of the reference genes, we analyzed the expression patterns of target genes related to anthocyanin regulation and synthesis.This work represents the initial systematic investigation of stable reference genes in A. kteniophylla, serving as a fundamental basis for future research on the molecular mechanisms underlying biological functions in A. kteniophylla, with a particular focus on the expression of genes related to anthocyanin synthesis.

Selecting reference genes for candidate analysis
An investigation was conducted to determine the gene expression patterns of A. kteniophylla utilizing RNA-seq data from various tissues, such as mainroot, fibrous root, root nodule, stem xylem, stem phloem, leaf, flower, and fruit, as well as roots of varying growth stages, including 1 month old, 3 months old, 6 months old, 1 year old, 2 years old, and 3 years old, which has yet to be published by our laboratory.A total of 42 data sets were used to calculate fragments per kilobase per million reads (FPKM), mean values (MV) and standard deviations (SD) for each gene.To identify genes with high expression levels and stability, we selected genes with MV values greater than 30 and ranked them based on their SD values from the smallest to the highest.The top 500 genes were further filtered on the basis of functional annotations to identify 29 candidates with commonly known functions.The heatmap were generated using TBtools software (Chen et al., 2020), and box plot using GraphPad Prism 7 software (Mitteer et al., 2018).Furthermore, the computation of the coefficient of variation (CV, FPKM/mean FPKM) and the maximum coefficient of variation (MFC; Maximum FPKM value/minimum FPKM value) was performed and utilized to generate dot plots for evaluating the expression stability of the 29 chosen genes, following the methodology described in the literature (Yin et al., 2021).

Plant materials and stress treatments
A. kteniophylla seeds were gathered from Nanxiong City, Guangdong Province, China.They were planted in pots (3 plants per pot) filled with soil mixture (peat: perlite = 1:1) and grown in the greenhouse of South China Botanical Garden, Chinese Academy of Sciences (Guangzhou, Guangdong Province, China, 112.3571°E, 23.1772°N).Healthy six-month-old seedlings were divided into different groups and subsequently receive abiotic stress and hormone treatments.

Total RNA extraction and cDNA synthesis
The Fastpure® Universal Plant Total RNA Isolation Kit (Vazyme, Nanjing, China) was utilized to extract total RNAs from each collected sample in accordance with the recommended protocols of the supplier.Subsequently, the integrity and quality of the total RNA samples were assessed through 1 g L −1 agarose gel electrophoresis, while the concentration and purity of RNAs were determined using a NanoDrop 2000c spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA).The RNAs were deemed acceptable if they exhibited an optical density of A260/280 ≈ 1.8-2.1 and A260/230 ≈ 1.7-2.2.The HiScript®II Reverse Reverse Transcriptase (Vazyme, Nanjing, China) was employed to synthesize single-strand cDNA from 1.0 µg of total RNA in a 20 µL volume system, following the manufacturer's protocols.Prior to use, all cDNA samples were diluted with nuclease-free water (1:20) and stored at -20°C.

Primer design and qRT-PCR
The gene-specific primers utilized for qRT-PCR (as presented in Table 1) were developed through the utilization of Primer Premier 5.0 (Premier, Inc., Toronto, Canada), with the following parameters: an amplicon length of 80-200 bp, primer length ranging from 18-22 bp, a melting temperature (Tm) of 58-61°C, and a GC content of 50-66%.The primers were subsequently synthesized by Tsingke Biotechnology Co., Ltd (Beijing, China).The primer pair specificity and amplicon size was tested by melting curve analysis and gel electrophoresis, PCR efficiency (E) and regression coefficient (R 2 ) were assessed using the diluted cDNA, as described earlier (Kudo et al., 2016).LightCycler® 480 System (Roche Molecular Systems, Germany) was used for qRT-PCR using LightCycler® 384 Plates.The reaction system were 10 mL in volume and contained 1 µL 20-fold diluted cDNA, 5 µL of 2 ×Vazyme qPCR SYBR Green Master Mix (Vazyme,Nanjing, China), 0.5µL of each forward and reverse primer and ddH 2 O 3.5 µL.The cycling conditions were initiated with an initial temperature of 95°C for 30 s, followed by 40 cycles of 95°C for 5 s and 60°C for 30 s, and concluded with a final melt-curve analysis.Analytical replicates of all qRT-PCR experiments were conducted three times.

Assessment of expression stability
To assess the stability of the 29 potential reference genes, five analytical tools, namely the Delta Ct method, BestKeeper, geNorm, Normfinder, and RefFinder, were utilized to determine the stabilities of each gene across varying experimental conditions (Andersen et al., 2004;Pfaffl et al., 2004;Rhinn et al., 2008;Silver et al., 2006;Vandesompele et al., 2002).Among them, the Delta Ct method enabled screening of appropriate references based on the SD value of Ct values, of which the smaller SD value indicate a higher gene stability (Andersen et al., 2004).The BestKeeper can determine the stabilization of the references via the coefficient of variation (CV) and standard deviation (SD) calculated based on the Ct values, with higher stability associated with lower CV and SD values (Pfaffl et al., 2004).GeNorm calculates an M-value based on the Ct value, with lower M-values represented greater stability (Rhinn et al., 2008).The Normfinder is an ANOVA-based model approach which calculates the SV value to assess the expression stability, where the lower SV value indicates a higher expression stability (Silver et al., 2006).As a final step, RefFinder, as a comprehensive tool (Vandesompele et al., 2002), computed the geometric average of the weights assigned to each gene in accordance with the four methods worked out previously.

Stability verification of reference genes
The leaves of A. kteniophylla seedling show different colors including red, semi-red and green for the different anthocyanin content.Here, to assess the stability of the potential suitable reference genes, three sample sets were prepared, including samples of red, semi-red, and green leaf of A. kteniophylla seedling.Each sample was immediately frozen in liquid nitrogen and stored at a temperature of −80°C after being collected.The 2 - △△Ct method was used to determine the relative expression of investigated genes (Kumar et al., 2011).For the quantification of double reference genes, △△Ct was calculated as [((Ct gene -Ct reference 1) experiment -(Ct gene -Ct reference1) control)/2 + ((Ct gene -Ct reference 2) experiment-(Ct gene -Ct reference2) control)/2] (Gou et al., 2020).

Results
RNA-seq data selection for candidate reference genes in A. kteniophylla In order to select suitable candidate references, RNA-seq was conducted to assess gene expression levels in A. kteniophylla, which were subsequently normalized by FPKM values.The FPKM value was then utilized to calculate the mean value (MV) and standard deviation (SD).The candidate genes were arranged in ascending order based on their SD, and 29 genes with MV exceeding 30 and possessing functional annotations were chosen (Figure 1A and Supplementary Table 1).In order to conduct a more comprehensive assessment of the expression stability of the genes in question, the Mean Fold Change (MFC) and Coefficient of Variation (CV) were computed.As depicted in Figure 1B and Supplementary Table 1, the MFC and CV values for the 29 candidate genes were found to be less than 6.3 and 0.5, respectively.Taken together, the aforementioned four evaluation metrics (i.e., MV, SD, MFC, and CV values) collectively suggest that the 29 selected candidate genes are suitable for subsequent screening as ideal reference genes.To further investigate the stability of expression of 29 chosen genes, heat maps (Figure 1C and Supplementary Table 2) and box plots (Figure 1D) were generated using RNA-seq data.The results demonstrate that the expression of these genes remained remarkably stable and consistent across various growth stages and tissues.Therefore, it can be inferred that the selected genes are dependable and trustworthy.

Primer specifcity and PCR efficiency analysis of reference gene
In order to identify optimal reference genes for gene expression profiling of A. kteniophylla, primers for 29 commonly utilized reference genes in qRT-PCR were designed using Primer Premier 5 (Supplementary Table 3).The specificity of these primer pairs was subsequently confirmed via electrophoresis on agarose gels, which demonstrated that each reference gene was specifically amplified and presented as a singular band, as anticipated (Supplementary Figure 1).Furthermore, melting curve analysis revealed that each amplicon exhibited only one peak (Supplementary Figure 2).The findings indicate that the primers exhibited a high degree of specificity and were suitable for gene expression analysis via qRT-PCR.The PCR amplification efficiency (E) was assessed and the results, presented in Supplementary Figure 3 and Table 1, revealed that the E values for each primer pair ranged from 94.08% (ADF) to 104.87% (ATP).Additionally, the regression coefficient (R 2 ) values varied between 0.9836 for a-TUB-3 and 0.9999 for EIF-1A.

Analyzing reference gene stability comprehensively
Genes that make up the ideal reference set should be capable of stably expressing themselves in different conditions, hormone treatments, abiotic stress, and organ types (Mughal et al., 2018).To evaluate the consistency of the expression of the 29 potential genes and determine the most suitable reference gene for normalizing expression in A. kteniophylla, their expression levels RNA-seq data selection for candidate reference genes in A. kteniophylla.The transcriptome's candidate genes are presented through dot plots that display their mean (FPKM) and standard deviation (SD) values (A), as well as their coefficient of variation (CV) and maximum fold change (MFC) (B).(C) The expression levels of 29 candidate genes across various tissues and growth years are depicted in the heat map.The bar charts represent the expression levels of corresponding genes across various stages and tissues.A higher bar chart indicates a higher level of expression, while a lower bar chart indicates a lower level of expression.(D) Overview of the expression levels of 29 potential reference genes is provided through box and whisker plot graphs that exhibit the log2 (FPKM) values of each selected gene in multiple tissues.The 25th and 75th percentiles and medians are represented by black boxes and lines, respectively.Whisker caps represent the maximum and minimum values.
were assessed across various treatment conditions(MeJA, ABA, SA, NaCl, Cu, Ag, Mn, Co,PEG-induced drought low temperature, heat, waterlogging, wounding, and oxidative stress) via qRT-PCR (Supplementary Table 4).Among candidate reference genes in plant species, five software tools, namely delta-CT, BestKeeper, NormFinder, geNorm, and RefFinder were usually selected to identify appropriate reference genes.Using the delta CT method, SD values were calculated to rank the stability of the candidates, and genes with lower SD values were selected as more stable reference genes (Duan et al., 2017).According to Figure 2B, the results of delta CT analysis showed that UBC2 (2.01), SAM5 (2.04), and UBA1 (2.04) were the most stable reference genes, while a-TUB-3 (3.64), HSP70 (3.83) and PP2A-2 (9.26) were the most unstable.BestKeeper evaluate the stability of candidates by calculating the CV and SD of the Ct (cycle threshold) values of the candidate genes, and the genes with the lowest CV and SD (CV ± SD) considered to be the most stable.As depicted in Figure 2C, the BestKeeper analysis suggested that UBA1 (3.49 ± 0.75), CYP1 (4.11 ± 0.80) and PP2A-1 (3.77 ± 0.81) were the most stable reference genes, while PP2A-2 (17.43 ± 4.78), RPL13 (12.94 ± 2.75) and a-TUB-4 (7.91 ± 2.36) were found to be the least stable.GeNorm ranked the candidates according to expression stability M values, which are the average paired variation of a certain gene with all other candidate genes, with lower M values indicating more stable expression (Vandesompele et al., 2002).As can be seen in Figure 2E, the UBC2 and UBC9 genes were found to be the most stable reference genes under various treatments, with an M value of 0.56.On the contrary, HSP70 and PP2A-2 were the unstable gene, with M values of 2.28 and2.76,respectively.NormFinder ranked the potential genes according to their minimal inter-and intra-group expression variation, with lower variation values indicating higher stability of the reference genes (Andersen et al., 2004).This study examined that UBC2 and UBA1 genes were the most stable with the lowest ranking values of 0.74 and 0.79, while PP2A-2 and HSP70 were the least stable, with the ranking values at 9.10 and 3.33, respectively (Figure 2D).These results indicated that NormFinder and geNorm gave almost identical stability orders for the candidate genes, although minor differences were observed.In conjunction with the aforementioned four methods, the all-encompassing analysis tool, Reffinder, computes the geometric mean value of the weights assigned to the potential references, culminating in a conclusive and comprehensive ranking of said references.Through this thorough analysis, UBC2 was identified as the most stable reference gene, with UBA1 and PP2A-1 following closely behind (refer to Figure 2A).
To further analyze the selection of ideal references in A. kteniophylla under different specific conditions, all the samples treated were categorized into three groups based on exposure to heavy metal stress (Cu, Ag, Mn and Co), hormone treatment (MeJA, ABA and SA) and abiotic stress (salt, drought, cold, heat, waterlogging, injury and oxidation stress).The results depicted in Figure 3 indicate that UBA1 and UBC2 exhibit the highest stability as reference genes for heavy metal stress, while PP2A-2 displays the least stability (Figure 3A).As for hormone treatment and abiotic stress treatment, UBC2 was shown to be the best gene (Figures 3B,  C), and PP2A-2 was confirmed as the least stable (Figures 3B, C).
The expression stability ranking of 29 selected reference genes across samples determined vis different evaluation algorithms: (A) RefFinder, (B) delta-CT, (C) BestKeeper, (D) NormFinder, and (E) GeNorm.Wang et al. 10.3389/fpls.2023.1284007Frontiers in Plant Science frontiersin.org Application of reference genes in the anthocyanin-related genes of A. kteniophylla As shown in Figures 4A-C, A. kteniophylla was currently used as the ornamental plant due to its bright leaf color (including red, semi red, and green) at seedling stage, which was mainly caused by the difference in anthocyanin content.At present, numerous studies have suggested that anthocyanins were widely distributed pigments in plants, and their biosynthetic pathways and regulatory mechanisms have been widely studied in model plants and ornamental plants (Zhao and Tao, 2015).Based on comprehensive analysis via RefFinder, the top two genes (UBC2 and UBA1) in A. kteniophylla were the most suitable for gene expression studies, while PP2A-2 was not recommended as the reference gene.To determine the reliability of the reference genes selected for normalization, we selected six key genes related to anthocyanin synthesis and regulation as target genes (Supplementary Tables 3 and 5) (including two transcription factors: PAP1, TT8 and four anthocyanin synthase genes: ANS1, ANS2, ANS3 and ANS4) for normalization with the most stable reference gene or combination as well as the least stable reference gene.
According to Figures 4D-I, the expression of these genes tested showed similar trends when the reference genes with the most stable expression (UBC2 and UBA1) were used alone or in combination under different color leaf samples, while the expression pattern of the most unstable expression reference gene (PP2A-2) obtained was significantly different from those of the stable references used alone or in combination.For example, when normalized with the most suitable reference genes UBC2 or UBA1 and their combination, the expression levels of PAP1 was highest in the red leaves, followed by semi red leaves and lowest in green leaves.However, the expression patterns of PAP1 were quite different from those of the optimal reference gene when normalized with the unstable gene PP2A-2, which PAP1 gene showed the highest expression level in the semi red leaves, followed by green leaves and lowest in red leaves (Figure 4D).Therefore, this results in biased results when inappropriate reference Stability rankings of 29 selected reference genes across (A) stress treatments, (B) heavy metal stress, and (C) hormone treatments.Wang et al. 10.3389/fpls.2023.1284007Frontiers in Plant Science frontiersin.orggenes are selected.In addition to confirming the significance of selecting reference genes prior qRT-PCR, these results demonstrated that our selected genes could serve as reliable and feasible references.Furthermore, this study found that among the six key genes that may be potentially related to anthocyanin synthesis in different leaf color samples, the transcription factor PAP1 and structural genes ANS2, ANS3 and ANS4 had the highest expression levels in red leaves (Figure 4).Regardless of which reference genes and their combination (UBC2, UBA1 and UBC2 & UBA1).Meanwhile, the expression patterns of PAP1 and ANS3 were similar: red leaf > semired leaf > green leaf, indicating that these two genes may play a significant role in the formation and change of leaf color in A. kteniophylla.In contrast, the transcription factor TT8 and structural genes ANS1 may not have performed.

Discussion
A. kteniophylla, as a significant medicinal and ornamental plant, has drawn more attention of researchers due to its showy leaf color, anticancer and other pharmacological activities.In recent years, high-throughput sequencing has become increasingly popular, RNA sequencing based on genome sequence has been widely used to detect gene expression, which provides a basis for the study of molecular biology studies of A. kteniophylla.Additionally, qRT-PCR was a method of high sensitivity, accuracy, and throughput for analyzing gene expression, required the ideal reference genes as a prerequisite to assure the accuracy and dependability of the results.So far, the development of stable reference genes in A. kteniophylla has received little attention, despite evaluating and validating reference genes in a variety of plants.Furthermore, the establishment of reference genes for expression in A. kteniophylla using RNA-seq datasets has not been previously reported.

Identification of potential reference genes using RNA seq data
Recently, RNA-seq has been used to identify novel and potential reference genes for several non-model plants using expression data and assembly sequences (Li et al., 2022;Yin et al., 2021).In this study, 42 RNA-Seq datasets (unpublished) were analyzed to measure the expression levels of different genes in A. kteniophylla various tissues and growth years.In the present study, appropriate reference genes were selected using expression profile data to quantify genes in A. kteniophylla, with special attention to the application of reference genes in the study of anthocyanin synthesis and regulation related genes during the color change process of A. kteniophylla leaves.Mughal et al. (Mughal et al., 2018).proposed that the ideal reference gene suitable for qRT-PCR must meet the following criteria: the gene must have medium to high expression level, generally higher than the basic cell background; Meanwhile, the gene must be stably expressed across different tissues, growth stages and stress conditions, and the expression level must maintain low variation.Therefore, in order to select potential candidate reference genes, the expression level of A. kteniophylla genes were standardized to FPKM values using RNA seq data.These values were subsequently utilized to compute four crucial indicators, namely SD (FPKM standard deviation), MV (FPKM mean), MFC (maximum multiple change of FPKM), and CV (coefficient of variation), for assessing candidate expression stability (Figure 1).Meanwhile, according to numerous previous studies (Wang et al., 2021;Yang et al., 2020;Zhong et al., 2022), we tried to select some commonly used housekeeping genes as candidate genes.Finally, the 29 selected candidate genes were suitable for further screening of ideal reference genes (Table 1).

Reference genes UBC2 and UBA1 were validated most stable
The standard qRT-PCR method was utilized to evaluate the amplification efficiency and primer specificity of all 29 candidate genes, as demonstrated in Supplementary Figure 3 and Table 1.All 29 pairs of primers were found to meet the criteria for qRT-PCR experiments.To identify the optimal genes for gene expression normalization in A. kteniophylla, the expression levels of the 29 candidates were examined under various treatment conditions using qRT-PCR.In this study, we conducted a systematic evaluation and relatively stable reference genes using four most frequently used algorithms (Delta CT, BestKeeper, geNorm, and NormFinder) to improve the reliability of reference gene evaluation (Figure 2).The ranking order of the 29 candidates determined by the results of all four algorithms was reasonable but not completely consistent, mainly caused by different statistical algorithms, which have been similar findings in previous studies (Abbas et al., 2021;Li J. et al., 2017;Zhao et al., 2016;Zhong et al., 2022).According to previous studies, these algorithms identified mostly the same most unstable genes, but the order in which the most stable genes were sorted was different (Duan et al., 2017).Similarly, PP2A-2 was considered the least stable genes among the four algorithms, with slightly different sequences of stable genes in this study.Therefore, we utilized the comprehensive tool RefFinder to perform final comprehensive sorting on the results of the four algorithms to improve the reliability of the results, and found that among the top stability reference genes according to RefFinder, UBC2 was the top reference gene, followed by UBA1 and PP2A-1 (Figure 2A).Furthermore, the samples were analyzed further to select ideal reference genes in A. kteniophylla under three groupings of different specific conditions.Figure 3 demonstrates that UBA1 and UBC2 was the most appropriate reference gene for heavy metal stress, and UBC2 was confirmed as the most stable gene for hormone treatment and abiotic stress treatments, while PP2A-2 was determined to be the least stable gene under all the above conditions.The finding was similar to previous research (Yin et al., 2021).The results showed that under different experimental conditions, UBC2 and UBA1, with a particular emphasis on UBC2, were the optimal reference genes for gene expression analysis in A. kteniophylla across diverse experimental conditions.
Our study indicated UBC2 was the most suitable reference gene across most samples except for heavy metal stress samples (Figures 2 and 3) among the 29 candidates, which may be helpful when analyzing target genes in A. kteniophylla using qRT-PCR.UBC2, as an essential component of the ubiquitin/proteasome pathway gene, has a crucial role in ubiquitination of proteins.At present, UBC2 performed as an the optimal reference gene in Cicer arietinum (Castro et al., 2012), Trifolium pratense (Mehdi Khanlou and Van Bockstaele, 2012), Cajanus cajan (Sinha et al., 2015), Corchorus capsularis (Lin et al., 2023), and Glycyrrhiza glabra L (Li et al., 2020) for its highly conservative character, while exhibited low stability in Fucus distichus (Linardićand Braybrook, 2021).Furthermore, UBA1 was the most stably expressed gene in heavy metal stress samples, and the second stably expressed gene in all samples.UBA1, as the first enzyme in the ubiquitination pathway, which encodes instructions for the production of the ubiquitinactivating enzyme E1 (Panchamia et al., 2020).It has been reported that this enzyme was necessary for the ubiquitin-proteasome system, which was one of the most important cellular processes for protein modification in plants, and involved in signaling relating to host defense (Hatfield and Vierstra, 1992;Shirsekar et al., 2010).In contrast, PP2A-2 exhibits particularly excellent instability in all samples, despite its stable expression in other species such as Nicotiana tabacum (Schmidt and Delaney, 2010), Elymus sibiricus (Zhang et al., 2019), and Oenanthe javanica (Jiang et al., 2014).This discovery was similar to the results of previous studies in Momordica charantia (Wang et al., 2019).Previous studies have found that homologous genes exhibit similar sequences and functions, and commonly used as reference genes for gene transcription analysis (Zhou et al., 2019).As an example, as appropriate reference gene for Fortunella crassifolia, the Actin gene ACT6, ACT7, and ACT8 were selected (Pfaffl et al., 2004).However, other research results have found opposite conclusions, such as UBC19, UBC22, and UBC29 as homologous genes in Isatis indigotica, with expression levels and relatively poor stability (Li T. et al., 2017).Here, We found that the homologous ubiquitin-conjugating enzyme genes (UBC2, UBC5, UBC9, and UBC35) possess completely distinct expression levels and stability characteristics (Figure 2).Due to this, homologous genes may show quite different expression levels and stability despite having similar sequences and functions.
UBC2 and UBA1 performed well in the quantification of anthocyanin related genes in A. kteniophylla In A. kteniophella, anthocyanins were the cause of the formation of different colors in seedling leaves (Figures 4A-C).Anthocyanins, as a class of flavonoids, impart color to different types of leaves, fruits, flowers, and many other tissues to prevent various biological and abiotic damage, as well as damage caused by ultraviolet or visible light (Noda et al., 2017).Additionally, anthocyanins are crucial to the cultivation of flowers, nutrition, healthcare, and the provision of medical care in human society (Yoshida et al., 2009).Molecular mechanisms underlying leaf color formation in A. kteniophella seedlings, however, remain poorly understood.The biosynthesis pathway of anthocyanins has been well studied in various plants, and found that anthocyanin synthase (ANS), as a key enzyme gene in anthocyanin synthesis pathway, catalyzed the conversion of leucoanthocyanidins to anthocyanins (Jun et al., 2018;Luo H. et al., 2018).Transcription factors bHLH and R2R3-MYB were typically involved in regulating structural genes in the flavonoid pathway, and it is particularly noteworthy that PAP1 (R2R3-MYB family) and TT8 (bHLH family) genes can positively regulate structural genes in the anthocyanin pathway and promote anthocyanin synthesis and accumulation (Fraser et al., 2013;Jaakola, 2013;Zhang et al., 2020).To ascertain the accuracy and reliability of the selected reference genes, a normalization procedure was conducted on the relative expression levels of target genes (namely, PAP1, TT8, ANS1, ANS2, ANS3, and ANS4) in three distinct colored leaves of A. kteniophella seedlings.This normalization involved the use of the two most stable reference genes, their combinations, and the least stable reference gene.The findings revealed that the employment of the most stably expressed reference genes (UBC2 and UBA1) either individually or in conjunction facilitated the clear expression patterns of all target genes across the various colored leaves samples.
Even though the expression trends of each stable reference genes were consistent with those of stable reference gene combinations, there were still slight differences in expression levels.Nevertheless, the expression profile of the inappropriate reference gene (PP2A-2) showed complete inconsistency with the results of the stable gene and significant changes occurred.This clearly demonstrated the impact of using unstable reference genes on qRT-PCR results, and also indicated that the relative expression level of target genes varies with the selection and quantity of references.Previously, similar results have been studied in multiple species, including Codonopsis pilosula (Yang et al., 2020), Scutellaria baicalensis (Wang et al., 2021), Aegilops tauschii (Abbas et al., 2021), and Dioscorea opposita (Zhao et al., 2016).In addition, we also found that among all six key genes, the expression levels of transcription factor PAP1 and structural genes ANS2, ANS3, and ANS4 were highest in red leaves (Figure 4), especially the expression profiles of PAP1 and ANS3 were red leaves>semi red leaves>green leaves, suggesting that these two genes may play a significant role in the anthocyanin biosynthesis pathway in embryonic leaves.Taking these findings together, they highlight the fundamental significance of selecting ideal reference genes before performing qRT-PCR.

Conclusion
Here, we systematically identified and screened ideal reference genes under diverse experimental conditions for the first time, which were used for qRT-PCR analysis of target gene expression in Ardisia genus plant, A. kteniophella.From our results, two stably expressing genes, UBC2 and UBA1, were identified as appropriate for accurately normalizing target gene expression in A. kteniophella across all experimental conditions.In addition, UBC2 and UBA1 have also been successfully used for gene expression analysis in the biological study of color difference in the A. kteniophella leaves, and found that PAP1 and ANS3 may play a critical role in the formation and change of leaf color in A. kteniphella for the first time.These findings not only make a significant contribution to functional genomics analysis for A. kteniphella, but they also lay a good foundation for the study of molecular mechanism of biological function in A. kteniphella in the future.Besides, the stable references of this study provide reference guidelines for other species that have not validated qRT-PCR analysis reference genes.

Funding
The author(s) declare financial support was received for research, authorship, and/or publication of this article.This work was supported by Guangdong Province Plant Ex-situ Protection System Planning and National Key Protected Plant Ex-situ Protection Expansion (Research on Key Techniques for Rescue Protection of Three Key Protected Plants) (E336050041), and Protection and construction of rare and endemic plants in highlevel special parks (Research on Key technologies for rescue protection of Ardisia kteniophylla) (E236120011).
FIGURE 4 Reference genes identification and application in the anthocyanin-related genes of A. kteniophylla.(A) The leaf phenotype of red, semi red and green.(Scale bar: 2cm).(B, C) Relative anthocyanin content in red, semi-red, and green leaves of seedlings.(D-I) Relative expression levels of target gene PAP1 (D), TT8 (E), ANS1 (F), ANS2 (G), ANS3 (H), and ANS3 (I) in different leaf colors were normalized with the most stable reference gene or combination and the most unstable gene.Statistically significant differences between samples are marked with different letters (P<0.05,Student's ttest).

TABLE 1
Details of candidate reference genes and validation genes for qRT-PCR normalization in A. kteniophylla.