Revista Internacional de Andrología. 2026; 24(3): 29-45. doi: 10.22514/j.androl.2026.028
Original Research

Selection and validation of optimal housekeeping gene through RNA-sequence and qRT-PCR analysis in rat testis tissues and Sertoli cells under mangancism and Coridius chinensis extracts rescue
Selección y validación del gen constitutivo óptimo mediante análisis de transcriptoma (RNA-seq) y qRT-PCR en tejidos testiculares y células de Sertoli de rata bajo manganismo y tratamiento con extractos de Coridius chinensis

Jingjia Wu1,†, Pan Tao1,†, Zhengyi Zhang1, Dengfeng Luo1, Huixian Lou1, Chuanbo Zhu1, Xiaolan Yan1, Yu Hu1,2, Xiaowei Sun3, Xiaohui Hou1,*,, Fei Gao4, Nanjian Luo1,*,

1School of Preclinical Medicine, Zunyi Medical University, 563000 Zunyi, Guizhou, China

2Department of Reproductive Medicine, Affiliated Hospital of Zunyi Medical University, 563000 Zunyi, Guizhou, China

3Yunnan Animal Husbandry Station, 650000 Kunming, Yunnan, China

4State Key Laboratory of Stem Cell and Reproductive Biology, Institute of Zoology, Chinese Academy of Sciences, 100101 Beijing, China

*Corresponding Author(s):houxh1@zmu.edu.cn (Xiaohui Hou); luonj@zmu.edu.cn (Nanjian Luo)

† These authors contributed equally.

