Semen sampling as a simple, noninvasive surrogate for prostate health screening

ABSTRACT The detection rates for prostate cancer (pCa) by invasive biopsy are high, fully justifying its use in confirmatory testing. False-positive results of prior, relatively insensitive screening tests, however, can lead to expensive and often unnecessary surgery. Several reports have suggested the potential use of the ejaculate to screen for prostate conditions. Hitherto, the potential impact of sterilization on the diagnostic potential of seminal plasma screening has not been examined. Herein, we report cellular and molecular comparisons of semen samples obtained from normal (N = 5), vasectomized (N = 5) and prostate pathology patients (N = 4; confirmed by a biopsy) that were centrifuged over 60% PureSperm cushions. Non-penetrating cells were washed prior to immunocytochemistry with prostatic epithelial cell markers including PSMA, NKX3.1 and CD24. KRT18 was used to highlight epithelial cells in these samples. RNA sequencing was then used to identify differentially expressed small RNAs associated with vasectomy and prostate pathology. Specific gene transcripts were confirmed by RT-qPCR. PMSA+/KRT18+, CD24+/KRT18+ and NKX3.1/+KRT18+ cells were observed, albeit infrequently in most processed semen samples by indirect immunocytochemistry. Targeted RT-qPCR supported their enrichment, along with their putative designation as prostatic luminal cells. Small RNAs in seminal plasma were highly heterogeneous, with tRNAs and miRNAs being the dominant forms. Hsa-miR-143 and hsa-miR-199 were among the most prominent of the differentially expressed miRNAs upregulated in samples with prostate pathology but not vasectomy. The targets of these small RNAs illustrate biological processes involved among others in transcription regulation and collagen metabolism. Our outcomes strongly support an appraisal of selected biologically meaningful small RNAs of ejaculate semen for prostate health screening. A long-term goal would be a simple, routine, noninvasive test for monitoring prostate health, potentially among younger men.


Introduction
Prostate cancer (pCa) is the most common noncutaneous cancer of males, responsible for almost 25,000 deaths annually in the US alone (Siegel et al. 2018). Diagnosis, screening and surveillance of the disease relies on the detection of prostate-specific antigen (PSA) in the blood (Hendriks et al. 2017), followed up, where indicated, by highly invasive prostate biopsies. The efficacy of PSA screening, however, is equivocal (Pinsky et al. 2017) and the results of a recently published randomized clinical trial showed that single PSA testing was a poor indicator of likely survival (Martin et al. 2018). Based on the discovery of alleles that may confer an increased susceptibility to the disease (Schumacher et al. 2018), alternative DNA screening modalities using the saliva have been proposed (Oxtoby 2018). DNA tests, however, cannot detect the presence of overt or indolent disease. As almost a third of seminal fluid originates from the prostate (Mann 1954), ejaculate semen could offer an indirect and more appropriate window into extant prostate function.
The existence of cell-free RNA in ejaculated semen has recently driven a renewed interest in the molecular investigation of prostate health (Babiker et al. 2007;Huang et al. 2015). These small RNAs (<200 nt) include miRNA (~22 nt) and piRNA (~30 nt) and as regulators of gene expression, they are closely involved in the differentiation of the epididymal epithelium and acquisition of male fertility (Papaioannou et al. 2010;Bjorkgren et al. 2012). Some miRNAs originate in the prostate (Liu and Tang 2011;Ostling et al. 2011;Huang et al. 2015) and can be indicative of prostate health (Selth et al. 2014). Any tests using the ejaculate to detect pCa should consider the higher frequency of vasectomy among older men. This study hypothesized that when vasectomy procedures are taken into account, small RNAs derived from ejaculate semen should be informative of prostate health.

