Male rats (6–9 weeks old, Sprague-Dawley) were purchased and allowed to acclimate for approximately one week before being used in the study. The animals were subjected to 12 h/12 h light/dark cycles in temperature-controlled rooms. Water and food were provided ad libitum. 4-Nonylphenol (100 mg/kg, Sigma-Aldrich) was dissolved in olive oil and then applied by (ip) injection, while the equivalent volume of vehicle alone was applied as a control. 4-Nonylphenol exposure was continued for 14 days. The rats were used for the experiments according to the animal care guidelines, and all animal procedures were approved by the Institutional Laboratory Animal Use and Care Committee of Xiamen University.

Testicular samples were collected and incubated in 4% paraformaldehyde (PFA) for fixation on a rotator overnight at 4°C and subsequently embedded and sectioned using a rotary microtome (Leica). For hematoxylin eosin (HE) staining, the sections were deparaffinized in xylene, rehydrated in a graded series of ethanol washes, washed with tap water, and then stained with hematoxylin and eosin.

For immunofluorescence staining, the sections were deparaffinized in xylene, rehydrated in a graded series of ethanol washes and washed with PBS. Then, antigen retrieval was performed in a microwave oven for 20 min in antigen-unmasking solution. The activity of endogenous peroxidase was blocked by 3% hydrogen peroxide in methanol, and nonspecific staining was blocked with a BSA solution for 30 min. The testicular sections were incubated with a primary antibody against androgen receptor (AR, Santa Cruz) overnight at 4°C and subsequently washed three times in PBS. Then, an IgG-Alexa 594 (Invitrogen-Molecular Probes) secondary antibody was applied for 1 h at 37°C. The sections were mounted in antifade fluorescence mounting medium with DAPI (Invitrogen). The slides were imaged using an Olympus microscope (BX51) with a digital camera.

Isolation and culture of primary sertoli cells were performed as described in a previous study (Hu et al., 2015). The testes were harvested from rats (5–10-day-old Sprague-Dawley rats) and processed for sertoli cell enrichment. Briefly, the testes were sliced using a razor blade into fine (∼1 mm wide) pieces and digested with trypsin (1 mg/ml, Gibco) and DNase I (0.8 mg/ml, Sigma-Aldrich) in F12/DMEM for 20 min at 37°C. Then, the dissociated cells were passed through a 100 µm filter and subjected to centrifugation at 800 × g for 5 min, and the pellets were resuspended in 40 ml of a solution containing 0.01% soybean trypsin inhibitor, 0.8 mg/ml dnase I, 2 mM EDTA, and 1 M glycine for 10 min. The samples were then centrifuged at 800 × g for 5 min and collected. The collected samples were then suspended and further digested in a solution of collagenase I (Sigma-Aldrich, 1 mg/ml), DNase I (Sigma-Aldrich, 0.8 mg/ml) and hyaluronidase III (Sigma-Aldrich, 1 mg/ml) at 37°C for 20 min. The cells were then centrifuged at 800 × g for 5 min and collected. The resulting cells were cultured in F12/DMEM (1:1) containing 10% fetal bovine serum (FBS). Primary sertoli cells (P2) were treated with 0, 0.2, 0.5, and 1 mM 4-nonylphenol (Sigma-Aldrich) separately in a humidified incubator under 5% CO₂ at 37°C.

After treatment, sertoli cells were collected, and the protein content of the supernatant was detected with a Thermo Scientific Pierce BCA Protein Assay Kit. The protein samples were run and separated by SDS polyacrylamide gel electrophoresis. After electrophoresis, a Trans-Blot system (Bio-Rad) was used to transfer the proteins to PVDF membranes at 90 V for approximately 1–2 h. The blots were blocked with 5% (w/v) nonfat skimmed milk for 60 min at room temperature. Then, the membranes were incubated with primary antibodies against androgen receptor (AR, Santa Cruz) overnight at 4°C. Subsequently, horseradish peroxidase (HRP)-conjugated secondary antibodies were added, and the blots were incubated for 60 min at room temperature. The results were detected with enhanced chemiluminescence (ECL) solution. The immunoblots were finally photographed using a Bio-Rad image analyzer.