History Submitted: 14 May 2025 | Accepted: 30 July 2025 | Published: 30 September 2026
Copyright:  ©2026 The Author(s). Published by MRE Press.
This is an open access article under the CC BY 4.0 license (https://creativecommons.org/licenses/by/4.0/).

Collapse table of contents

Abstract

Background: Testicular tissues and cells undergoing injury or therapeutic interventions often exhibit significant changes in functional and housekeeping gene (HKG) expression. Using unstable HKGs can cause inaccurate quantification of target gene expression in quantitative real-time polymerase chain reaction (qRT-PCR) analyses. Therefore, identifying reliable HKGs is critical for accurately assessing functional gene expression in studies of testicular injury and therapeutic interventions. Methods: In this study, 17 candidate HKGs were selected based on RNA-seq data and previous reports. The cycle threshold (Ct) values of these genes were determined by qRT-PCR in testicular tissue and Sertoli cells from the control (Ctrl) group, manganese ion (Mn2+) group, and Mn2+ combined with Coridius chinensis extract (Mn2+ + CcE) group. Expression stability was evaluated using four algorithms (geNorm, Normfinder, Bestkeeper, and Delta Ct), and comprehensively assessed by RefFinder and ComprFinder. The three most stable HKGs were validated by normalization to specific marker genes for spermatogonia and Sertoli cells (CASPASE3 (Cysteine-dependent aspartate-specific protease-3), DDX4 (DEAD box polypeptide 4), SOX9 (SRY-box transcription factor 9), and OCCLUDIN (Occludin)). Results: PPIA (Peptidylprolyl Isomerase A) exhibited the lowest Ct values in both testicular tissues and related cells. In testicular tissues, YBX1 (Y-box binding protein 1) and PPIA had the lowest Ct values across different treatment groups, while UBB (Ubiquitin B) and TUBULIN (tubulin) exhibited the highest Ct values. In Sertoli cells, PPIA and ACTB (Actin Beta) were identified as the most stable HKGs under various culture conditions. However, classical HKGs, such as TUBULIN and GAPDH (Glyceraldehyde-3-Phosphate Dehydrogenase) displayed unstable expression, indicating their unsuitability for target gene normalization. Conclusions: Our study demonstrates that PPIA may serve as the optimal HKG in rat testicular tissues and Sertoli cells under various pathological conditions, providing valuable reference data for functional gene expression studies related to testicular diseases in mammals.

Resumen
Antecedentes: Los tejidos y células testiculares sometidos a lesiones o intervenciones terapéuticas suelen presentar cambios significativos en la expresión de genes clave, incluidos los genes funcionales y los genes constitutivos (HKGs, por sus siglas en inglés). El uso de HKGs inestables puede conducir a una cuantificación errónea de la expresión génica objetivo en los análisis de PCR cuantitativa en tiempo real (qRT-PCR). Por lo tanto, la identificación de HKGs confiables es crucial para evaluar con precisión la expresión génica funcional en estudios relacionados con lesiones testiculares e intervenciones terapéuticas. Métodos: En este estudio se seleccionaron inicialmente 17 genes constitutivos candidatos basándose en datos de secuenciación de ARN (RNA-seq) y en informes previos. Los valores del umbral de ciclo (Ct) de estos genes se determinaron mediante qRT-PCR en tejido testicular y células de Sertoli de los grupos control (Ctrl), con ion manganeso (Mn2+) y con Mn2+ combinado con extracto de Coridius chinensis (Mn2+ + CcE). La estabilidad de expresión de estos genes candidatos fue evaluada mediante cuatro algoritmos (geNorm, NormFinder, BestKeeper y Delta Ct), y analizada de forma integral con RefFinder y ComprFinder. Los tres HKGs más estables fueron validados mediante su normalización con respecto a genes marcadores específicos de espermatogonias y células de Sertoli (CASPASE3 (Cysteine-dependent aspartate-specific protease- 3), DDX4 (DEAD box polypeptide 4), SOX9 (SRY-box transcription factor 9) y OCCLUDIN (Occludin)). Resultados: PPIA (Peptidylprolyl Isomerase A) mostró los valores de Ct más bajos tanto en los tejidos testiculares como en las células asociadas. En los tejidos testiculares, YBX1 (Y-box binding protein 1) y PPIA como los valores de Ct más bajos entre los distintos grupos de tratamiento, mientras que UBB (Ubiquitin B) y TUBULINA (tubulin) presentaron los valores de Ct más elevados. En las células de Sertoli, PPIA y ACTB (Actin Beta) fueron identificados como los HKGs más estables bajo diversas condiciones de cultivo según la evaluación integral. Sin embargo, HKGs clásicos como TUBULINA y GAPDH (Glyceraldehyde-3-Phosphate Dehydrogenase) mostraron una expresión inestable, lo que indica que no son adecuados para la normalización de genes diana. Conclusiones: Nuestro estudio demuestra que PPIA podría considerarse el único HKG óptimo en tejidos testiculares y células de Sertoli de rata bajo diversas condiciones patológicas, proporcionando datos de referencia valiosos para la selección de HKGs en estudios de expresión génica funcional relacionados con enfermedades testiculares y tratamientos en mamíferos.

Keywords:Testis;Sertoli cell;Transcriptome;Housekeeping gene;Manganism;Coridius chinensis
Palabras Clave:Testículo;Célula de Sertoli;Transcriptoma;Gen constitutivo;Manganismo;Coridius chinensis
PDF(16.26 MB)|EndNote (RIS)|BibTeX|RefMan|RefWorks

Cite this article

Jingjia Wu, Pan Tao, Zhengyi Zhang, Dengfeng Luo, Huixian Lou, Chuanbo Zhu, Xiaolan Yan, Yu Hu, Xiaowei Sun, Xiaohui Hou, Fei Gao, Nanjian Luo. Selection and validation of optimal housekeeping gene through RNA-sequence and qRT-PCR analysis in rat testis tissues and Sertoli cells under mangancism and Coridius chinensis extracts rescue. Revista Internacional de Andrología. 2026; 24(3): 29-45. doi: 10.22514/j.androl.2026.028

1. Introduction

The testis is a crucial organ in the mammalian male reproductive system, playing vital roles in spermatogenesis as well as androgen synthesis and secretion [1, 2, 3]. Spermatogenesis primarily depends on meiotic differentiation of spermatogonial stem cells (SSCs) located in the spermatogenic epithelium of the testis [3]. Sertoli cells, positioned between the SSCs and meiotic spermatogenic cells, provide nutritional support to various stages of sperm cells [4], and contribute to the formation of the blood-testis barrier through tight junctions between adjacent Sertoli cells [5, 6]. Additionally, Leydig cells, located outside the seminiferous tubules, produce androgens that maintain secondary sexual characteristics and positively regulate spermatogenesis [4, 7]. However, with global urbanization and industrialization, heavy metals, such as lead (Pb), cadmium (Cd), and manganese (Mn) increasingly, accumulate in groundwater, soil, surface water and ecosystems, causing testicular tissue damage and dyszoospermia in both humans and animals [8, 9, 10]. Several studies from our group have shown that Coridius chinensis (C. chinensis) extract (CcE) significantly ameliorates testicular disorders induced by divalent manganese ion (Mn2+), acting through antioxidant, antitoxic, and antiapoptotic mechanisms [11, 12]. These physiological and cellular recoveries were attributed to changes in the expression of various functional genes regulated by CcE. Therefore, accurately determining target gene expression levels is essential for understanding how CcE alleviates spermatogenic disorders induced by Mn2+ exposure.

Quantitative real-time polymerase chain reaction (qRT-PCR) is a commonly used technique for assessing gene expression across diverse biological samples. To accurately quantify relative gene expression using qRT-PCR data, stable housekeeping genes (HKGs) are employed for normalization. Stable expression of HKGs under varying physiological and pathological conditions is crucial. Numerous software programs, including geNorm [13], Normfinder [14], BestKeeper [15], and Delta Ct [16], have been developed to evaluate HKG stability through statistical analyses of cycle threshold (Ct) values. Additionally, RefFinder and ComprFinder integrate these four algorithms to comprehensively evaluate HKG stability [17, 18]. Meanwhile, transcriptome sequencing provides a comprehensive mRNA database, wherein metrics such as fragments per kilobase million (FPKM), coefficient of variation (CV), maximum fold change (MFC), and dispersion measure (DPM) are utilized to identify stable or novel HKGs exhibiting low variability [17, 19].

Currently, several HKGs have been identified in testicular tissues and associated cells. In humans, the variability of ACTB, GAPDH and HSPCB (Heat Shock Protein 90 Alpha Family Class B Member 1) was not statistically significant in ejaculates, cryopreserved samples, or formalin-fixed, paraffin-embedded testicular tissues [20]. Meanwhile, the mRNA levels of B2M (beta-2-microglobulin gene) and PPIA in Sertoli cells were unaffected by transcription-disruptive toxins [21]. In animals, B2M and GAPDH minimal the least variation in testes and associated tissues [22]. Moreover, PPIA, GAPDH and ACTB were identified as suitable HKGs in mouse and rat testes during embryonic and postnatal development [23, 24]. B2M and HSP90 (Heat Shock Protein 90) were the most stable in sheep testis tissues in both uninfected and sheep pox virus-infected groups [25]. For yak testicular development, TBP (TATA box (Hogness box)-binding protein) was ranked highest in stability across four developmental stages from 0.5 to 6 years of age [26]. However, appropriate HKGs for testicular tissues under heavy metal toxicity combined with traditional Chinese medicine (TCM) treatment have not yet been investigated, potentially affecting the accuracy of functional gene expression analyses in studies of male reproductive disorders.

In this study, optimal HKG selection was performed in rat testicular tissues and Sertoli cells exposed to manganese and subsequently treated with TCM. Seventeen candidate HKGs were initially identified based on transcriptomic data and existing literature related to testicular tissues and cells. The Ct values of these candidate HKGs were determined by qRT-PCR, and their stability was assessed using four analytical programs (geNorm, NormFinder, BestKeeper, and Delta Ct) and two comprehensive indices (RefFinder and ComprFinder) across three experimental groups: Control (Ctrl), manganese (Mn2+), and manganese combined with CcE (Mn2+ + CcE). The top three most stable HKGs were further validated by normalizing the expression levels of four target genes in both testicular tissues and Sertoli cells. This rigorous selection and validation process ensures reliable normalization of gene expression data, thereby facilitating an accurate investigation of molecular mechanisms underlying testicular spermatogenic disorders and their amelioration through TCM. These findings provide crucial references for selecting appropriate HKGs, ensuring accurate analyses of functional gene expression in studies of testicular diseases and therapeutic interventions.

2. Material and methods

2.1 Animals and sample collection

A total of 27 male Sprague Dawley rats, aged 12 weeks with an average body weight of 230 ± 18 g, were used in this study. They were randomly and equally assigned to three groups: control (Ctrl), manganese exposed (Mn2+), and manganese exposure with Coridius chinensis extract treatment (Mn2+ + CcE). Rats in the Ctrl group were intraperitoneally injected with 0.5 mL/day of 0.9% normal saline and administered 1 mL/day of distilled water by gavage. The Mn2+ and Mn2+ + CcE groups received daily intraperitoneal injections of MnCl2 (Manganese chloride,15 mg/kg body weight, 0.5 mL), with concurrent intragastric administration of distilled water (1 mL/day) or CcE (100 mg/kg, 1 mL/day), respectively [12]. Rats were maintained under a 12-hour light/dark cycle with ad libitum access to water and standard rodent chow. Following a 30-day treatment period, all rats were simultaneously euthanized by carbon dioxide anesthesia. Testicular tissues were promptly excised, flash-frozen in liquid nitrogen, and stored at −80 °C to preserve RNA integrity for subsequent extraction.

2.2 Primary Sertoli cells isolation and culture

Following a previously described method [27], testes were obtained from 7-day-old postpartum Sprague Dawley rats, rinsed three times in phosphate-buffered saline (PBS) supplemented with antibiotics (1% penicillin and 1% streptomycin (P1400, Solarbio, Beijing, China)). For cell isolation, rat testes were minced finely and digested in a solution containing 0.1% collagenase IV (C5138, Sigma-Aldrich, Shanghai, China), 0.1% hyaluronidase (H3506, Sigma-Aldrich, Shanghai, China), and 0.25% trypsin (25300, GIBCO, Shanghai, China) at 37 °C for 20 min. Following plating for 48 h, residual germ cells were removed by treatment with 20 mM Tris buffer (pH 7.4) for 3 min. Purified Sertoli cells were subsequently cultured in Dulbecco’s Modified Eagle Medium (DMEM/F12) supplemented with 10% fetal bovine serum (FBS) and maintained in a humidified incubator containing 5% CO2 at 34 °C. Sertoli cells from testicular tissues of Ctrl, Mn2+, and Mn2+ + CcE groups were isolated separately. After 12 h of cultivation, cells were used for subsequent experiments.

2.3 Statistical analysis of candidate housekeeping genes (HKGs) from RNA-Seq data

RNA-Seq data for testicular tissues and Sertoli cells from Ctrl, Mn2+, and Mn2+ + CcE groups were obtained from our previous study. Transcriptomic sequencing was performed on nine testicular tissue samples and nine Sertoli cell samples using paired-end sequencing technology on an Illumina NovaSeq 6000 platform (Illumina, San Diego, CA, USA). Following assembly and annotation, gene expression levels were quantified as fragments per kilobase of transcript per million mapped reads (FPKM) using the BMKCloud platform (www.biocloud.net). Based on the FPKM values, the coefficient of variation (CV), maximal fold change (MFC), and detection p-value measure (DPM) were calculated for each gene using Microsoft Excel. According to prior research [17], candidate HKGs were selected based on the following criteria: FPKM >100, CV <15%, MFC <1.5, and DPM <0.3.

2.4 RNA extraction and complementary DNA (cDNA) synthesis

Total RNA was isolated from testicular tissues and Sertoli cells using Trizol® reagent according to the manufacturer’s protocols (15596026CN, Thermo Fisher Scientific, Waltham, MA, USA). To eliminate genomic DNA contamination, the extracted RNA was treated with DNase I (2270A, Takara Biotechnology, Dalian, China). RNA concentrations were quantified utilizing the Qubit® RNA Assay Kit (Q10210, Waltham, MA, USA) on a Qubit® 2.0 Fluorometer (Thermo Fisher Scientific, Waltham, MA, USA). RNA purity was assessed using a NanoPhotometer® spectrophotometer (N50, IMPLEN, Munich, BY, Germany), by measuring optical densities at 260 nm (OD260) and 280 nm (OD280). Only RNA samples exhibiting an OD260/280 ratio between 1.8 and 2.0 were considered suitable for subsequent analyses. For cDNA synthesis, 1 μg of RNA from each sample was reverse transcribed using the All-In-One 5X RT MasterMix (Abm, Cat#G592, Shanghai, China). The reverse transcription protocol consisted of incubation at 37 °C for 15 minutes, followed by 60 °C for 10 minutes, and a final denaturation at 95 °C for 3 minutes.

2.5 Quantitative real-time polymerase chain reaction (qRT-PCR)

To ensure the reliability of qRT-PCR analysis, raw gene expression data were normalized to minimize variations among samples and experimental conditions according to the Minimum Information for Publication of Real-Time Quantitative PCR Experiments (MIQE) guidelines [28]. The SYBR green qRT-PCR assay was performed using a Bio-Rad CFX 96 PCR system (Bio-Rad, Hercules, CA, USA) with BlasTaq 2× qPCR MasterMix (Abm, Cat#G891, Shanghai, China). Reaction mixtures were prepared in a total volume of 10 μL, consisting of 5.0 μL of 2× qPCR MasterMix, 4.0 μL of cDNA template, and 0.5 μL of each primer. Thermocycling was initiated with denaturation at 95 °C for 3 minutes, followed by 40 amplification cycles comprising 15 seconds at 95 °C for denaturation and 1 minute at 60 °C for annealing and extension. Subsequently, gel electrophoresis of PCR products was performed following each PCR run to monitor the formation of non-specific products (Supplementary Fig. 1). Primers were designed utilizing the Primer-BLAST web tool (https://www.ncbi.nlm.nih.gov/tools/primer-blast/) and are detailed in Supplementary Table 1.

2.6 Stability analysis of candidate housekeeping genes (HKGs)

The expression stability of selected HKGs was evaluated employing four widely recognized algorithms. geNorm assesses gene expression stability by calculating the gene stability measure (M value), wherein lower M values correspond to higher stability [13]. NormFinder determines HKG stability by accounting for intra- and inter-group variations critical for accurate normalization [14]. BestKeeper analyzes variability in expression of candidate reference genes by computing the standard deviation (SD) from Ct data [15]. Additionally, the Delta Ct method compares relative gene expression between pairs of genes within each sample [16]. After applying these standard algorithms, we conducted an integrated ranking by employing supplementary computational tools, specifically RefFinder [18] and ComprFinder [17], to consolidate and refine the stability evaluations.

2.7 Validation of selected housekeeping genes (HKGs) and their combinations

To corroborate the screening results, we selected candidate HKGs from the extremes of the stability ranking, the three most stable and the three least stable, for normalizing the expression of target genes. Specifically, CASPASE3, DDX4, SOX9, and OCCLUDIN were chosen as target genes in both testicular tissues and Sertoli cells. Total mRNA samples previously analyzed from 47 samples were re-extracted for this validation. When multiple genes were combined, their average Ct values were computed using the geometric mean approach [13]. The relative target gene expression levels were normalized to these designated HKGs and quantified using the 2−ΔΔCt method. Correlation analyses were subsequently performed to explore relationships among these candidate HKGs and their combinations, assessing the relative expression of target genes normalized to different HKG sets.

2.8 Statistical analysis

Cycle threshold (Ct) values and relative mRNA expression data are presented as mean ± standard error of the mean (SEM), whereas indices calculated from algorithms are reported as individual values. All statistical analyses were conducted utilizing SPSS version 23.0 (SPSS Inc., Armonk, NY, USA), and graphical representations were prepared with GraphPad Prism version 6.0 (GraphPad Software, Santiago, CA, USA). Comparisons among multiple groups were performed using one-way analysis of variance (ANOVA) followed by the least significant difference (LSD) post hoc test. Pearson’s correlation coefficient was applied to analyze correlations in relative target gene expression normalized by different HKGs. A p-value of less than 0.05 was regarded as statistically significant.

3. Results

3.1 Selection of candidate HKGs in testicular tissues and Sertoli cells following Mn2+ exposure and CcE treatment using transcriptomic data

Histological examination (HE) staining revealed pronounced detachment of spermatogenic cells from seminiferous tubule basal membranes following Mn2+ exposure (Fig. 1A). Notably, this damage was alleviated after subsequent treatment with CcE (Fig. 1A). Similarly, Sertoli cells cultured in vitro displayed morphological shrinkage and increased cell death after Mn2+ treatment, while CcE administration ameliorated these deleterious effects (Fig. 1B). Profiling of testicular tissue and Sertoli cells exposed to Mn2+ and then treated with CcE identified 761 and 195 overlapping genes, respectively, according to stringent criteria including FPKM, CV, MFC, and DPM (Fig. 1C). A four-color Venn diagram illustrated the intersections of these gene sets, revealing 41 shared genes between testicular tissue and Sertoli cell samples. Gene Ontology (GO) enrichment analysis indicated these candidate genes primarily participated in biological processes involving protein folding, cytoplasmic cellular components, and unfolded protein-binding molecular functions (Supplementary Fig. 2). Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis annotated these genes to pathways such as protein processing in the endoplasmic reticulum and antigen processing and presentation (Supplementary Fig. 3). Subsequently, 12 genes (UBB, ACTB, SOD1 (Superoxide dismutase 1), PPIA, EEF1A1 (Eukaryotic Translation Elongation Factor 1 Alpha 1), TPT1 (Tumor Protein Translationally Controlled 1), CALM2 (Calmodulin 2), CALR (Calreticulin), YWHAE (Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta), CCT5 (Chaperonin Containing TCP1 Subunit 5), DDX5 (DEAD-Box Helicase 5), and YBX1 (Y-box binding protein 1)) were selected for further evaluation as candidate HKGs (Fig. 1D). Additionally, five previously validated HKGs (B2M, TBP, HSP90AA1 (Heat Shock Protein 90 Alpha Family Class A Member 1), GAPDH, and TUBULIN) from testicular tissues and Sertoli cells were included. The transcriptome derived FPKM, CV, MFC, and DPM values for these candidate HKGs are detailed in Table 1. Under conditions of manganese and CcE treatment, UBB exhibited the highest average FPKM values in both testicular tissue and Sertoli cells, with values of 616.08 and 1173.49, respectively. The lowest CV and MFC values were observed in YWHAE (3.11%, 1.13) in Sertoli cells and EEF1A1 (2.42%, 1.06) in testicular tissue, indicating high expression stability. Additionally, CALR and YWHAE demonstrated the lowest DPM values in testicular tissue and Sertoli cells, respectively, with values of 0.03 and 0.04. These findings suggest that UBB, YWHAE, and EEF1A1 are robust candidates for use as HKGs in subsequent quantitative analyses. Among them, UBBdisplayed the highest average FPKM expression levels under conditions involving Mn2+ and CcE, with values of 616.08 in testicular tissue and 1173.49 in Sertoli cells. The lowest CV and MFC values, indicative of stable expression, were observed for YWHAE in Sertoli cells (CV = 3.11%, MFC = 1.13) and EEF1A1 in testicular tissue (CV = 2.42%, MFC = 1.06). Furthermore, CALR (0.03) in testicular tissue and YWHAE (0.04) in Sertoli cells exhibited the lowest DPM values. These collective findings highlight UBB, YWHAE, and EEF1A1 as robust HKG candidates suitable for subsequent quantitative gene expression analyses.

Histological and genetic alterations in testicular tissues and Sertoli cells following heavy metal exposure and rescue treatment with traditional Chinese medicine. (A) Representative histological sections (H&amp;E staining) illustrating the seminiferous tubules and adjacent interstitial tissue of testes from the control (Ctrl), manganese-exposure (Mn2+), and manganese-exposed treated with Coridius chinensis extract (Mn2+ + CcE) groups. Scale bars = 100 μm. (B) Representative microscopic images depicting Sertoli cells in the Ctrl, Mn2+, and Mn2+ + CcE groups. Scale bars = 50 μm. (C) Venn diagram illustrating the number of overlapping genes identified through mRNA sequencing based on stringent filtering criteria, including FPKM (&gt;100, blue area), CV (&lt;15%, yellow area), MFC (&lt;1.5, green area), and DPM (&lt;0.3, red area). Red circles indicate the overlapping gene sets from testicular tissue and Sertoli cells, respectively. (D) Combined Venn diagram analysis of overlapping genes from both testicular tissues and Sertoli cells. Twelve candidate housekeeping genes (HKGs) highlighted within the dotted box were selected from the overlapping gene set. FPKM: fragments per kilobase million; CV: coefficient of variation; MFC: maximum fold change; DPM: dispersion measure; ACTB: Actin Beta; CALM2: Calmodulin 2; CALR: Calreticulin; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; PPIA: Peptidylprolyl Isomerase A; SOD1: Superoxide dismutase 1; TPT1: Tumor Protein Translationally Controlled 1; UBB: Ubiquitin B; YBX1: Y-box binding protein 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta.

Fig. 1.Histological and genetic alterations in testicular tissues and Sertoli cells following heavy metal exposure and rescue treatment with traditional Chinese medicine. (A) Representative histological sections (H&E staining) illustrating the seminiferous tubules and adjacent interstitial tissue of testes from the control (Ctrl), manganese-exposure (Mn2+), and manganese-exposed treated with Coridius chinensis extract (Mn2+ + CcE) groups. Scale bars = 100 μm. (B) Representative microscopic images depicting Sertoli cells in the Ctrl, Mn2+, and Mn2+ + CcE groups. Scale bars = 50 μm. (C) Venn diagram illustrating the number of overlapping genes identified through mRNA sequencing based on stringent filtering criteria, including FPKM (>100, blue area), CV (<15%, yellow area), MFC (<1.5, green area), and DPM (<0.3, red area). Red circles indicate the overlapping gene sets from testicular tissue and Sertoli cells, respectively. (D) Combined Venn diagram analysis of overlapping genes from both testicular tissues and Sertoli cells. Twelve candidate housekeeping genes (HKGs) highlighted within the dotted box were selected from the overlapping gene set. FPKM: fragments per kilobase million; CV: coefficient of variation; MFC: maximum fold change; DPM: dispersion measure; ACTB: Actin Beta; CALM2: Calmodulin 2; CALR: Calreticulin; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; PPIA: Peptidylprolyl Isomerase A; SOD1: Superoxide dismutase 1; TPT1: Tumor Protein Translationally Controlled 1; UBB: Ubiquitin B; YBX1: Y-box binding protein 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta.

Table 1.The information of 17 HKGs in transcriptome data of Sertoli cells and testis tissues.
TypeGene symbolTestis tissuesSertoli cells
Mean FPKM (N)CV (%)MFCDPMMean FPKM (N)CV (%)MFCDPM
Predicted candidate HKGs
UBB616.082.521.070.061719.168.401.350.14
ACTB564.826.261.170.12415.585.641.240.11
SOD1413.789.841.300.07158.479.201.300.05
PPIA399.778.381.240.14504.966.761.270.14
EEF1A1313.602.421.060.06536.928.051.290.12
TPT1290.677.521.210.09288.485.531.250.08
CALM2252.355.911.170.05408.666.971.310.10
CALR246.038.021.230.03220.854.991.230.09
YWHAE238.088.311.260.09187.023.111.130.04
CCT5219.877.401.200.09389.644.141.160.09
DDX5210.395.071.120.11537.555.381.190.12
YBX1191.144.541.130.04663.875.631.200.11
Additional reported HKGs
B2M31.6129.063.530.28517.2962.385.950.13
TBP83.186.071.220.095.775.571.160.15
HSP90AA1402.156.691.240.13452.4312.911.370.15
GAPDH1087.545.811.240.0915.7021.901.850.10
TUBULIN873.803.991.140.15314.9723.041.920.17
FPKM: fragments per kilobase million; CV: coefficient of variation; MFC: maximum fold change; DPM: dispersion measure; HKGs: housekeeping genes; UBB: Ubiquitin B; ACTB: Actin Beta; SOD1: Superoxide dismutase 1; PPIA: Peptidylprolyl Isomerase A; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; TPT1: Tumor Protein Translationally Controlled 1; CALM2: Calmodulin 2; CALR: Calreticulin; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; YBX1: Y-box binding protein 1; B2M: Beta-2-Microglobulin; TBP: TATA-Box Binding Protein; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase; TUBULIN: Tubulin.

3.2 Expression profiling of candidate HKGs via qRT-PCR

For initial selection of appropriate HKGs, Ct values and associated fold changes of candidate genes were measured via qRT-PCR. In testicular tissue, Ct values for candidate HKGs ranged from 21.23 ± 1.02 to 28.80 ± 1.91 (Fig. 2A). Among these, EEF1A1, PPIA, and YBX1 showed the lowest Ct values (21.23 ± 1.02, 21.44 ± 0.84, and 22.45 ± 0.85, respectively). Conversely, B2M, UBB, and TUBULIN exhibited the highest Ct values (28.80 ± 1.91, 27.92 ± 2.43, and 27.75 ± 2.50, respectively) (Fig. 2A). After normalization against the mean Ct values of all 17 candidate genes, fold-change values ranged from 2.14 ± 1.99 to 198.58 ± 54.80 (Fig. 2B). Notably, despite their initially low Ct values, EEF1A1, PPIA, and YBX1 displayed the greatest fold changes (198.58 ± 54.80, 158.82 ± 41.90, and 80.65 ± 20.86, respectively). In contrast, B2M, UBB, and TUBULIN exhibited the lowest fold-change values, ranging from 2.45 ± 1.22 to 4.97 ± 2.41 (Fig. 2B). In Sertoli cells, candidate HKG Ct values varied from 21.81 ± 0.78 to 29.21 ± 1.17 (Fig. 2C). The lowest Ct values were observed for PPIA, HSP90AA1, and GAPDH, ranging from 21.81 ± 0.78 to 24.27 ± 1.74 (Fig. 2C), while CALM2, ACTB, and TBP showed the highest Ct values (29.21 ± 1.17, 27.17 ± 0.78, and 26.91 ± 1.15, respectively). Following normalization to the average Ct value, fold changes in Sertoli cells ranged from 2.14 ± 1.99 to 198.58 ± 54.80 (Fig. 2D), with PPIA exhibiting the highest and CALM2 the lowest fold-change values (Fig. 2D). Integrated analysis of testicular tissues and Sertoli cell data revealed that PPIA maintained the lowest Ct value, while B2M exhibited the highest Ct values (Fig. 2E). After normalization to the average Ct across all candidate genes, PPIA demonstrated the highest fold change, while B2M exhibited the lowest fold-change values (Fig. 2F).

Analysis of Ct values and relative expression levels of 17 candidate HKGs in testicular tissues and Sertoli cells by qRT-PCR. (A,C,E) Ct values for candidate HKGs from testicular tissues, Sertoli cells, and combined samples in the control (Ctrl), Mn2+ exposure, and Mn2+ + CcE treatment groups (N = 4–6 samples per group). The dotted line denotes a Ct threshold value of 28, used as the selection criterion for further analysis based on high gene abundance. (B,D,F) Relative gene expression levels of HKGs normalized to the geometric mean of all candidate HKGs within each sample. Data are presented as means ± SEM. ACTB: Actin Beta; B2M: Beta-2-Microglobulin; CALM2: Calmodulin 2; CALR: Calreticulin; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; PPIA: Peptidylprolyl Isomerase A; SOD1: Superoxide dismutase 1; TBP: TATA-Box Binding Protein; TPT1: Tumor Protein Translationally Controlled 1; TUBULIN: Tubulin; UBB: Ubiquitin B; YBX1: Y-box binding protein 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta.

Fig. 2.Analysis of Ct values and relative expression levels of 17 candidate HKGs in testicular tissues and Sertoli cells by qRT-PCR. (A,C,E) Ct values for candidate HKGs from testicular tissues, Sertoli cells, and combined samples in the control (Ctrl), Mn2+ exposure, and Mn2+ + CcE treatment groups (N = 4–6 samples per group). The dotted line denotes a Ct threshold value of 28, used as the selection criterion for further analysis based on high gene abundance. (B,D,F) Relative gene expression levels of HKGs normalized to the geometric mean of all candidate HKGs within each sample. Data are presented as means ± SEM. ACTB: Actin Beta; B2M: Beta-2-Microglobulin; CALM2: Calmodulin 2; CALR: Calreticulin; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; PPIA: Peptidylprolyl Isomerase A; SOD1: Superoxide dismutase 1; TBP: TATA-Box Binding Protein; TPT1: Tumor Protein Translationally Controlled 1; TUBULIN: Tubulin; UBB: Ubiquitin B; YBX1: Y-box binding protein 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta.

3.3 Stability assessment of HKGs using four analytical algorithms under different experimental conditions

The stability of candidate HKGs was first evaluated using the geNorm algorithm, which employs the M value as a stability indicator, with lower values signifying higher stability (Supplementary Table 2). In testicular tissues, EEF1A1 exhibited the greatest stability (lowest M value), whereas TUBULIN displayed the least stability (Fig. 3A). For Sertoli cells, ACTB was identified as the most stable, and GAPDH the least (Fig. 3B). In combined analysis of testicular tissues and Sertoli cells, PPIA and TBP showed identical M values, jointly ranking as the most stable, while HSP90AA1 exhibited the lowest stability (Fig. 3C). Additionally, geNorm calculated the pairwise variation values (V2/V3) as 0.010 in testicular tissues, and 0.013 both in Sertoli cells and combined samples (Supplementary Fig. 4). These results suggest that employing a single stable HKG would be adequate for accurate normalization of qRT-PCR data.

The NormFinder method (Supplementary Table 3) indicated that PPIA had the highest stability across testicular tissues, Sertoli cells, and combined groups (Fig. 3D–F). Conversely, TUBULIN, GAPDH, and HSP90AA1 were identified as least stable in testicular tissues, Sertoli cells, and combined groups, respectively (Fig. 3D–F).

The BestKeeper algorithm assessed HKG stability by calculating SD values, with lower SD values indicating greater stability. In testicular tissues and the combined dataset, PPIA consistently ranked as the most stable HKG (Fig. 3G,H). In contrast, UBB and HSP90AA1 showed the highest variability in testicular tissues and combined samples, respectively (Fig. 3G,H). In Sertoli cells, UBB demonstrated the lowest SD value (highest stability), followed closely by PPIA, whereas TUBULIN exhibited the highest SD, indicating the lowest stability (Fig. 3I). Comprehensive SD values for each HKG are provided in Supplementary Table 4.

The Delta Ct approach further provided comprehensive stability evaluations for candidate HKGs (Supplementary Table 5), with lower Delta Ct values indicating superior stability. In testicular tissues, YBX1 and PPIA emerged as most stable, whereas TUBULIN and HSP90AA1 were the least stable (Fig. 3J). In Sertoli cells, PPIA was again ranked highest, and GAPDH lowest, consistent with NormFinder results (Fig. 3K). In combined analysis, PPIA retained the highest stability, while HSP90AA1 was the least stable (Fig. 3L).

Stability rankings of candidate HKGs assessed using geNorm, NormFinder, Bestkeeper and Delta Ct algorithms. (A,D,G,J) Stability rankings of candidate HKGs in testicular tissues according to the four algorithms. (B,E,I,K) Stability rankings of candidate HKGs in Sertoli cells. (C,F,H,L) Stability rankings of candidate HKGs in combined testicular tissues and Sertoli cells. Stability decreases sequentially from left (highest) to right (lowest). Green dots indicate stability values in testicular tissues, red dots indicate Sertoli cells, and blue dots represent combined samples. Ct: Cycle threshold; std: Standard Deviation Value; STDEV: Standard Deviation; ACTB: Actin Beta; B2M: Beta-2-Microglobulin; CALM2: Calmodulin 2; CALR: Calreticulin; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; PPIA: Peptidylprolyl Isomerase A; SOD1: Superoxide dismutase 1; TBP: TATA-Box Binding Protein; TPT1: Tumor Protein Translationally Controlled 1; TUBULIN: Tubulin; UBB: Ubiquitin B; YBX1: Y-box binding protein 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta.

Fig. 3.Stability rankings of candidate HKGs assessed using geNorm, NormFinder, Bestkeeper and Delta Ct algorithms. (A,D,G,J) Stability rankings of candidate HKGs in testicular tissues according to the four algorithms. (B,E,I,K) Stability rankings of candidate HKGs in Sertoli cells. (C,F,H,L) Stability rankings of candidate HKGs in combined testicular tissues and Sertoli cells. Stability decreases sequentially from left (highest) to right (lowest). Green dots indicate stability values in testicular tissues, red dots indicate Sertoli cells, and blue dots represent combined samples. Ct: Cycle threshold; std: Standard Deviation Value; STDEV: Standard Deviation; ACTB: Actin Beta; B2M: Beta-2-Microglobulin; CALM2: Calmodulin 2; CALR: Calreticulin; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; PPIA: Peptidylprolyl Isomerase A; SOD1: Superoxide dismutase 1; TBP: TATA-Box Binding Protein; TPT1: Tumor Protein Translationally Controlled 1; TUBULIN: Tubulin; UBB: Ubiquitin B; YBX1: Y-box binding protein 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta.

3.4 Comprehensive ranking of candidate HKGs Using RefFinder and ComprFinder

Although rankings of gene stability within specific tissues or cell types were generally consistent across the four previously applied algorithms, each algorithm identified different optimal HKGs. To integrate the individual rankings and produce a unified stability assessment, the RefFinder and ComprFinder algorithms were applied to the results obtained from the earlier methods (Tables 2 and 3). Consistent with previous algorithms, lower comprehensive scores indicate higher stability. Notably, both RefFinder and ComprFinder consistently identified the same three most stable candidate HKGs, regardless of tissue or cell type. In testicular tissues, both comprehensive analyses indicated YBX1, PPIA, and EEF1A1 as the most stable HKGs, while HSP90AA1, UBB, and TUBULIN were consistently the least stable (Tables 2 and 3). In Sertoli cells, the top three stable HKGs were identified as PPIA, ACTB, and UBB, with PPIA consistently ranking first. Conversely, B2M, TUBULIN, and GAPDH were identified as least stable. When considering combined testicular tissues and Sertoli cells, the most stable HKGs were PPIA, TBP, and ACTB, while UBB, TUBULIN, and HSP90AA1 were ranked as least stable, receiving the highest scores in both RefFinder (Table 2) and ComprFinder (Table 3).

Table 2.Comprehensive rankings calculated using the RefFinder method.
Ranking numberTestis tissuesSertoli cellsTotal
GeneScoreGeneScoreGeneScore
1YBX11.57PPIA1.57PPIA1.00
2PPIA1.57ACTB1.86ACTB1.68
3EEF1A12.78UBB2.99TBP3.22
4DDX54.12DDX53.03TPT14.16
5TBP4.90CCT54.61SOD15.42
6ACTB4.95HSP90AA15.18CCT56.82
7CALM26.65TBP6.74DDX56.96
8GAPDH8.65YBX18.00CALR7.09
9SOD19.07CALM29.00GAPDH8.97
10TPT19.24TPT110.49EEF1A19.80
11CALR10.44SOD110.49YBX19.92
12CCT511.74EEF1A112.72YWHAE11.70
13B2M13.00CALR13.21B2M13.49
14YWHAE14.00YWHAE13.47CALM213.88
15UBB15.48B2M14.47UBB13.98
16TUBULIN16.23TUBULIN16.00TUBULIN16.00
17HSP90AA117.00GAPDH17.00HSP90AA117.00
ACTB: Actin Beta; B2M: Beta-2-Microglobulin; CALM2: Calmodulin 2; CALR: Calreticulin; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; PPIA: Peptidylprolyl Isomerase A; SOD1: Superoxide dismutase 1; TBP: TATA-Box Binding Protein; TPT1: Tumor Protein Translationally Controlled 1; TUBULIN: Tubulin; UBB: Ubiquitin B; YBX1: Y-box binding protein 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta.
Table 3.Comprehensive rankings calculated using the ComprFinder method.
Ranking numberTestis tissuesSertoli cellsTotal
GeneScoreGeneScoreGeneScore
1YBX10.0197PPIA0.0456PPIA0.0021
2PPIA0.0385ACTB0.0466ACTB0.0231
3EEF1A10.0718UBB0.0761TBP0.1155
4DDX50.1392DDX50.0800TPT10.2289
5ACTB0.1753CCT50.1536SOD10.3158
6TBP0.2062HSP90AA10.1700CALR0.3368
7CALM20.2582EEF1A10.2713DDX50.3447
8GAPDH0.3476YBX10.2980CCT50.4011
9SOD10.3842TBP0.3131GAPDH0.4404
10TPT10.4189CALM20.3648EEF1A10.5327
11CALR0.4574TPT10.4138YBX10.5367
12CCT50.5030SOD10.4226YWHAE0.6663
13B2M0.5781CALR0.5331B2M0.7358
14YWHAE0.6649YWHAE0.5381CALM20.7420
15UBB0.9123B2M0.5741UBB0.7578
16HSP90AA10.9623GAPDH0.7537TUBULIN0.8692
17TUBULIN1.0000TUBULIN1.0000HSP90AA11.0000
ACTB: Actin Beta; B2M: Beta-2-Microglobulin; CALM2: Calmodulin 2; CALR: Calreticulin; CCT5: Chaperonin Containing TCP1 Subunit 5; DDX5: DEAD-Box Helicase 5; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; PPIA: Peptidylprolyl Isomerase A; SOD1: Superoxide dismutase 1; TBP: TATA-Box Binding Protein; TPT1: Tumor Protein Translationally Controlled 1; TUBULIN: Tubulin; UBB: Ubiquitin B; YBX1: Y-box binding protein 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta.

3.5 Validation of optimal HKGs through target gene expression

To validate the normalization accuracy provided by the optimal HKGs, four target genes (CASPASE3, DDX4, SOX9, and OCCLUDIN), specifically expressed in seminiferous tubules or Sertoli cells, were selected for assessment. Relative expression levels of these target genes were normalized to the most stable and least stable HKGs. RNA sequencing data provided the average FPKM values for these genes under manganese exposure and subsequent CcE treatment in testicular tissues and Sertoli cells (Fig. 4A–D). Average FPKM values of DDX4, SOX9, and OCCLUDIN significantly decreased in the Mn2+ group relative to the control, but were restored to control levels after CcE treatment (Fig. 4B–D). Particularly, CASPASE3 expression in the control and Mn2+ + CcE groups was significantly lower and nearly identical, compared with elevated expression in the Mn2+ group (Fig. 4A). qRT-PCR analysis revealed similar expression patterns for these four target genes when normalized against the most stable and unstable HKGs, consistent with the RNA sequencing data (Fig. 4E–H). In testicular tissues, relative expressions of CASPASE3 and DDX4, normalized using the three most stable HKGs (YBX1, PPIA, and EEF1A1), were consistent across treatment groups (Fig. 4E,F). However, normalization with the three least stable HKGs (YWHAE, UBB, HSP90AA1) resulted in significant distortion of expression patterns for these genes (Fig. 4E,F). In Sertoli cells, normalization with the two most stable HKGs (PPIA and ACTB) produced comparable expression levels of SOX9 and OCCLUDIN across groups (Fig. 4G,H). Unfortunately, despite previous stability indications, UBB provided unreliable normalization for SOX9 and OCCLUDIN (Fig. 4G,H). Furthermore, normalization using the three least stable HKGs (YWHAE, B2M, GAPDH) led to inconsistent and chaotic expression profiles across treatments (Fig. 4G,H).

Relative FPKM values of target genes and their normalized relative expression using the three most and least stable HKGs. (A,B) Relative FPKM values of CASPASE3 and DDX4 derived from transcriptomic data in testicular tissues across treatment groups. (C,D) Relative FPKM values of SOX9 and OCCLUDIN from transcriptomic analysis of Sertoli cells in three treatment groups. (E,F) Relative expression levels of CASPASE3 and DDX4 in testicular tissues, normalized using either the three most or three least stable HKGs under different treatments. (G,H) Relative expression levels of SOX9 and OCCLUDIN in Sertoli cells, normalized to the three most and least stable HKGs across treatments. ns, no significant. *p &lt; 0.05. **p &lt; 0.01. ns: no significance; FPKM: fragments per kilobase million; Ctrl: control; Mn2+: manganese ion; CcE: Coridius chinensis extract; CASPASE3: Cysteine-dependent aspartate-specific protease-3; DDX4: DEAD-Box Helicase 4; SOX9: SRY-box transcription factor 9; OCCLUDIN: Occludin; YBX1: Y-box binding protein 1; PPIA: Peptidylprolyl Isomerase A; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta; UBB: Ubiquitin B; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; ACTB: Actin Beta; B2M: Beta-2-Microglobulin; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase.

Fig. 4.Relative FPKM values of target genes and their normalized relative expression using the three most and least stable HKGs. (A,B) Relative FPKM values of CASPASE3 and DDX4 derived from transcriptomic data in testicular tissues across treatment groups. (C,D) Relative FPKM values of SOX9 and OCCLUDIN from transcriptomic analysis of Sertoli cells in three treatment groups. (E,F) Relative expression levels of CASPASE3 and DDX4 in testicular tissues, normalized using either the three most or three least stable HKGs under different treatments. (G,H) Relative expression levels of SOX9 and OCCLUDIN in Sertoli cells, normalized to the three most and least stable HKGs across treatments. ns, no significant. *p < 0.05. **p < 0.01. ns: no significance; FPKM: fragments per kilobase million; Ctrl: control; Mn2+: manganese ion; CcE: Coridius chinensis extract; CASPASE3: Cysteine-dependent aspartate-specific protease-3; DDX4: DEAD-Box Helicase 4; SOX9: SRY-box transcription factor 9; OCCLUDIN: Occludin; YBX1: Y-box binding protein 1; PPIA: Peptidylprolyl Isomerase A; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta; UBB: Ubiquitin B; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; ACTB: Actin Beta; B2M: Beta-2-Microglobulin; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase.

3.6 Correlation analysis of different HKGs by normalizing target gene expression

To investigate relationships between stable and unstable HKGs, correlation analyses were conducted for the relative expression of four target genes normalized to different sets of HKGs. In testicular tissues, relative CASPASE3 expression normalized by stable HKGs (YBX1, PPIA, EEF1A1) showed strong correlations (R = 0.81–0.88, p < 0.001; Fig. 5A). Conversely, expression normalized by unstable HKGs (YWHAE, UBB, HSP90AA1) exhibited no significant inter-correlations (p > 0.05), though YWHAE and HSP90AA1 did correlate significantly with stable HKGs (p < 0.01). Similarly, DDX4 expression normalized by stable HKGs (YBX1, PPIA, EEF1A1) was significantly correlated (R = 0.62–0.77, p < 0.001; Fig. 5B). Notably, unstable HKGs (UBB, HSP90AA1) also showed significant correlations with stable HKGs (R = 0.53–0.73, p < 0.01; Fig. 5B). In Sertoli cells, relative SOX9 expression normalized by stable HKGs (PPIA, ACTB, UBB) displayed significant correlations (R = 0.48–0.88, p < 0.01; Fig. 5C). Interestingly, normalization by the three least stable HKGs (YWHAE, B2M, GAPDH) also yielded significant correlations (p < 0.001). However, no significant correlations were found between stable genes (e.g., PPIA, ACTB) and unstable HKG sets (p > 0.05; Fig. 5C). For the OCCLUDIN, significant correlations were evident within stable and unstable HKG normalization sets (R = 0.49–0.92, p < 0.05; Fig. 5D). Specifically, the two most stable HKGs (PPIA and ACTB) demonstrated negligible correlations with the least stable HKGs, highlighting significant differences in their normalization performance (p > 0.05; Fig. 5D).

Correlation analysis of target gene expression normalized to different housekeeping genes (HKGs) in testicular tissues and Sertoli cells. (A,B) Correlation heatmaps illustrating relative expression of CASPASE3 and DDX4 normalized to YBX1, PPIA, EEF1A1, YWHAE, UBB, and HSP90AA1 in testicular tissues. (C,D) Correlation heatmaps depicting relative expression of SOX9 and OCCLUDIN normalized to PPIA, ACTB, UBB, YWHAE, B2M, and GAPDH in Sertoli cells. Statistical significance: *p &lt; 0.05, **p &lt; 0.01, ***p &lt; 0.001. CASPASE3: Cysteine-dependent aspartate-specific protease-3; DDX4: DEAD-Box Helicase 4; SOX9: SRY-box transcription factor 9; OCCLUDIN: Occludin; YBX1: Y-box binding protein 1; PPIA: Peptidylprolyl Isomerase A; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta; UBB: Ubiquitin B; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; ACTB: Actin Beta; B2M: Beta-2-Microglobulin; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase.

Fig. 5.Correlation analysis of target gene expression normalized to different housekeeping genes (HKGs) in testicular tissues and Sertoli cells. (A,B) Correlation heatmaps illustrating relative expression of CASPASE3 and DDX4 normalized to YBX1, PPIA, EEF1A1, YWHAE, UBB, and HSP90AA1 in testicular tissues. (C,D) Correlation heatmaps depicting relative expression of SOX9 and OCCLUDIN normalized to PPIA, ACTB, UBB, YWHAE, B2M, and GAPDH in Sertoli cells. Statistical significance: *p < 0.05, **p < 0.01, ***p < 0.001. CASPASE3: Cysteine-dependent aspartate-specific protease-3; DDX4: DEAD-Box Helicase 4; SOX9: SRY-box transcription factor 9; OCCLUDIN: Occludin; YBX1: Y-box binding protein 1; PPIA: Peptidylprolyl Isomerase A; EEF1A1: Eukaryotic Translation Elongation Factor 1 Alpha 1; YWHAE: Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta; UBB: Ubiquitin B; HSP90AA1: Heat Shock Protein 90 Alpha Family Class A Member 1; ACTB: Actin Beta; B2M: Beta-2-Microglobulin; GAPDH: Glyceraldehyde-3-Phosphate Dehydrogenase.

4. Discussion

The testis is an essential organ in the male reproductive system, critical for sperm production and maintaining secondary sexual characteristics in both humans and animals [1, 2, 3]. Exposure to heavy metals or treatments involving TCM can modulate testicular histology and induce cellular physiological alterations, often accompanied by substantial changes in the expression of various genes, including functional genes and HKGs. Therefore, identifying suitable HKGs is crucial for accurately assessing target gene expression and elucidating the molecular mechanisms underpinning testicular damage and recovery. In this study, we demonstrated that CcE effectively repaired damage to seminiferous tubule epithelia and restored Sertoli cell viability and tight junction integrity disrupted by Mn2+ exposure. Identifying stable HKGs is necessary to reliably quantify functional gene expression and investigate molecular pathways underlying testicular injury and repair. Hence, we utilized transcriptomic data and qRT-PCR results from rat testicular tissues and Sertoli cells exposed to Mn2+ and treated with CcE, selecting candidate HKGs under these conditions. Four separate algorithms and two integrated comprehensive indices were employed to identify the most stable HKGs. We subsequently normalized target gene expression data to these optimal HKGs and conducted correlation analyses to validate HKG stability.

Recent research suggests that transcriptomic sequencing not only facilitates the identification of differentially expressed genes, but also serves as an extensive resource for discovering novel or highly stable HKGs. Such HKG identification methods have been employed across various tissues, such as skin [17], placenta [29], and adipose tissue [30]. In this study, transcriptomic analysis of rat testicular tissues and Sertoli cells under control conditions, Mn2+ exposure, and subsequent CcE treatment provided FPKM values for each gene. To effectively narrow down candidate HKGs from the vast transcriptome dataset, we applied stringent selection criteria based on multiple metrics, specifically the CV, MFC, DPM, and FPKM values. These parameters align with approaches from previous studies [17, 19]. Rigorous thresholds on these four parameters substantially reduced candidate HKG numbers, thereby enhancing reliability. Comparison between testicular tissues and Sertoli cells yielded 41 shared genes, enriched significantly in the protein processing pathway within the endoplasmic reticulum, according to KEGG analysis. Interestingly, this pathway was similarly identified in transcriptomic screening for HKGs in placental development studies [29], suggesting minimal influence from heavy metals, drugs, or developmental processes. Consequently, these identified genes might serve as novel candidate HKGs for testicular studies.

The accuracy of Ct values is essential for reliable quantification of mRNA expression using qRT-PCR. Higher Ct values indicate the requirement for more amplification cycles before detection, reflecting lower initial transcript abundance. For candidate HKGs, lower transcript abundance (high Ct values) frequently corresponds to greater variability due to stochastic fluctuations, especially under varied experimental conditions. Thus, elevated Ct values might compromise measurement reliability and accuracy [31]. Therefore, selected HKGs should ideally exhibit relatively low Ct values to ensure robust and reproducible quantification. Under standardized experimental conditions, including consistent cDNA synthesis methods and dilution protocols—we selected HKGs with Ct values below 28, in line with established recommendations [31, 32], to ensure high reliability. Surprisingly, average Ct values exceeded 28 for B2M, TUBULIN, and UBB in testicular tissues and for TUBULIN and CALM2 in Sertoli cells. Surprisingly, B2M, TUBULIN, and UBB exhibited mean Ct values exceeding 28 in testicular tissue, while TUBULIN and CALM2 demonstrated mean Ct values above 28 in Sertoli cells. Although these genes were excluded as potential HKGs for testicular tissues and related cell types, they were retained in subsequent screening analyses. Notably, PPIA showed the lowest Ct values consistently in both testicular tissues and Sertoli cells, indicating its superior precision and stability in subsequent quantification.

During Mn2+ exposure and CcE treatment, PPIA emerged as the optimal HKG candidate in both testicular tissues and Sertoli cells. PPIA encodes peptidylprolyl isomerase A, and previous studies have identified it among the most stable HKGs in brain tissue [33], polycystic kidneys [34], and osteosarcoma cells [35, 36], demonstrating stable expression across normal and pathological conditions. Additionally, PPIA maintained stability across the hypothalamic-pituitary-gonadal (HPG) axis in pigs [37] and birds [38], and was considered suitable in mouse testes spanning developmental stages, from embryonic day 18 to 4 months postnatally [23]. Our results corroborate the consistent stability of PPIA expression during testicular injury and therapeutic intervention.

We specifically observed stable YBX1 expression in testicular tissues following Mn2+ exposure and CcE treatment. YBX1 encodes a nucleic acid-binding protein implicated in numerous cellular processes, including transcriptional and translational translation [39]. YBX1 was previously identified as a stable internal reference gene in lung tissues infected with ovine pulmonary adenocarcinoma virus [40], with FPKM values >200 and SEM <20, demonstrating stability according to NormFinder analysis [40]. Moreover, YBX1 was abundantly expressed across diverse testicular tissues and cells, including seminiferous tubules, Leydig cells, and spermatogonial stem cells [41]. Our findings further validate the stability of YBX1expression unaffected by heavy metal exposure and TCM treatment, suggesting YBX1 as an ideal HKG candidate for testicular research.

In Sertoli cells, ACTB ranked second in stability after PPIA. ACTB, a widely adopted HKG, has previously shown consistent expression across developmental stages in testicular tissues [23, 42]. Moreover, ACTB expression remained stable in immature porcine Sertoli cells exposed to acute heat stress [43], suggesting that nonspecific physiological stress does not significantly impact its expression. Our findings confirmed ACTB’s robust stability in Sertoli cells subjected to heavy metal exposure and pharmacological intervention.

Conversely, other commonly used HKGs, including GAPDH and TUBULIN, did not demonstrate high stability. GAPDH exhibited only moderate stability in testicular tissues, whereas TUBULIN consistently ranked lowest among the 17 HKG candidates evaluated. Although GAPDH has shown stability during mouse testicular development [23], its expression is vulnerable to perturbation by viral infections and environmental stressors. GAPDH instability was previously observed in lamb testicular cells infected with sheep pox virus [25] and in porcine immature Sertoli cells subjected heat stress [43]. Thus, fluctuations [32] in GAPDH expression under heavy metal exposure or TCM intervention were expected. Similarly, TUBULIN expression instability has been documented in testicular tissues from rodents, poultry [42], and fish [44]. Environmental conditions, such as temperature or photoperiod disturbances, can alter TUBULIN expression [42, 44], indicating that its expression is also unstable in testicular cells after heavy metal poisoning or drug treatment.

The requisite number of HKGs for precise normalization was determined by calculating the average pairwise variation (Vn/Vn+1) using the geNorm algorithm [13]. This method employs a cut-off threshold of Vn/Vn+1 below 0.15 to identify the appropriate number (n) of HKGs needed for reliable normalization [17, 32]. In our study, all Vn/Vn+1 values for both testicular tissues and Sertoli cells were below the 0.15 threshold (Supplementary Fig. 4), indicating that inclusion of additional HKGs is unnecessary for accurate normalization. Our previous research also indicated negligible differences between using single versus multiple stable HKGs for normalization [29, 32]. Moreover, utilizing a single HKG is cost-effective and simplifies experimental workflows. Therefore, employing a single highly stable HKG is sufficient for accurate quantification of target gene expression in qRT-PCR experiments. This conclusion is beneficial for researchers, facilitating convenient and economical selection of a stable reference gene for evaluating gene expression in testicular injury and therapeutic intervention studies.

The target genes have been reported previously as effective markers to validate the stability of candidate HKGs [29, 32, 45]. In this study, four commonly employed target genes (CASPASE3, DDX4, SOX9 and OCCLUDIN) were selected, and their expression levels normalized by candidate HKGs were evaluated in testicular tissues and Sertoli cells subjected to various treatments. Specifically, the relative expression trends of these target genes corresponded closely with their FPKM values obtained from transcriptomic analyses (Fig. 4A–D). CASPASE3 encodes a cysteine-aspartic acid protease involved prominently in the apoptotic signaling pathway [12, 46], with prior reports demonstrating notably lower expression in manganese-exposed testicular tissues treated with CcE [12]. Another gene, DDX4 (also termed DEAD-box Helicase 4), is expressed predominantly in male germ cells and is essential for regulating germ cell proliferation and differentiation across mammalian groups [47, 48, 49]. Our study identified significantly reduced DDX4 expression in damaged seminiferous epithelium following manganese exposure, which was subsequently restored by CcE treatment. In Sertoli cells, SOX9 and OCCLUDIN served as representative biomarkers for assessing the stability of candidate HKGs. SOX9 plays an essential role in mammalian sex determination during fetal development and is persistently expressed in Sertoli cells to regulate their proliferation and differentiation throughout testicular maturation [50, 51]. Our results demonstrated that SOX9 expression was markedly decreased in Mn2+-exposed Sertoli cells but could be restored to control levels by CcE. Additionally, the tight junction is a key structure in Sertoli cells involved in the blood-testis barrier, and its function relies on the expression of the OCCLUDIN gene [11, 52]. Here, OCCLUDIN expression in Sertoli cells was disrupted by Mn2+ exposure but restored by CcE [11]. Thus, these four target genes with varying expression levels can serve as reliable standards for evaluating HKG stability in testicular tissues and cells under disease and treatment conditions.

5. Conclusions

In this study, RNA sequencing and qRT-PCR data were utilized to screen optimal HKGs for testis tissue and Sertoli cells under conditions of manganese exposure and Coridius chinensis extract treatment in rats. After stability analysis using four algorithms and two comprehensive indices, YBX1 and PPIA were identified as the two most stable HKGs in testis tissues under different treatments. Meanwhile, PPIA and ACTB were identified as the top two HKGs in Sertoli cells across all groups. PPIA was comprehensively evaluated as the single optimal HKG for both testis tissues and cells under conditions of heavy metal poisoning and Chinese herbal treatment (Fig. 6). Our results provide robust and valuable recommendations for selecting HKGs in testicular diseases and related cellular physiological processes.

Experimental procedure and conclusions of the study. (A) Sample information for transcriptome data from testis tissues and Sertoli cells in three treatment groups. (B) Candidate HKGs were preliminarily screened using four criteria, including FPKM, CV, DPM, and MFC, and further overlapped by Venn diagram analysis. (C) qRT-PCR was performed with candidate HKG primers in testicular tissues and Sertoli cells. (D) Stable HKGs were determined using four algorithms (geNorm, NormFinder, BestKeeper, and Delta Ct) and two comprehensive analytical methods (ComprFinder and RefFinder). (E) Stable HKGs were validated using four target genes (DDX4, OCCLUDIN, SOX9, and CASPASE3) in testicular tissues and cells. (F) Optimal HKGs were ultimately selected by considering both testis tissues and cells. Ctrl: control; Mn2+: manganese ion; CcE: Coridius chinensis extract; HKGs: Housekeeping Genes; ACTB: Actin Beta; PPIA: Peptidylprolyl Isomerase A; YBX1: Y-box binding protein 1; DDX4: DEAD-Box Helicase 4; OCCLUDIN: Occludin; SOX9: SRY-box transcription factor 9; CASPASE3: Cysteine-dependent aspartate-specific protease-3; CV: coefficient of variation; FPKM: fragments per kilobase million; MFC: maximum fold change; DPM: dispersion measure; qRT-PCR: Quantitative Reverse Transcription Polymerase Chain Reaction; Ct: Cycle threshold.

Fig. 6.Experimental procedure and conclusions of the study. (A) Sample information for transcriptome data from testis tissues and Sertoli cells in three treatment groups. (B) Candidate HKGs were preliminarily screened using four criteria, including FPKM, CV, DPM, and MFC, and further overlapped by Venn diagram analysis. (C) qRT-PCR was performed with candidate HKG primers in testicular tissues and Sertoli cells. (D) Stable HKGs were determined using four algorithms (geNorm, NormFinder, BestKeeper, and Delta Ct) and two comprehensive analytical methods (ComprFinder and RefFinder). (E) Stable HKGs were validated using four target genes (DDX4, OCCLUDIN, SOX9, and CASPASE3) in testicular tissues and cells. (F) Optimal HKGs were ultimately selected by considering both testis tissues and cells. Ctrl: control; Mn2+: manganese ion; CcE: Coridius chinensis extract; HKGs: Housekeeping Genes; ACTB: Actin Beta; PPIA: Peptidylprolyl Isomerase A; YBX1: Y-box binding protein 1; DDX4: DEAD-Box Helicase 4; OCCLUDIN: Occludin; SOX9: SRY-box transcription factor 9; CASPASE3: Cysteine-dependent aspartate-specific protease-3; CV: coefficient of variation; FPKM: fragments per kilobase million; MFC: maximum fold change; DPM: dispersion measure; qRT-PCR: Quantitative Reverse Transcription Polymerase Chain Reaction; Ct: Cycle threshold.

This study has several limitations that should be acknowledged. Our sample was restricted to testicular tissue and Sertoli cells, which may reduce the generalizability of results to the process of spermatogenesis in vitro. Future studies should include all types of spermatogenic cells and Leydig cells. Moreover, the group design for manganism and CcE rescue may limit the applicable scope of the optimal HKG. Adding other clinical-related testicular injury models are needed to confirm the stability of HKG. Despite these limitations, our findings provide a valuable reference for the selection of HKGs in testicular tissues and cells.

Abbreviations

CcE, Coridius chinensis extracts; Cd, Cadmium; Ct, Cycle threshold; Ctrl, Control; CV, Coefficient of variation; DPM, Dispersion measure; HKG, Housekeeping gene; MFC, Maximum fold change; Mn2+, Manganism; Pb, Lead; PFKM, Fragments per kilobase of transcript sequence per million; qRT-PCR, Quantitative real-time polymerase chain reaction; SSCs, Spermatogonia stem cells; ACTB, Actin Beta; B2M, Beta-2-Microglobulin; CALM2, Calmodulin 2; CALR, Calreticulin; CASPASE3, Caspase3; CCT5, Chaperonin Containing TCP1 Subunit 5; DDX4, DEAD-Box Helicase 4; DDX5, DEAD-Box Helicase 5; EEF1A1, Eukaryotic Translation Elongation Factor 1 Alpha 1; GAPDH, Glyceraldehyde-3-Phosphate Dehydrogenase; HSP90AA1, Heat Shock Protein 90 Alpha Family Class A Member 1; HSPCB, Heat Shock Protein 90 Alpha Family Class B Member 1; OCCLUDIN, Occludin; PPIA, Peptidylprolyl Isomerase A; TBP, TATA-Box Binding Protein; TPT1, Tumor Protein Translationally Controlled 1; TUBULIN, Tubulin; SOD1, Superoxide dismutase 1; SOX9, SRY-box transcription factor 9; UBB, Ubiquitin B; YWHAZ, Tyrosine 3-Monooxygenase/Tryptophan 5-Monooxygenase Activation Protein Zeta; YBX1, Y-box binding protein 1.

Availability of data and materials

The datasets used and analyzed during the current study are Included in the article and Supplementary material. Raw and processed RNA-Seq data have been deposited in GEO under accession number GSE339089. Additional data supporting the findings of this study are available from the corresponding author upon reasonable request.

Author contributions

NL, FG and XH—conceived the project. NL, XH and JW—designed the experiments and wrote the manuscript. JW, PT, ZZ and DL—performed the experiments, analyzed the data, and images the graphs. NL, XS and YH—reviewed and edited the manuscript. JW, PT, HL, CZ and XY—revised the manuscript. All authors contributed to the conception of the manuscript. All authors read and approved the final manuscript.

Ethics approval and consent to participate

All experimental procedures in this study were approved by the Ethics Committee of Zunyi Medical University (permit number SYXK2013068). The Sprague Dawley rats were purchased from Chongqing Tengxin Biotechnology Co. Ltd. (license number: Chongqing ICP No. 16006497).

Acknowledgment

We appreciate the Department of Histology and Embryology, School of Preclinical Medicine, Zunyi Medical University for providing experimental samples. And we also thank Guofa Zhang (School of Medical Information Engineering, Zunyi Medical University) for providing statistical analysis.

Funding

This work was financially supported by the Joint Projects of City and School of Science and Technology Bureau of Zunyi City (grant No. Zunshikehe HZ [2019] 6 and Zunshikehe HZ [2023] 180), the Doctoral Research Initiation Project of Zunyi Medical University (grant No. [2020]-36), National Undergraduate Training Program for Innovation and Entrepreneurship (No. ZYDC2021045, ZYDC202202236, and ZYDC202302018), 2021 Academic New Seedling Cultivation and Innovation Exploration Special Project [Grant No. Qian Ke He Basic [2021]1350-016], the Project of Basic Research Program of Guizhou Province [Grant No. Qian Ke He Basic-ZK[2023]-547].

Conflict of interest

The authors declare that they have no competing interests.

Supplementary material

Supplementary material associated with this article can be found, in the online version, at https://files.intandro.com/files/article/2105189194416439296/attachment/Supplementary%20material.docx.

References

Li L, Lin W, Wang Z, Huang R, Xia H, Li Z, et al. Hormone regulation in testicular development and function. International Journal of Molecular Sciences. 2024; 25: 5805.

[Google Scholar]

Barth A, Perry VEA, Hamilton LE, Sutovsky P, Oko R. Bovine spermatogenesis. Advances in Anatomy, Embryology and Cell Biology. 2025; 240: 65–136.

[Google Scholar]

Ishiguro KI. Mechanisms of meiosis initiation and meiotic prophase progression during spermatogenesis. Molecular Aspects of Medicine. 2024; 97: 101282.

[Google Scholar]

O’Donnell L, Smith LB, Rebourcet D. Sertoli cells as key drivers of testis function. Seminars in Cell & Developmental Biology. 2022; 121: 2–9.

[Google Scholar]

Xiao Y, Zhang J, Guan Y, Wang M, Liu D, Xiong S, et al. Research progress on Sertoli cell secretion during spermatogenesis. Frontiers in Endocrinology. 2024; 15: 1456410.

[Google Scholar]

Wanjari UR, Gopalakrishnan AV. Blood-testis barrier: a review on regulators in maintaining cell junction integrity between Sertoli cells. Cell and Tissue Research. 2024; 396: 157–175.

[Google Scholar]

Liu C, Peng H, Yu J, Luo P, Xiong C, Chen H, et al. Impaired ketogenesis in Leydig cells drives testicular aging. Nature Communications. 2025; 16: 4224.

[Google Scholar]

Mitra S, Patra T, Saha D, Ghosh P, Mustafi SM, Varghese AC, et al. Sub-chronic cadmium and lead compound exposure induces reproductive toxicity and development of testicular germ cell neoplasia in situ in murine model: attenuative effects of resveratrol. Journal of Biochemical and Molecular Toxicology. 2022; 36: e23058.

[Google Scholar]

Souza TL, Batschauer AR, Brito PM, Leao-Buchir J, Spercoski KM, Neto FF, et al. Evaluation of Mn exposure in the male reproductive system and its relationship with reproductive dysfunction in mice. Toxicology. 2020; 441: 152504.

[Google Scholar]

Ogunro OB, Agboola AO, Asejeje FO, Gyebi GA. Ethnopharmacological insights: hyperoside from Marsdenia latifolia (Benth.) K.Schum. mitigates reproductive toxicity in male Wistar rats induced by manganese exposure. Natural Product Research. 2026; 40: 1265–1276.

[Google Scholar]

Cen C, Wang F, Xiong K, Jiang L, Hou X. Protective effects of Coridius chinensis extracts on rat reproductive damage induced by manganese. Andrologia. 2022; 54: e14326.

[Google Scholar]

Liu Q, Cen C, Fu H, Wang F, Wang Y, Xu T, et al. Antioxidant activity of Coridius chinensis extracts on manganese-induced testicular damage in rats. Environmental Toxicology. 2019; 34: 1067–1073.

[Google Scholar]

Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, et al. Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biology. 2002; 3: RESEARCH0034.

[Google Scholar]

Andersen CL, Jensen JL, Ørntoft TF. Normalization of real-time quantitative reverse transcription-PCR data: a model-based variance estimation approach to identify genes suited for normalization, applied to bladder and colon cancer data sets. Cancer Research. 2004; 64: 5245–5250.

[Google Scholar]

Pfaffl MW, Tichopad A, Prgomet C, Neuvians TP. Determination of stable housekeeping genes, differentially regulated target genes and sample integrity: BestKeeper—excel-based tool using pair-wise correlations. Biotechnology Letters. 2004; 26: 509–515.

[Google Scholar]

Silver N, Best S, Jiang J, Thein SL. Selection of housekeeping genes for gene expression studies in human reticulocytes using real-time PCR. BMC Molecular Biology. 2006; 7: 33.

[Google Scholar]

Zhang J, Deng C, Li J, Zhao Y. Transcriptome-based selection and validation of optimal house-keeping genes for skin research in goats (Capra hircus). BMC Genomics. 2020; 21: 493.

[Google Scholar]

Xie F, Xiao P, Chen D, Xu L, Zhang B. miRDeepFinder: a miRNA analysis tool for deep sequencing of plant small RNAs. Plant Molecular Biology. 2012; 80: 75–84.

[Google Scholar]

Wang L, Chen X, Song T, Zhang X, Zhan S, Cao J, et al. Using RNA-seq to identify reference genes of the transition from brown to white adipose tissue in goats. Animals. 2020; 10: 1626.

[Google Scholar]

Cavalcanti MC, Failling K, Schuppe HC, Bergmann M, Stalf T, Weidner W, et al. Validation of reference genes in human testis and ejaculate. Andrologia. 2011; 43: 361–367.

[Google Scholar]

Ribeiro MA, dos Reis MB, de Moraes LN, Briton-Jones C, Rainho CA, Scarano WR. Defining suitable reference genes for RT-qPCR analysis on human Sertoli cells after 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) exposure. Molecular Biology Reports. 2014; 41: 7063–7066.

[Google Scholar]

Klein C, Rutllant J, Troedsson MH. Expression stability of putative reference genes in equine endometrial, testicular, and conceptus tissues. BMC Research Notes. 2011; 4: 120.

[Google Scholar]

Gong ZK, Wang SJ, Huang YQ, Zhao RQ, Zhu QF, Lin WZ. Identification and validation of suitable reference genes for RT-qPCR analysis in mouse testis development. Molecular Genetics and Genomics. 2014; 289: 1157–1169.

[Google Scholar]

Liu Y, Liu M, Zhang Q, Chen G, Wang H. Selection and identification of the panel of reference genes for quantitative real-time PCR normalization in rat testis at different development periods. Toxicology and Applied Pharmacology. 2020; 406: 115243.

[Google Scholar]

Sonowal J, Lal Patel C, Dev K, Singh R, Barkathullah N, Akram Malla W, et al. Selection and validation of suitable reference gene for qPCR gene expression analysis in lamb testis cells under Sheep pox virus infection. Gene. 2022; 831: 146561.

[Google Scholar]

Zhou X, Wu X, Chu M, Liang C, Ding X, Pei J, et al. Validation of suitable reference genes for gene expression studies on yak testis development. Animals. 2020; 10: 182.

[Google Scholar]

Wu W, Hu Y, Zhang Q, Xu Y, Su W. TNFα stimulates the proliferation of immature Sertoli cells by attenuating UPS-degradation of cyclin D1 and leads to the delay of BTB maturation in pubertal rats. Andrology. 2023; 11: 575–590.

[Google Scholar]

Bustin SA, Ruijter JM, van den Hoff MJB, Kubista M, Pfaffl MW, Shipley GL, et al. MIQE 2.0: revision of the minimum information for publication of quantitative real-time PCR experiments guidelines. Clinical Chemistry. 2025; 71: 634–651.

[Google Scholar]

Luo N, Zhou Y, Chen X, Zhao Y, Hu Y. Screening the optimal housekeeping genes (HKGs) of placenta tissues by RNA-sequence and qRT-PCR throughout gestation in goat (Capra hircus). Gene. 2024; 895: 147966.

[Google Scholar]

Zhao L, Yang H, Li X, Zhou Y, Liu T, Zhao Y. Transcriptome-based selection and validation of optimal reference genes in perirenal adipose developing of goat (Capra hircus). Frontiers in Veterinary Science. 2022; 9: 1055866.

[Google Scholar]

Yang G, Yu X, Zhang Y, Luo J, Li X, Zhu L, et al. Screening and validation of stable reference genes for qRT-PCR analysis in epicauta gorhami (Coleoptera: Meloidae). Insects. 2024; 15: 942.

[Google Scholar]

Zhou Y, Li X, Zhang X, Li M, Luo N, Zhao Y. Screening of candidate housekeeping genes in uterus caruncle by RNA-sequence and qPCR analyses in different stages of goat (Capra hircus). Animals. 2023; 13: 1897.

[Google Scholar]

Schwarz AP, Nikitina VA, Krytskaya DU, Shcherbakova KP, Trofimov AN. Reference gene expression stability within the rat brain under mild intermittent ketosis induced by supplementation with medium-chain triglycerides. PLOS ONE. 2023; 18: e0273224.

[Google Scholar]

Munoz JJ, Anauate AC, Amaral AG, Ferreira FM, Watanabe EH, Meca R, et al. Ppia is the most stable housekeeping gene for qRT-PCR normalization in kidneys of three Pkd1-deficient mouse models. Scientific Reports. 2021; 11: 19798.

[Google Scholar]

Barta N, Ördög N, Pantazi V, Berzsenyi I, Borsos BN, Majoros H, et al. Identifying suitable reference gene candidates for quantification of DNA damage-induced cellular responses in human U2OS cell culture system. Biomolecules. 2023; 13: 1523.

[Google Scholar]

Franko N, Vrščaj LA, Zore T, Ostanek B, Marc J, Lojk J. TBP, PPIA, YWHAZ and EF1A1 are the most stably expressed genes during osteogenic differentiation. International Journal of Molecular Sciences. 2022; 23: 4257.

[Google Scholar]

Kim HD, Jo CH, Choe YH, Lee HJ, Jang M, Bae SG, et al. PPIA, HPRT1, and YWHAZ are suitable reference genes for quantitative polymerase chain reaction assay of the hypothalamic-pituitary-gonadal axis in sows. Animal Bioscience. 2022; 35: 1850–1859.

[Google Scholar]

Zinzow-Kramer WM, Horton BM, Maney DL. Evaluation of reference genes for quantitative real-time PCR in the brain, pituitary, and gonads of songbirds. Hormones and Behavior. 2014; 66: 267–275.

[Google Scholar]

Lyabin DN, Eliseeva IA, Ovchinnikov LP. YB-1 protein: functions and regulation. Wiley Interdisciplinary Reviews: RNA. 2014; 5: 95–110.

[Google Scholar]

Yang H, Zhang L, Liu S. Determination of reference genes for ovine pulmonary adenocarcinoma infected lung tissues using RNA-seq transcriptome profiling. Journal of Virological Methods. 2020; 284: 113923.

[Google Scholar]

Liang J, Zheng Y, Zeng W, Chen L, Yang S, Du P, et al. Comparison of proteomic profiles from the testicular tissue of males with impaired and normal spermatogenesis. Systems Biology in Reproductive Medicine. 2021; 67: 127–136.

[Google Scholar]

Zhang BB, Shen X, Li XJ, Tian YB, Ouyang HJ, Huang YM. Reference gene selection for expression studies in the reproductive axis tissues of Magang geese at different reproductive stages under light treatment. Scientific Reports. 2021; 11: 7573.

[Google Scholar]

Wang Y, Wu ZW, Fang T, Zhang YQ, Chen L, Du ZQ, et al. Identification of internal reference genes for porcine immature Sertoli cells under heat stress. Reproduction in Domestic Animals. 2022; 57: 1344–1352.

[Google Scholar]

Mahanty A, Purohit GK, Mohanty S, Nayak NR, Mohanty BP. Suitable reference gene for quantitative real-time PCR analysis of gene expression in gonadal tissues of minnow Puntius sophore under high-temperature stress. BMC Genomics. 2017; 18: 617.

[Google Scholar]

Zhao L, Yang Y, Yang H, Luo N, Li X, Zheng J, et al. Screening genes related to embryo implantation in Dazu black goats (Capra hircus) by morphological and transcriptome analyses. Journal of Animal Science. 2023; 101: skac401.

[Google Scholar]

Xu D, Liu L, Zhao Y, Yang L, Cheng J, Hua R, et al. Melatonin protects mouse testes from palmitic acid-induced lipotoxicity by attenuating oxidative stress and DNA damage in a SIRT1-dependent manner. Journal of Pineal Research. 2020; 69: e12690.

[Google Scholar]

Long XB, Shi WJ, Yao CR, Li SY, Zhang JG, Lu ZJ, et al. Norethindrone suppress the germ cell development via androgen receptor resulting in male bias. Aquatic Toxicology. 2023; 261: 106604.

[Google Scholar]

Yang R, Zhang B, Zhu W, Zhu C, Chen L, Zhao Y, et al. Expression of phospholipase D family member 6 in bovine testes and its molecular characteristics. International Journal of Molecular Sciences. 2023; 24: 12172.

[Google Scholar]

Amirian M, Azizi H, Hashemi Karoii D, Skutella T. VASA protein and gene expression analysis of human non-obstructive azoospermia and normal by immunohistochemistry, immunocytochemistry, and bioinformatics analysis. Scientific Reports. 2022; 12: 17259.

[Google Scholar]

Gao Y, Wang Z, Long Y, Yang L, Jiang Y, Ding D, et al. Unveiling the roles of Sertoli cells lineage differentiation in reproductive development and disorders: a review. Frontiers in Endocrinology. 2024; 15: 1357594.

[Google Scholar]

Yang M, Lin SY, Yang CQ, Chen YS, He ZY. Progress on SOX9 and its enhancers in mammalian sex determination. Hereditas. 2024; 46: 677–689.

[Google Scholar]

Wang S, Ma Y, Li W, Zhao Y, Gao Y, Wang S. LncRNA SNHG5/IGF2BP1/Occludin axis regulates Nd2O3 induced blood-testis barrier disruption. Environmental Pollution. 2024; 358: 124527.

[Google Scholar]