Semen sample characteristics
Our study employed samples from three different groups, including five normozoospermic (N) controls with no known prostate pathologies, five samples from sterile vasectomized men (S) and four men with confirmed prostate pathologies (P). All control samples were derived from men below 40 years old, with an average age of 28 years. The average ages for the vasectomized and prostate pathology groups were 55 and 64 years, respectively. The volume of semen and sperm concentration were on average approximately 3 mL and 80 million sperm/mL, respectively, for the N group. All samples from the S group had an average volume of 2.5 mL and fewer than 100,000 non-motile sperm/mL, indicative of a successful vasectomy (Sharlip et al. 2012). The average volume and sperm concentration for the P group were 2.5 mL and 26.4 millions/ml, respectively.
We first examined the round cell content of ejaculates from N and S samples, recovered from the 60% density fractions using antibodies to PSMA, NKX3.1 and CD24. Co-labeling of these markers with antibodies to KRT18 enabled us to distinguish between epithelial cells of likely prostatic origin (PSMA + / KRT18 + , NKX3.1 + /KRT18 + and CD24 + /KRT18 + ) and other cells, such as CD24 + /KRT18and CD24 -/ KRT18 + , present in the ejaculate of less determinate origin ( Figure 1). PSMA + /KRT18 + and NKX3.1 + / KRT18 + cells displaying an angular, central nucleus with two prominent cytoplasmic lobes, suggestive of prostatic epithelial cells, were observed in these samples. These cells were often seen alongside two types of CD24 + cells, which were either KRT18or KRT18 + , the former likely to be neutrophils and/or round cells of indeterminate origin. PSMA + cells were also detected in the P samples examined. For example, the cells illustrated in Figure 1E, F that may have been undergoing mitosis prior to isolation have two diffuse nuclei bounding a cytoplasmic space. We estimated that between 100 and 500 PSMA + cells were present in each sample. The virtual absence of sperm in the S group was noted, with DAPI-only and/or immunesignal labeling, suggesting the presence of other unidentified cell types possibly including necrotic cells. Single-channel images of the field are presented (Supplementary Figure 1).

Efficacy of cell enrichment in semen samples using density gradient centrifugation (DGC)
The effectiveness of cell enrichment by DGC was tested by extracting total RNA from the recovered 60% fractions (N and S groups) and pellets (N samples only) using RT-qPCR for PSMA and CD24 transcripts, normalized against endogenous control RNAs ( Figure 2). PCR amplification for CD24 and PSMA generated from the 60% fractions required statistically significantly fewer cycles (higher quantities) than those generated from the pellet fractions (one-way ANOVA), suggesting enrichment of putative prostatic cells in this fraction. Together with the immunocytochemical results from unfractionated semen versus the 60% fractions, our molecular data suggest that DGC enriched the round cell content of ejaculate semen, including potential cells of prostatic origin.