Total RNA was extracted from 4-NP-treated primary sertoli cells for bulk RNA sequencing using an RNA extraction kit (TaKaRa, Dalian, China), and the A260/A280 values were quantified in the range of 1.8–2.1. Quality control of the RNA was performed, and the RNA met the sequencing requirements. Raw sequencing data have been submitted to the NCBI Sequence Read Archive under accession number GSE167912.

Quality control of the raw reads was performed with a standard procedure and verified with FastQC (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/) in this study. We used the R package edgeR to analyze the differentially expressed ceRNAs, and ceRNAs with log2 values > 1 or log2 values < −1 were selected as significantly differentially expressed (p < 0.05).

To further perform bioinformatics analysis, the Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases were used to analyze the differentially expressed circRNA cognate mRNAs, lncRNA targets and mRNAs. GO terms and KEGG pathways with corrected p values <0.05 were significantly enriched.

Tissue was immediately saturated in RNAlater until analysis. Total RNA from sertoli cells was reverse-transcribed into cDNA with a HiScript Reverse transcriptase Kit (Vazyme). Quantitative PCR was performed with Fast Start universal SYBR Green Master Mix (Roche) in a Real-Time PCR System (7,500, ABI). Gapdh transcript levels were used to normalize the relative expression of the target genes. The primer pairs (circRNAs and mRNAs) for quantitative PCR are shown in Supplementary Table S1.

The 3’ untranslated regions (UTRs) of mRNAs, circRNAs and lncRNAs were used to predict miRNA targets with miRanda and TargetScan. As circRNAs and lncRNAs contain multiple miRNA binding sites, both are considered ceRNAs that bind to miRNA competitively and indirectly influence mRNA expression. A network of ceRNA interactions, such as lncRNA/miRNA/mRNA and circRNA/miRNA/mRNA interactions, was constructed with local Perl scripts.

The data analyses were performed with GraphPad Prism 7. Probability values of < 0.05 were considered to indicate statistical significance.

The morphology of 100 mg/kg/day 4-NP-treated (for 14 days) rat testes groups was examined by HE staining. Histological examinations of 4-NP-treated rats showed severe damage and arrest of spermatogenesis (Figure 1A) and testes that were significantly smaller than those in the control group (Figure 1B). Compared with the seminiferous tubules in the control group, which contained the full complement of germ cells, those in the 4-NP-treated group had smaller diameters, and only morphologically normal sertoli cells and undifferentiated spermatogonia cells were observed in the 4-NP-treated group (Figure 1A) and the seminiferous epithelium was also completely lacking in differentiated germ cells; at the basal membrane, only sertoli cells and undifferentiated spermatogonia were observed (Figure 1A). Moreover, although immunofluorescence analysis revealed that the signal intensity of AR in 4-NP-treated groups was weaker than that of the controls, the numbers of AR-positive cells (sertoli cells) in tubules were not significantly lower in 4-NP-treated rat testes than in control testes in vivo (Figure 1C). This result suggested that 4-NP maybe decreased the expression of AR and impair the function of sertoli cells.

Then, to evaluate whether 4-NP impaired the function of sertoli cells and reduced the expression of AR, we investigated the effect of 4-NP on rat primary sertoli cells in vitro and treated the cells with differing concentrations of 4-NP (0, 0.2, 0.5, and 1 mM) at 24 h. When compared with no 4-NP-treated rat primary sertoli cells, there were no morphological alterations in sertoli cells exposed to 4-NP (0.2 and 0.5 mM) and no apparent change in the numbers of sertoli cells exposed to 4-NP (0.2 and 0.5 mM), but cell death and decreased numbers of cells were observed in 1 mM 4-NP-treated group (Figure 2A). Furthermore, we investigated the changes in the protein levels of AR in primary sertoli cells associated with 4-NP exposure. Sertoli cells were exposed to 4-NP (0, 0.2, 0.5, and 1 mM), and the AR levels in the sertoli cells were determined by western blot analysis (Figure 2B). The protein levels of AR were significantly lower in the 0.2, 0.5, and 1 mM 4-NP-treated cells than in the control cells (Figure 2B). These findings revealed that 4-NP reduced AR protein levels in rat sertoli cells and impaired sertoli cell function.

