Identification of Appropriate Endogenous Controls for Circulating miRNA Quantification in Working Dogs under Physiological Stress Conditions

Simple Summary Circulating miRNAs are not only present in cells but also in the extracellular environment, especially in different biofluids including blood, and can act in a paracrine manner by facilitating a diversity of signaling mechanisms between cells. In qPCR gene expression profiling analysis, the endogenous control must be a stable gene to allow an accurate cross-sample gene expression comparison. The appropriate selection of an endogenous control is a crucial step. This research aims to select, in the miRNome, appropriate circulating miRNAs that can serve as an endogenous control. The study model are working dogs used in the search and rescue of missing persons after natural disasters. Abstract Cell-free miRNAs, called circulating miRNAs (cmiRNAs), can act in a paracrine manner by facilitating a diversity of signaling mechanisms between cells. Real-time qPCR is the most accepted method for quantifying miRNA expression levels. The use of stable miRNA endogenous control (EC) for qPCR data normalization allows an accurate cross-sample gene expression comparison. The appropriate selection of EC is a crucial step because qPCR data can change drastically when normalization is performed using an unstable versus a stable EC. To find EC cmiRNA with stable expression in search and rescue (SAR) working dogs, we explored the serum miRNome by Next-Generation Sequencing (NGS) at T0 (resting state) and T1 immediately after SAR performance (state of physiologically recovered stress). The cmiRNAs selected in the NGS circulating miRNome as probable ECs were validated by qPCR, and miRNA stability was evaluated using the Delta Ct, BestKeeper, NormFinder, and GeNorm algorithms. Finally, RefFinder was used to rank the stability orders at both T0 and T1 by establishing miR-320 and miR-191 as the best-circulating ECs. We are confident that this study not only provides a helpful result in itself but also an experimental design for selecting the best endogenous controls to normalize gene expression for genes beyond circulating miRNAs.


Introduction
MicroRNAs (miRNAs) are small non-coding RNAs that regulate target protein expression through two primary mechanisms: translational blocking in the initiation step or the elongation phase and removal of the polyA tail, promoting deadenylase activity followed by mRNA degradation [1][2][3]. Recently, it has been found that miRNAs are not only present in cells but also in the extracellular environment, particularly in various biofluids [4,5].

Animal Enrolment
The dogs included in this research work are all experts in avalanche or rubble search because they were trained at the GdF Dog Breeding and Training Center [14] (Castiglione del Lago, Perugia, Italy) and the SAR Alpine School (Passo Rolle, Trento, Italy). The dogs were physically and behaviorally tested [14] to certify their suitability for SAR training and work. Table 1 shows the dogs enrolled in the study. Table 1. SAR dog characteristics. The table shows the distribution and characteristics of the SAR dogs (n = 22); dog identification (ID) has its letter, and the gender is indicated in brackets (M)/(F).