Determining seminal plasma small RNA biomarkers by RNA sequencing
A multidimensional scaling (MDS) plot for cell-free seminal small RNAs identified by RNA-Seq showed that samples could be satisfactorily resolved into three distinct groups/phenotypes (N, S and P; Figure 3). Approximately, half (48%) of the sequencing reads mapped to the human genome, of which 78% were assigned to small RNAs (Supplementary Figure 2; Supplementary Table 1). The remaining reads corresponded to repetitive elements, including LINEs, HERVs, MERs and 7SL RNA, while fragments of protein-coding mRNAs were also present. TRNAs, miRNAs, YRNAs and RNA fragments deriving from the 28S large ribosomal subunit (LSU) dominated the three groups.
The presence of piRNA and tRNA subclusters distinguished between P and N samples, respectively .1 (C) and CD24 (D). All except the P samples (E, F) were co-labeled (green channel) for KRT18 and DAPI (blue channel) was included throughout. Panel A shows two PSMA + /KRT18 + features each with two nuclei indicating either individual, closely juxtaposed cells or possibly two bi-nucleate cells surrounded by numerous, much smaller sperm nuclei. Panel B shows a similarly PSMA + /KRT18 + sized cell, with triangular-shaped nucleus and no sign of sperm. Panel C shows an example of an NKX3.1 + / KRT18 + cell and panel D shows a low powered field with CD24 + /KRT18 + cells (short barbed arrowheads) and CD24 − /KRT18 + cells (long barbed arrowheads). The latter may be neutrophils. The panel also shows an example of a CD24 + /KRT18 − /DAPI − 'entity' that may be a cytoplasmic fragment (broken arrow). The remaining panels show examples of PSMA + cells from P samples with (E/E') and without (F/F') signs of sperm. Note the presence of bright PSMA + vacuoles in the cytoplasm of these cells, more clearly visible in the higher-powered images shown in E' and F'. Green signals in fields obtained from co-labeled cells were sometimes overwhelmed by the red signals but as can be seen in Supplementary Figure 1 (single channels), were usually present. DAPI-only signals indicated by asterisks may represent cellular fragments remaining after necrosis or apoptosis. The scale bars throughout are 20 µm.  Table 1 are included in the more extensive heatmap (Supplementary Figure 3).
The 20 and 100 most significantly differentially expressed (DE) RNAs between the N and P groups (>4-fold change (FC); False Discovery Rate (FDR) <0.05) are shown, respectively, in Figure 4 and Supplementary Figure 3. Among others, we looked for potentially informative features that were not DE exclusively in N vs. S comparison, but in N vs. P and/ or S vs. P comparisons. Several features met this requirement (Table 1), including four miRNAs (mir-199a, mir-143, mir-127, mir-4800) and tRNAs such as tRNA-Gly-GGY and tRNA-Ile-ATT. Candidate seminal biomarkers potentially informative of the underlying prostate pathology (P group), regardless of sterilization, were apparent. The interplay between the significantly downregulated small RNAs in the P group, such as miR-127 and miR-4800, and the significantly upregulated small RNAs, including mir-199a and mir-143, along with tRNAs, such as the tRNA-Ile-ATT variants, could indicate, in seminal plasma, the transition from a normal to an abnormal prostate.
Our study also reports the ontologies (biological processes) of the target genes potentially regulated by the significantly differentially expressed miRNAs of the P and N groups. The more highly expressed miRNAs in the P group include mir-187, mir-199a and mir-143. The corresponding gene targets are involved in protein phosphorylation (FDR 1.4E-4), positive regulation of transcription (FDR 2.9E-4) and WNT signaling pathway (FDR 3.1E-2). The more highly expressed miRNAs of the N group include among others mir-152, mir-127 and mir-  Table 2). Their associated gene targets are involved in the negative regulation of transcription (FDR 1.3E-6), extracellular matrix organization (FDR 4.5E-6), collagen catabolic process (FDR 3.3E-5) and cell adhesion (FDR 4.2E-4) ( Table   2), illustrating the importance of biological processes distinguishing normal from abnormal prostate pathology in seminal plasma. All but one (miR-7153) of the above microRNAs have been linked with prostate origin/pathologies (Casanova-Salas Log 2 fold-change (>4-fold difference) relating to mean N vs. S (vs P) small RNA expression is shown. Although piRNAs were detectable in S and P samples (see Supplementary Figure 2), their expression was lower compared with N samples. Note the relative absence of piRNAs in S samples and the strong expression of miRNA-143 and miRNA199a-2 in P samples. The panel (upper right) indicates the log 2 FC (cpm/ mean(cpm)), ranging from dark blue (−12.51) to purple (12.51). Hierarchical clustering analysis was employed to classify the small RNAs in each category (top right). et al. 2014; Theodore et al. 2014;Ma et al. 2017;Arai et al. 2019;Fan et al. 2020) and the presence of miR-7153 in sperm has been described in relation to unexplained asthenozoospermia (Heidary et al. 2019).