To uncover the ceRNA mechanism underlying the effect of 4-NP on sertoli cells in vitro, systematic ceRNA/small RNA sequencing (ceRNA/small RNA-Seq) was performed. Raw sequencing data have been submitted to the NCBI Sequence Read Archive under accession number GSE167912. We included three controls and three 4-NP (0.5 mM)-treated sertoli cell samples in the miRNA-Seq experiment, and we included two controls and two 4-NP (0.5 mM)-treated sertoli cell samples in the ceRNA-Seq experiment. The ceRNA/miRNA-Seq data were analyzed by bioinformatics methods. Strict filtering methods were used to compare the 4-NP group and the control group (fold change: >2 or < 0.5; p < 0.05) (Supplementary Figures S1A, S2A, S3A, S4A). We identified 54 upregulated and 82 downregulated miRNAs (Supplementary Figure S1B), 465 upregulated and 419 downregulated circRNAs (Supplementary Figure S2B), 74 upregulated and 191 downregulated lncRNAs (Supplementary Figure S3B), and 1,747 upregulated and 2,659 downregulated mRNAs (Supplementary Figure S4B).

FastQC quality control was performed. A Q-score higher than 10 (error rate < 10%) was required for each replicate for circRNAs (Supplementary Figure S2D) and lncRNAs (Supplementary Figure S3E). The variability of the ceRNA sequences in 4-NP-treated samples relative to controls is displayed in box plots (Supplementary Figures S2C, S3C, S4C). In addition, the correlation coefficients of the replicated samples had good consistency (Supplementary Figures S1C, S3D).

For miRNAs and ceRNAs. The length distribution of small RNAs is shown in Supplementary Figure S1D, and the precision of the RPKM analysis for the mRNA-Seq experiments is shown in Supplementary Figure S4D. Moreover, based on the differentially expressed miRNAs (Figure 3A), circRNAs, lncRNAs, and mRNAs, a tree was generated by MeV_4_9_0 cluster analysis. The tree clearly separated 4-NP-treated samples from controls (Figures 3B–D). These results indicated that ceRNA expression could be used to robustly separate control and 4-NP-treated samples.

We determined the functional categories of differentially expressed miRNAs, circRNAs, lncRNAs and mRNAs between 4-NP-treated rat primary sertoli cells (0.5 mM, in vitro) and without 4-NP treatment via GO analysis to further understand the biological processes that may be mediated by these differentially expressed molecules in 4-NP-treated cells. For differentially expressed miRNAs, 827 GO terms in the category of biological processes were significantly enriched.

As shown in Figure 4A, the top 10 biological process terms with significant enrichment were the developmental process, regulation of developmental process, cell differentiation, cellular developmental process, tissue development, phosphate-containing compound metabolic process, protein phosphorylation, regulation of cell communication, macromolecule modification, and transport terms.

A total of 100 differentially expressed circRNAs were annotated with GO terms. The top 10 biological process terms with significant enrichment are shown in Figure 4B. They were the insulin secretion, response to glucocorticoid, axon extension, I-kappaB kinase/NF-kappaB signaling, negative regulation of neuron differentiation, cell death, negative regulation of bicellular tight junction assembly, response to testosterone, regulation of cell morphogenesis, and negative regulation of mitogen-activated protein kinases (MAPK) cascade terms. Many of the biological process terms were related to sertoli cell development, such as the response to testosterone and negative regulation of the MAPK cascade terms.

A total of 100 differentially expressed lncRNAs were annotated with GO terms. As shown in Figure 4C, the top 10 biological process terms with significant enrichment were the RNA phosphodiester bond hydrolysis, response to drug, phosphatidylinositol 3-kinase signaling, hormone activity, cytokine receptor binding, protein tyrosine/serine/threonine phosphatase activity, ER-associated ubiquitin-dependent protein catabolic process, negative regulation of cell aging, positive regulation of protein phosphorylation, and G-protein coupled receptor signaling pathway terms. The most abundant terms were those associated with the functions of sertoli cells.