Experimental Flowchart
SAR trial. The simulated SAR trial required dogs to find a target odor (a hidden operator and his breath, simulating a missing person) on a rubble field (30 × 35 m), within a maximum time of 15 min. Blood was collected at T0 immediately before the SAR (resting baseline) and at T1 later in the SAR (physiological stress). T0 versus T1 enables a relation between two different physiological states, where T1 comprises a multifactorial combination of stresses physiologically recovered in trained dogs [15].
The experimental framework of the study is divided into two phases.
In the first step, dog blood was drawn at T0 (n = 22) and T1 (n = 22) (dogs reported in Table 1), and total serum RNA was isolated. The total RNA of six dogs listed in At T0 and T1, 3 mL of the blood is collected as described by Guelf ysis was controlled in all serum samples because a significant source of comes from contamination by cellular-derived miRNAs resulting from rum was first visually inspected, and the presence of a pink color ind of free hemoglobin. Later, hemoglobin concentrations are evaluated i rum samples by the optical density at 414 nm (absorbance peak of free a NanoDrop™ 1000 spectrophotometer (Thermo Scientific, Scoresby, V If the absorbance of oxyhemoglobin at 414 nm exceeded a value of 0.2 classified as hemolyzed and were excluded from the experiment. T tracted from 200 µL of serum, and quantitative analysis of miRNAs described in Guelfi et al. [15]. The RNA was divided into two aliquo step and another for the second step, and both aliquots were stored at survey RNA extraction, cDNA production, and qPCR, control quality ( used. During RNA extraction, UniSp2 and UniSp4 spike-ins were add ficiency and yield of the RNA isolation, while UniSp6 was included in cDNA synthesis efficiency and control PCR inhibitor presence [15]. UniSp4, and UniSp6 spike-ins were amplified and quantified with the pair.

First Step: NGS Circulating miRNome and Stable miRNAs Selection
T0 and T1 total RNA (5 µL) were converted into miRNA NGS libr aration and next-generation sequencing to select candidate EC miRNA as described in Guelfi et al. [15]. NGS output data were utilized to ident expressed and most stable miRNAs in the T0 and T1 miRNomes. The T0 and T1 miRNome allowed us to screen, within a broad physiolog stable miRNAs. As a preliminary approach, NGS miRNAs are ordered ber to choose only those miRNAs abundantly expressed in the miRN the dogma that EC miRNAs should be highly abundant in cells and tiss are shown in the Supplementary Material (SM) Table S1. Then, only th NAs, i.e., with T0 versus T1 p > 0.05 significance value, FC value < 1.1, Figure 1. Research flowchart. In the first step, in a sample size of six dogs, a subset of miRNAs is identified by NGS as stable miRNAs. The white box shows all the phases performed in the first step: blood sampling at T0 and T1, serum RNA extraction, NGS circulating miRNA profiling, and stable miRNA selection. In the second step (gray boxes), the phases are stable miRNA qPCR validation (22 dogs) and bioinformatics analysis of qPCR output data to determine EC miRNAs' stability ranking.

First
Step: T0 and T1 Blood Sampling and Serum RNA Extraction At T0 and T1, 3 mL of the blood is collected as described by Guelfi et al. [15]. Hemolysis was controlled in all serum samples because a significant source of variation in serum comes from contamination by cellular-derived miRNAs resulting from hemolysis. The serum was first visually inspected, and the presence of a pink color indicated the presence of free hemoglobin. Later, hemoglobin concentrations are evaluated in all processed serum samples by the optical density at 414 nm (absorbance peak of free hemoglobin) using a NanoDrop™ 1000 spectrophotometer (Thermo Scientific, Scoresby, Victoria, Australia). If the absorbance of oxyhemoglobin at 414 nm exceeded a value of 0.2, the samples were classified as hemolyzed and were excluded from the experiment. Total RNA was extracted from 200 µL of serum, and quantitative analysis of miRNAs was performed as described in Guelfi et al. [15]. The RNA was divided into two aliquots, one for the first step and another for the second step, and both aliquots were stored at −80 • C until use. To survey RNA extraction, cDNA production, and qPCR, control quality (QC) spike-ins were used. During RNA extraction, UniSp2 and UniSp4 spike-ins were added to assess the efficiency and yield of the RNA isolation, while UniSp6 was included in the RT to observe cDNA synthesis efficiency and control PCR inhibitor presence [15]. In qPCR, UniSp2, UniSp4, and UniSp6 spike-ins were amplified and quantified with the respective primer pair.

First
Step: NGS Circulating miRNome and Stable miRNAs Selection T0 and T1 total RNA (5 µL) were converted into miRNA NGS libraries. Library preparation and next-generation sequencing to select candidate EC miRNAs were performed as described in Guelfi et al. [15]. NGS output data were utilized to identify the most highly expressed and most stable miRNAs in the T0 and T1 miRNomes. The comparison of the T0 and T1 miRNome allowed us to screen, within a broad physiological state, the most stable miRNAs. As a preliminary approach, NGS miRNAs are ordered by the count number to choose only those miRNAs abundantly expressed in the miRNome, according to the dogma that EC miRNAs should be highly abundant in cells and tissues [16]. NGS data are shown in the Supplementary Material (SM) Table S1. Then, only the most stable miRNAs, i.e., with T0 versus T1 p > 0.05 significance value, FC value < 1.1, and Bonferroni test value > 1, were selected. In the panel of candidate miRNAs for EC, cfa-miR-320 was included because, in previous research, it was elected as the best EC miRNA in dog serum [15]. NGS data are shown in the Supplementary Materials Table S2.

Second Step: Probable EC miRNA qPCR Validation
MiRNAs selected as probable ECs were validated by qPCR in an aliquot of serum RNA, extracted in the first step, at T0 (n = 22) and T1 (n = 22). Total RNA (10 ng) was reverse transcribed as described in Guelfi et al. [15]. QPCR amplification was carried out using 3 µL of cDNA (diluted 1:50) and appropriate primers as reported in Guelfi et al. [10] ( Table 2). Table 2. QPCR-specific pair primers. The resuspended primer mix contains both forward and reverse sequences. * Previous studies have identified miR-320 as the best EC miRNA in SAR dog serum.

Primer
Recommended for: GeneGlobe ID

UniSp2
Spike-in (RNA isolation efficiency evaluation) YP00203950 UniSp4 Spike-in (RNA Isolation efficiency evaluation) YP00203953 UniSp6 Spike-in ( PCR cycling conditions included initial denaturation for 15 s at 95 • C; followed by 40 cycles of amplification: 15 s at 94 • C for denaturing double-stranded DNA, 15 s at 95 • C for annealing, and 15 s at 70 • C extension steps. Data acquisition should be performed during the annealing/extension step. Melting curve analysis at 60-95 • C was performed to assess amplification specificity. The amplification and the results of interpreting spike-ins were performed, referring to Guelfi et al. [15]. Only reactions with efficiencies ranging from 95 to 100% were included in the subsequent analysis.

Statistics Analyses
The statistical tests utilized are described in the figure legends and methods. The graphs and statistical tests were executed employing GraphPad Prism 8 (GraphPad Software Inc., San Diego, CA, USA).

Serum RNA Quality/Quantity
The mean RIN (RNA Integrity Number) was 8.9, and the mean 260/280 ratio was 1.9 (range 1.88-2.00), which indicates high RNA purity. The amount of extracted RNA ranged from 9 to 16 ng/µL and was corrected before RT [15].

Stable EC miRNA Profiling Based on NGS miRNome Outputs
The comparison between the T0 miRNome and the T1 miRNome allowed us to identify eight stable miRNAs based on the mean stability value of count per million reads (CPM average), T0 versus T1 significance value (p-value) p > 0.05, fold change (FC) value < 1.1, and Bonferroni test value > 1. The additional information is reported in Supplementary Materials-NGS data (Tables S1 and S2). The results reported in Table 3 indicate the best EC miRNAs were selected in the first step of the experimental framework.

Bioinformatics Results of EC Expression Stability and EC Final Ranking
According to MIQE guidelines, the correct way to normalize miRNA qPCR data is by using more stably expressed EC miRNAs. The most stable candidate EC miRNAs were evaluated using four mathematical approaches: comparative Delta-Ct, BestKeeper, NormFinder, and GeNorm. As a result of each of the four bioinformatics methods, an EC miRNA (at T0 and T1) stability value (SV) is obtained (Figure 2). Additional information on the comparative ∆Ct method, BestKeeper, Normfinder, GeNorm, and RefFinder are reported in the Supplementary Materials Table S2. EC: expression stability.
Finally, RefFinder, based on the SV results of Delta-Ct, BestKeeper, NormFinder, and GeNorm, assigns an appropriate weight to the individual EC miRNAs by calculating the geometric mean of their weights for the overall T0 and T1 final ranking. RefFinder rank identified the top five miRNAs (gray box) at T0: miR-16, miR-191, miR-320, miR-93, and let-7b, while in T1, rank identified: miR-320, miR-191, let-7b, miR-16, and miR-93. MiR-191 and miR-320 (written in white in the dark gray box) maintain the best stability performances at both T0 and T1 (Table 4).

QPCR Data Quality Control and qPCR Data Normalization
To achieve an accurate identification of miRNA expression level, as a first hypothesis, the sources of technical variability were reduced by controlling the sample quality by inspecting the UniSp2, UniSp4, and UniSp6 amplification curves to remove outliers (Cq > 35). All tests did not show Cq values > 35. As the second step, following MIQE guidelines, qPCR data normalization with the most stable endogenous control miRNA should be performed to decrease the variability of experiments related to different amounts of starting RNA. The RefFinder algorithm ranked as the best EC normalizer: miR-16, miR-191, miR-320, miR-93, and let-7b (Table 4).

Discussion
This research evaluated different normalization strategies used for quantifying circulating miRNAs in the serum of working dogs trained to search for lost people after a natural catastrophe. As far as the GdF-trained dogs are concerned, SAR simulation represents a state of stress that dogs are able to regain physiologically as a result of training and exercise [15]. The SAR-learning-course includes acquisitions related to the enriched living environment in which dogs improve their problem-solving abilities, improve their physical performance, and develop better relationships by listening to and obeying the handler. In dogs, as in humans, physical activity triggers epigenetic changes in pathways associated with energy metabolism, insulin sensitivity, and more, regulating body homeostasis during exercise [22]. This adaptive physiological response is based on epigenetic regulation of gene expression mediated by miRNAs [23]. Changes in circulating miRNA levels in response to exercise allow adaptive regulation of the physiological processes including regulation of contraction and calcium signaling [24,25], bone metabolism [26,27], myocardial, skeletal muscle metabolism, repair and remodeling, mitochondrial physiology, inflammation, and angiogenesis [28]. Taking into consideration the importance of these epigenetic studies on the expression level of circulating miRNAs, a rigorous definition of the most stable normalizer is needed to reduce the number of confounding factors that influence the analytical outcome of the qPCR data.
To give relevance to the post-analytical normalization-dependent variability, this study paid special attention to stable circulating miRNA characterization by NGS profiling in a cohort of SAR dogs during two physiological states: resting and after simulated SAR performance. The NGS output data allowed simultaneous exploration and subsequent comparison of circulating miRNAs at time points T0 and T1. The stringent selection to improve the NGS data quality allowed us to identify a subset of cmiRNAs with stable expression in the T0 and T1 miRNomes that could be considered probable ECs. The selected cmiRNAs that were not differentially expressed were then validated by qPCR in a larger cohort of SAR dogs. The NGS-selected cmiRNA panel comprising cfa-miR-486, cfa-miR-16, cfa-miR-92, cfa-miR-191, cfa-miR-223, cfa-miR-423, cfa-let-7b, cfa-miR-25, and cfa-miR-93 was enriched with cfa-miR-320 because the previous studies showed that cfa-miR-320 is remarkably stable in serum [15]. In a panel of 179 circulating miRNAs, Falardi et al. [29] identified hsa-miR-320d (which has 100% homology with cfa-miR-320) as the most appropriate EC reference for circulating miRNA qPCR data normalization. The research of Falardi et al. emphasizes the stability of hsa-miR-320d in the basic physiological and physically active states due to the low association between miR-320 and the circulating miRNAs involved in exercise remodeling. The relevance of the claim is supported by the finding that the physical activity state, as opposed to the basic physiological state, undergoes epigenetic remodeling of skeletal muscle and bone [29]. The let-7 miRNA family, defined as one of the most conserved miRNA families in different species [30], is frequently deregulated in a wide variety of diseases in the dog as well as in other mammals [31]. The elevated presence of let-7b is considered strongly predictive of canine myxomatous mitral disease [32] or the poor outcome following ischemic stroke [33], and this probably supports the fact that let-7b in physiologically healthy dogs is a stable and reliable EC.
Research over two decades has demonstrated that qPCR data need to be normalized against an unregulated endogenous reference control gene to minimize variation that can mask or exaggerate biologically significant differences in gene expression levels. Extensive research on gene expression shows that using endogenous reference genes selected with bioinformatics criteria prior to conducting gene expression studies increases the accuracy of relative quantification. Conscious of this, to ensure the accurate selection of reference ECs, the putative ECs selected in the T0 and T1 miRNome were validated by qPCR in an extended sample cohort. Then the most acclaimed bioinformatics approaches, such as the Delta-Ct method, BestKeeper, NormFinder, and GeNorm, were used to choose the best EC reference gene [34,35]. In the final analysis, the RefFinder algorithm, based on the rankings of Delta CT, BestKeeper, NormFinder, and GeNorm, assigned an appropriate weight to the individual miRNAs by calculating the geometric mean of their weights and getting the T0 final EC ranking: miR-16, miR-191, miR-320, miR-93, and let-7b; and the T1 final EC ranking: miR-320, miR-191, let-7b, miR-16, and miR-93. Taking into consideration the aim of this study, to prevent an experimental bias that invalidates output data during extraction before reverse transcription and qPCR amplification, the reference exogenous miRNAs spike-in were added to control RNA quality and RT-qPCR reaction efficiency.
The optimal EC should be unrelated to biological variations, stages of disease, or treatments, and have storage stability, extraction properties, and quantification efficiency similar to the target miRNA. Therefore, before carrying out an experiment, it is imperative to first determine the best EC, keeping in mind that the EC should be expressed approximately at the same level as the target gene. Although many reports describe the ideal EC miRNAs in different species, tissues, or diseases, it is to be pointed out that, to date, there has been a lack of research describing the selection of circulating EC miRNAs in working dog serum.
A number of normalisation approaches have been suggested, but the use of one EC, selected as the best by Delta-Ct, BestKeeper, NormFinder, and GeNorm bioinformatic methods, is currently the favorite approach for the selection of the best EC to normalize qPCR data. Although some authors prefer to normalize with the average Cq value (called the mean centering method) derived from the means of different endogenous ECs [36]. In our opinion, normalization by the mean centering method should be limited to studies that analyze a large number of miRNAs. A variant of this method is the restricted mean centering method, developed for experiments in which a substantial proportion of miRNA data values are lacking. This method minimizes the technical variance in normalized expression values by utilizing a small number of normalizing ECs [20,37,38]. Other authors also suggest normalization with exogenous synthetic miRNAs added during the RNA extraction process. In our opinion, normalization with exogenous synthetic miRNAs (such as spike-ins) only correct variables linked to the extraction method, RT, or qPCR, but not for intrinsic analytical variables related to the initial sample quantity, collection method, and storage conditions. At present, a universally accepted normalization strategy is still lacking, and, unfortunately, the diverse approaches in use provide different results, with a high risk of confounding.
The differing results for the selection of the best EC, derived from Delta-Ct, Best-Keeper, NormFinder, and GeNorm bioinformatic approaches show that each method has its strengths and weaknesses, so a good strategy might be to let all of them contribute to the best EC identification [27]. We are confident that this study not only provides a helpful result in itself but also gives an experimental design for selecting the best endogenous controls to normalize gene expression, not only in circulating miRNAs.

Conclusions
Detecting gold-standard EC is crucial for producing accurate qPCR results. Our study showed the criteria for identifying appropriate endogenous controls to normalize the qPCR expression data of circulating miRNAs. The research was performed on search and rescue dogs under resting conditions and in response to SAR physiological stress. Possible circulating miRNAs as endogenous controls were selected through NGS followed by qPCR validation. The best-known algorithms Delta Ct, BestKeeper, NormFinder, GeNorm, and RefFinder were used to rank cmiRNAs in order of stability. miR-191 and miR-320 appeared to be the most stable serum endogenous controls. We are confident that this study provides a helpful result on the best-circulating EC miRNAs and provides agile and reliable criteria for identifying suitable endogenous controls to normalize qPCR expression data.
Author Contributions: S.D., G.G. and M.M.S. made substantial contributions to conception and design; S.D. and G.G. carried out the experiment and contributed to sample preparation. Acquisition of data and bioinformatics analysis were undertaken by G.G. and C.C. The manuscript was drafted by G.G. and S.D. All authors have read and agreed to the published version of the manuscript.