Discussion
Prostate cancer accounts for more than a quarter of all new male cancer diagnoses, equating to almost 48,500 cases annually in the UK and predicted to increase by  The corresponding biological processes are listed alongside their respective P-values and FDRs and the number of genes involved per ontology along with their proportion (%).
~12% by 2035 (Smittenaar et al. 2016;Ridgway and Aning 2021). Clinical detection is mostly reliant on PSA levels, which even in conjunction with digital rectal examination can be equivocal (Halpern et al. 2018). In this study, we provided evidence for the presence of prostatic cells in the human ejaculate and their enrichment following DGC was demonstrated by RT-qPCR. Small RNAs in seminal plasma were investigated using RNA sequencing. Despite the age differences between N and P groups, none of the differentially expressed seminal plasma microRNAs overlapped with the age-related seminal plasma microRNAs reported by Paoli et al. (2019). This supports the view that the RNA differences reported in the current study were prostate rather than age related. Gardiner et al. (1996) followed by Barren et al. (1998) reported cells of likely prostatic origin in the ejaculate, respectively, by immunocytochemistry and flow cytometry, suggesting that semen sampling could be useful for prostate screening. Moreover, higher levels of PSA mRNA were reported from the ejaculates of pCa patients, presumed to have been derived from higher loads of prostatic cells in pCa samples (Clements et al. 1999). Barren et al. (1998) did not find PSA to be a particularly informative marker but reported that PSMA + cells were more common in pCa than control samples. In the current study, we took an immunocytochemical approach to characterize the presence of nonsperm cells in the N, P and S samples, with the vasectomized samples setting a background for an ejaculate 'free' of sperm-/testis-derived RNA. Our data based mainly on PSMA labeling suggested that while present, prostatic cells in the ejaculate were relatively rare. The success of RNA isolation from such uncommon cells depends on their enrichment, their PSMA gene expression levels and the variation of gene expression among samples, which together made it difficult to define thresholds for the PSMA gene expression based on simple cell isolation technologies, such as DGC, that could assess the development of prostate pathologies with confidence. We therefore turned our attention to the small RNAs of seminal plasma, which appeared to have greater potential for detecting prostate pathologies. Moreover, with more cost-efficient molecular-based technologies, screening for prostate pathologies would be easier to apply and scale. Assuming that a distinct small noncoding RNA (sncRNA) signature arising from sterilization might interfere or complicate the detection of prostate pathology, we also included samples from recently vasectomized men where small RNAs of testicular and epididymal origin would be reduced or absent.
In agreement with previous reports (Selth et al. 2014;Vojtech et al. 2014), the majority of our mapped sequences belonged to tRNA and rRNA fragments, while the remainder comprised miRNAs, piRNAs, YRNAs and U2 RNAs. U2 RNA is essential for RNA processing through the spliceosome (Sashital et al. 2007). The 28S RNA subunit of the 80S complex was strongly represented in our S and P samples. Small RNAs of this subunit were most likely breakdown products arising from cell turnover in the reproductive tract, which may be elevated in inflammatory responses that follow vasectomy or accompany underlying but unidentified pathology of the prostate (Bhindi et al. 2017). Seminal tRNAs could have originated from 5ʹ/ 3ʹ terminal fragments or endonuclease-derived tRNA halves, although as previously reported (Olvedy et al. 2016;Magee et al. 2018) more likely the former. Our analysis revealed a number of tRNAs strongly upregulated in the P samples, although many of the same tRNAs were also upregulated in S samples making their suitability as markers for pCa screening less appealing. Focusing on the seminal plasma miRNAs elevated in P, but not in N or S samples, two species (miR-143 and miR-199a) met this requirement. Both are reported to be tumor suppressing in cancer pathologies, with ectopic expression inducing apoptosis He et al. 2016). Hence, normal expression of these molecules in healthy cells could be acting in a protective capacity to help prevent cancer progression. Alternatively, prostate cancer cells may secrete these small RNAs into the extracellular space, increasing their abundance in the seminal plasma. Either process may underlie their relatively greater abundance in the seminal plasma of P patients.
Other studies have taken different approaches to investigate seminal plasma at the mRNA, small RNA and protein levels in prostate cancer/pathologies (Drabovich et al. 2014). A review by Sapre and Selth (2013) reported mainly RT-qPCR studies of seminal plasma and blood serum, revealing potential prostate cancer microRNA biomarkers. There is a little correspondence, however, among studies reporting on seminal plasma RNA. For example, two studies (Selth et al. 2014;Vojtech et al. 2014), both of which employed RNA-Seq, reported one common and possibly informative marker (miR-375). This biomarker, however, was not reported in our study. Following on from their previous work, Roberts et al. (2015) reported a combination of molecular biomarkers, including serum PSA, microRNAs and protein biomarkers for the detection of prostate cancer (Roberts et al. 2015). A more recent study employed a larger cohort of patients using an RT-qPCR panel to reveal seminal plasma biomarkers that could be associated with the aggressiveness of prostate cancer (Farran et al. 2018), with Barceló (2020) following a similar approach testing certain microRNAs in men with moderately altered PSA (Barceló et al. 2020). The divergent results of the candidate seminal plasma microRNA biomarkers may be explained by the relatively low numbers of prostate cancer/pathology samples screened by high throughput technologies that limit our ability to resolve differences in their respective microRNA expression profiles, distinguishing healthy from prostate pathologies. This divergence could also be explained by the different experimental designs and methodological approaches employed in each study, including small RNA isolation approaches, library construction and sequencing methods.
In-depth analysis of the gene ontologies of the miRNA gene targets from P and N samples supported the findings of the current study, indicating their likely participation in important biological processes that could be deregulated in the transition from a normal to an abnormal prostate pathology. These processes involve protein phosphorylation, positive transcription regulation and WNT signaling (Bisson and Prowse 2009). The downregulated seminal plasma miRNAs of P patients suggest an association with the negative regulation of transcription, extracellular matrix organization (Walker et al. 2018), collagen catabolism and cell adhesion (Burns-Cox et al. 2001) in their prostates. These processes potentially reflect the extracellular microenvironment of the seminal plasma in prostate cancer patients. A changing microenvironment may indicate a deregulated extracellular matrix with potential gain/loss of cell adhesion properties (Janiszewska et al. 2020). Collagen catabolism may lead to collagen accumulation, which is an important process in cancer progression and invasion (Xu et al. 2019). Small RNAs regulating these biological processes may be informative in our understanding and useful as potential biomarkers of prostate pathology progression.
Our study showed a number of highly expressed microRNAs, including miR-143 and miR-199a among P samples, which corresponds with the ontological description of their target genes. Moreover, tRNAs, such as tRNA-Gly-GGY and tRNA-Ile-ATT, which have not been reported previously, could also meet the requirements for a biomarker of prostate health. Our study did not consider the possibility that these biomarkers arose from other extracellular vesicles or other tissues including the epididymis and Cowper's (bulbourethral) gland. If not of prostatic origin, however, prostate pathology may still have led to either their appearance or elevation in the ejaculate where they were identified by our RNA-Seq approach accordingly.
Changes in the small RNA biomarkers of ejaculate semen could alert practitioners of early pathological changes in prostate pathology well before overt and possibly more clinically serious disease has developed. The validation of these biomarkers on a large sample cohort will be necessary to confirm their utility and functionality, however. Semen is normally obtained without an invasive procedure and can be repeatedly sampled without harm and relatively inexpensively. Such a scenario would be more difficult in the case of elderly men with potential prostate pathologies. However, the kind of test envisaged here would be more applicable to younger men who can ejaculate normally and who may be at greater risk based on family history or PSA screening but who have not yet developed overt signs of pathological change or are currently asymptomatic.