There were also one hundred GO terms that were significantly enriched for differentially expressed genes. As shown in Figure 4D, the top 10 biological process terms with significant enrichment were the cell adhesion, negative regulation of cell proliferation, response to estradiol, regulation of cell growth, aging, single organismal cell-cell adhesion, response to steroid hormone, DNA methylation involved in gamete generation, response to drug, and response to estrogen terms. The most abundant terms were those associated with the functions of sertoli cells and the medicinal properties of 4-NP.

One hundred fifty-six KEGG pathways (FDR < 0.05) were associated with differentially expressed miRNAs. As shown in Figure 5A, the top 10 enriched pathway terms were the Focal adhesion, Wnt signaling pathway, Hippo signaling pathway, MAPK signaling pathway, Gap junction, EGFR tyrosine kinase inhibitor resistance, tight junction, mTOR signaling pathway, Notch signaling pathway, and Adherens junction terms. These results indicated that 4-NP injured the blood-testis barrier (BTB) in sertoli cells in vitro.

Sixteen KEGG pathways (FDR < 0.05) were associated with host genes of differentially expressed circRNAs. As shown in Figure 5B, the top Ten enriched pathway terms were the Chemokine signaling pathway, Insulin signaling pathway, ABC transporters, Ether lipid metabolism, FoxO signaling pathway, TGF-beta signaling pathway, Insulin resistance, mTOR signaling pathway, VEGF signaling pathway, and cGMP-PKG signaling pathway terms.

Seven KEGG pathways (FDR < 0.05) were associated with differentially expressed lncRNAs. As shown in Figure 5C, the top 10 enriched pathway terms were the Metabolism of xenobiotics by cytochrome P450, Biosynthesis of unsaturated fatty acids, Fatty acid elongation, Drug metabolism - cytochrome P450, Glutathione metabolism, Asthma, Circadian rhythm, Inflammatory bowel disease (IBD), Tyrosine metabolism, and Phenylalanine metabolism terms.

Fifty-one KEGG pathways (FDR < 0.05) were associated with differentially expressed mRNAs. As shown in Figure 5D, the top 10 enriched pathway terms were the Focal adhesion, ECM-receptor interaction, Arrhythmogenic right ventricular cardiomyopathy (ARVC), Malaria, Hypertrophic cardiomyopathy (HCM), Complement and coagulation cascades, Amebiasis, Dilated cardiomyopathy (DCM), PI3K-Akt signaling pathway, and Lysosome terms.

A core set of mRNAs (Ar, Atf6 and Cbp) and circRNAs (circ673, circ1377, circ1789 and circPTEN) identified as differentially expressed by ceRNA/miRNA-Seq (Figure 6A) were validated by real‐time quantitative polymerase chain reaction (RT‐qPCR). As shown in Figure 6A, the RT-qPCR results were highly consistent with the RNA-sequence data. In addition, we performed sanger sequencing of RT-qPCR products in PMD-19T plasmids, and the results shown in Figure 6B confirmed that the products were circ673, circ1377, circ1789 and circPTEN. Melting curve analysis and agarose gel electrophoresis were used to check for the specificity of the RT-qPCR products.

Based on our sequencing data, we used Cytoscape_V2_8_3 to predict a global ceRNA network. As shown in Supplementary Table S2 and Figure 7, we conducted further characterization of circRNA-associated ceRNA networks and performed GO analysis of the biological processes. Notably, we observed that the target genes of this module were related to biological processes and molecular functions of sertoli cells, such as cell adhesion and bicellular tight junction formation. In addition, we selected circ673, circ1377, circ1789, and circPTEN from the global ceRNA network and predicted their ceRNA interactions based on our RNA-Seq data by using Cytoscape_V2_8_3. As shown in Supplementary Table S3, both circ673 and circ1789 were predicted to regulate the expression of AR through ceRNA interactions. Both circ1377 and circPTEN were predicted to regulate the expression of ATF6 through ceRNA interactions.