Subjects and specimens
All semen samples were obtained between 2014 and 2015 by masturbation after 3-5 days of sexual abstinence. Control (normal; N) semen samples (n = 5) were obtained from healthy young volunteers with sperm concentrations and progressive motility of at least 20 × 10 6 and 50%, respectively, and from post-vasectomy (sterile; S) semen samples (n = 5) from patients undergoing checks for sperm clearance (no earlier than 8 weeks following the procedure).
Vasectomized samples provided a background view of the ejaculate, free of sperm and testis-derived RNAs and to check for potential interference with candidate molecular biomarkers of prostate health. Semen samples (n = 7) from patients with prostate pathology (P) were obtained. These samples were derived from two patients diagnosed with prostate cancer (Gleason 7), one patient with a diagnosis of premalignant lesions (prostate intraepithelial neoplasia) and one patient with rare atypical glands. Three samples were unsuitable for further analysis and were excluded. All samples were permitted to liquefy at 37°C and processed appropriately. Samples were treated either with Quinn's Advantage Sperm freezing medium (SAGE) or SpermFreeze (Origio, Denmark) prior to freezing in LN2 or shipping.

Immunocytochemistry
Semen samples (N, S and P) were thawed and incubated at 37°C for at least 30 min prior to centrifugation on a 60% cushion of Puresperm equilibrated with PureCeption Isotonic Solution (SAGEMedia, UK) mixed with Quinn's Sperm Washing Medium (SAGE, UK) at 400 × g for 20 min. The 60% DGC layers and the pellets were collected and washed twice with Quinn's Sperm Washing Medium at 500xg for 10 min.
For immunocytochemistry, cells recovered from 60% cushions were cytospun at 500 × g for 15 minutes on to poly-l-lysine coated slides (VWR, UK). Cells were fixed in 3:1 methanol/acetone for 15 min, washed with Phosphate Buffered Saline (PBS) and incubated in blocking solution (2% goat serum, 2% serum albumin, 0.1% Triton X-100 in PBS) for 30 min at room temperature (RT). Antibodies were applied to slides in a humidified chamber for 60 min at RT (Supplementary Table 3). Slides were subsequently washed in PBS and incubated for 30 min with appropriate fluorescence-conjugated secondary antibodies in blocking solution at RT. Slides were finally washed in PBS and mounted with an antifade-containing polyvinyl alcohol-mounting medium (DABCO, UK). Fluorescence microscopy was undertaken with a Leica LEITZ DMRB fluorescence microscope. The images were taken with a cooled CCD camera (Orca-03 R, Hamamatsu) controlled by Smart Capture 3 software (DSUK, UK).

Cellular RNA isolation and quantitative reverse transcription PCR (RT-qPCR)
RNA extraction was performed as described by Goodrich (Goodrich et al. 2013). Semen samples were centrifuged by DGC and washed twice in PBS. Cells from the 60% layer and pellets were lysed in RLT buffer (QIAGEN, UK) with 1:70 β-mercaptoethanol. Trizol (Ambion, UK) and chloroform (200 µL; VWR, PA, USA) were added to the mixture, followed by RT incubation after each step. The mixture was centrifuged in a benchtop microfuge at 13200 × g for 15 min (4°C). The upper layer was transferred into a fresh tube and mixed thoroughly with 100% ethanol. The mixture was then processed using the RNeasy mini kit (QIAGEN, Germany). RNA was treated with the Turbo DNase (Ambion, UK), followed by quality control and quantification using a Nanodrop ND-1000 (Thermo Fisher Scientific, MA, USA) and Agilent 2100 bioanalyzer (Agilent, CA, USA), respectively. RNA was then reverse transcribed using a mix of random primers and oligodT (Ambion, UK) and SuperTaq DNA polymerase according to the manufacturer's instructions (Thermo Fisher Scientific, MA, USA).
Primers for PSMA and CD24 were used to test for round cell enrichment (including prostatic cells) using RT-qPCR. A TaqMan probe-based assay (Primer Design, UK) was employed on an LC-480 thermal cycler (Roche, Switzerland). The same baseline was employed across all samples, using a maximum cycle cutoff of 35. GAPDH, RPLP0 and β-actin were used as endogenous controls and normalization was conducted using the GeNorm algorithm (Vandesompele et al. 2002). Oligo sequences are available (Supplementary  Table 4).

Small RNA isolation and RNA sequencing (RNA-Seq)
Spermatozoa, round cells and cell debris were cleared from liquefied semen samples by 10 min centrifugation at 1000xg followed by 20,000 × g for 30 min (4°C) to clarify the plasma. Cleared seminal plasma was stored at −80°C. Frozen seminal plasma was thawed at 37°C and mixed with TRIzol and an equal volume of isopropanol (Sigma, MI, USA) was added to and mixed thoroughly with the upper layer. The mixture was then centrifuged at full speed on a benchtop microcentrifuge for 15 min (4°C). The RNA pellets were washed twice with 70% ethanol and dissolved in 50 µl of RNase-free water containing RNase-out (20 U). RNA purity and quality were checked using the high sensitivity 6000 RNA Pico Kit (Agilent, CA, USA), following the manufacturer's instructions.
Isolated RNA was subjected to small RNA-Seq. Briefly, 150 ng of total RNA underwent NEBNext library construction (New England Biolabs, MA, USA), according to the manufacturer's protocol. We mainly focused on the miRNA, tRNA (18-23 nt) and piRNA (30-34 nt), which were size selected by gel electrophoresis (LabChip XT; Marshall Scientific, USA) and converted to cDNA, followed by adaptor ligation and purification. The cDNA was PCR enriched, using index PCR primers. After 15 PCR cycles, libraries were purified on QIAquick columns and the average fragment length and quality of the insert was evaluated on a Bioanalyzer 2100 (Agilent, USA). The cDNA libraries were quantified by RT-qPCR alongside a standard of known concentration. Equimolar library concentrations were pooled and sequenced on an Illumina NextSeq 500 (Illumina, CA, USA).

Bioinformatics, statistics and ontological analysis
Reads underwent quality control using FastQC (Andrews 2010) and adapter removal with TrimGalore (Krueger 2015). The reads were aligned to the human genome (hg38) using Subread v.1.6 (Liao et al. 2013). FeatureCounts function assigned reads into features using refseq (Pruitt et al. 2005), piRNA (Sai Lakshmi and Agrawal 2008) miRNA (Kozomara and Griffiths-Jones 2011) and repeatmasker (Smit et al. 2013(Smit et al. -2015. EdgeR performed normalization and differential expression analysis using Fisher's exact test (Robinson et al. 2010), revealing significant intergroup differences of small RNAs. One-way ANOVA with Tukey's post-hoc test was employed for statistical analysis of the RT-qPCR using GraphPad Prism 7.0. Small-RNA targets were predicted by the miR database, where 50% complementarity was considered (Chen and Wang 2020). DAVID was employed for ontological analysis of these small RNA targets (Huang et al. 2009).