Saliva microRNA Differentiates Children With Autism From Peers With Typical and Atypical Development

Objective: Clinical diagnosis of autism spectrum disorder (ASD) relies on time-consuming subjective assessments. The primary purpose of this study was to investigate the utility of salivary microRNAs for differentiating children with ASD from peers with typical development (TD), and non-autism developmental delay (DD). The secondary purpose was to explore microRNA patterns among ASD phenotypes. Method: This multi-center, prospective, case-control study enrolled 443 children (2–6 years). ASD diagnoses were based on DSM-5 criteria. Children with ASD or DD were assessed with the Autism Diagnostic Observation Schedule II and Vineland Adaptive Behavior Scales-II. MicroRNAs were measured with high throughput sequencing. Differential expression of microRNAs was compared among ASD (n=187), TD (n=125), and DD (n=69) children in the training set (N=381). Multivariate logistic regression was defined a panel of microRNAs that differentiated ASD and non-ASD children. The algorithm was tested in a prospectively collected, Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. Preliminary results from this study were presented as an abstract at the SUNY Upstate Medical University Autism Research Symposium, April 19, 2017, Syracuse, NY; the Pediatric Academic Societies Meeting, May 6, 2017, San Francisco, CA; and at the North American Saliva Symposium, September 15, 2017, Portland, OR. A sister-study involving characterization of the oral microbiome in this cohort (though alignment of these high throughput RNA sequencing results to the microbial database in K-SLAM) was recently submitted to Autism Research for publication. Dr. Middleton and Dongliang Wang, PhD, served as the statistical experts for this research. naïve set of 62 samples (ASD=37; TD=8; DD=17). Relationships between microRNA levels and ASD phenotypes were explored. Result: Fourteen microRNAs displayed differential expression (FDR<0.05) between ASD, TD, and DD groups. A panel of four microRNAs (controlling for medical/demographic covariates) best differentiated children with ASD from children without ASD in training (AUC=0.725) and validation (AUC=0.694) sets. Eight microRNAs were associated (R> [0.25], FDR<0.05) with social affect, and 10 microRNAs were associated with restricted/repetitive behavior. Conclusion: Salivary microRNAs are “altered” in children with ASD, and associated with levels of ASD behaviors. Salivary microRNA collection is non-invasive, identifying ASD-status with moderate accuracy. A multi-”omic” approach employing additional RNA families may improve accuracy, leading to clinical application.


INTRODUCTION
Autism spectrum disorder (ASD) represents a continuum of deficits in communication and social interaction, as well as restrictive, repetitive interests and behaviors. Healthcare providers have an opportunity to improve outcomes for children with ASD through early diagnosis and referral for evidence-based behavioral therapy. [1][2] Studies suggest earlier treatment contributes to improved social and behavioral outcomes.
An important barrier in the evaluation and treatment of ASD is the lack of objective assessment tools. [3][4][5] Recognition of ASD symptoms generally occurs no earlier than 18-24 months, when deficits in communication emerge. 6 Screening at this stage typically relies on the Modified Checklist for Autism in Toddlers -Revised (MCHAT-R). This parental survey is less than 50% specific. 7 In 2017 the U.S. Preventive Services Task Force determined that insufficient evidence existed to recommend ASD screening. 8 Nonetheless, the American Academy of Pediatrics continues to advocate for universal ASD screening, and pediatricians, faced with no alternative, continue to use subjective, non-specific tools. Clearly, a more accurate and objective toolset would improve ASD evaluation and therapy. extracellular biofluids such as saliva, 11 and renders them attractive biomarker candidates for disorders of the central nervous system (CNS). 12 Studies of miRNA in children with ASD have demonstrated differential expression patterns in post-mortem brain tissue, [13][14] serum, and cultured peripheral lymphoblasts. [15][16] Several of the miRNAs identified in these studies target genes known to be involved in ASD pathogenesis. 17 Brain biopsy is clearly too invasive to be suitable for ASD screening and the physiologic relevance of miRNA expression in cultured lymphoblasts introduces methodological concerns. Given the robust cranial nerve innervation of the oropharynx, its proximity to glymphatic structures, and the sensorimotor pathology observed in children with ASD (food texture sensitivity, 18 taste aversions, and speech apraxia 19 ) we previously explored the potential of salivary miRNA to differentiate children with ASD from typically developing peers. 20 A pilot study of 24 children with ASD demonstrated that salivary miRNAs are altered in ASD and broadly correlate with miRNAs reported to be altered in the brain of children with ASD.
Together, these studies support the potential utility of miRNA measurement in ASD screening. However, the clinical applicability of miRNA studies in persons with ASD has been limited by several factors: 1) No miRNA study has employed more than 55 participants with ASD, 21 despite the broad, heterogeneous nature of the disorder; 2) No miRNA study has enrolled children at the ages (2-6 years) when ASD diagnosis first occurs (i.e. when a diagnostic biomarker panel would have the most clinical utility); 3) No miRNA study has compared children with ASD to peers with non-autism developmental delay (DD) -a comparison required to develop a robust diagnostic toolset; and 4) No study has examined the ability of miRNA signatures to differentiate ASD phenotypes -a priority for the autism community.
The present study sought to address these deficiencies in the literature and establish the diagnostic utility of salivary miRNAs in ASD. We hypothesized that characterization of salivary miRNA concentrations in children with ASD, DD, and typical development (TD) would identify a panel of miRNAs with diagnostic potential. We posited that these miRNAs would exhibit brain-related targets on functional pathway analyses and display associations with specific autism phenotypes (assessed through standard measures of communication, socialization, and repetitive behavior). chosen to include children at the earliest ages of ASD diagnosis, when screening and diagnostic biomarkers would be of most clinical benefit. Recruitment occurred at academic, outpatient, primary and tertiary care clinics between October, 2015 and April, 2018. In the training set (used for miRNA exploration and creation of the regression algorithm) there were 187 children with ASD, 125 children with TD, and 69 children with DD. In the prospective test set (used for validation of the regression algorithm) there were 37 children with ASD, 8 children with TD, and 17 children with DD. Nearly equal numbers of ASD, TD, and DD participants were recruited from each site. An a-priori analysis using Power Analysis and Sample Size Software (v15; NCSS, LLC; Kaysville, Utah, USA), and setting the null area under the curve (AUC) to 0.7, determined that the sample size used in the training set provided 85% power to detect an area AUC = 0.77 (based on a one-sided z test, with an alpha = 0.05), and 99% power to detect an AUC>0.8. Similarly, the replication cohort (N=62) had 85.6% power to detect an AUC=0.78 comparing ASD to non-ASD children. ASD status was defined by DSM-5 diagnosis, confirmed by physician assessment within the previous 12 months, and supported by evaluation with the Autism Diagnostic Observation Schedule (ADOS)-II (or other standardized assessment tool such as the Checklist for Autism Spectrum Disorder, the Autism Diagnostic Interview -Revised, or the Childhood Autism Rating Scale). TD status was defined by history of negative ASD screening on the MCHAT-R and documentation of typical development at a pediatric well child visit within the previous 12 months. DD status was defined by a clinical deficit in gross motor, fine-motor, expressive communication, receptive communication, or socialization that was identified by standardized screening (Survey of Wellbeing in Young Children, MCHAT-R, or Parents Evaluation of Developmental Status) at a regularly scheduled visit, but not meeting DSM-5 criteria for ASD. Targeted recruitment was used to match age and sex across ASD, DD, and TD groups. Exclusion criteria for all groups included feeding-tube dependence, active periodontal disease, upper respiratory infection, fever, confounding neurological (i.e. cerebral palsy, epilepsy) or sensory (i.e. blindness, deafness) impairment, and wards of the state. TD participants with a medical condition requiring daily medication or pediatric specialist care were also excluded.

Participant Characterization
For all participants, extensive medical and demographic characterization was performed, including: age, sex, ethnicity, birth age, birth weight, perinatal complications, current weight, body mass index, oropharyngeal status (e.g. allergic rhinitis), dietary restrictions, medications, chronic medical issues, immunization status, medical allergies, early intervention services, surgical history, and family psychiatric history. Given the prevalence of attention deficit hyperactivity disorder (ADHD) 22 and gastrointestinal (GI) disturbance 23 among children with ASD, survey questions were included to identify these two common medical co-morbidities. GI disturbance was defined by presence of constipation, diarrhea, abdominal pain, or reflux on parental report, International Statistical Classification of Diseases and Related Health Problems (ICD-10) chart review, or use of stool softeners/ laxatives in the child's medication list. ADHD was defined by parental report, or ICD-10 chart review. Adaptive skills in communication, socialization, and daily living activities were measured in all participants using the Vineland Adaptive Behavior Scale (VABS)-II and standardized scores were reported. Evaluation of ASD symptomology (ADOS-II) was completed when possible for ASD and DD participants (n= 164). Social affect (SA), restricted repetitive behavior (RRB) and total ADOS-II scores were recorded.

Saliva Collection and RNA Processing
Saliva was collected from all children in a non-fasting state using a P-157 Nucleic Acid Stabilizing Swab (DNA Genotek; Ottawa, ON, Canada). Saliva was obtained from the sublingual and parotid regions of the oral cavity over a 5-10 second period, taking care to avoid the teeth when possible (https://www.youtube.com/watch? v=AzCpHWqhRQs&feature=youtu.be). Time of saliva collection was recorded, and swabs were kept at room temperature in stabilization solution for up to four weeks prior to storage at −20 °C. Salivary miRNA was purified using a standard Trizol method, followed by a second purification with an RNeasy mini column (Qiagen, Germantown, Maryland). The yield and quality of the RNA samples was assessed using the Agilent Bioanalyzer prior to library construction. RNA was sequenced at the SUNY Molecular Analysis Core at Upstate Medical University with an Illumina TruSeq Small RNA Sample Prep protocol (Illumina; San Diego, California). The targeted read depth for each sample was ten million reads, using 50 base pair single end reads on a NextSeq500 instrument (Illumina, San Diego, CA). Reads for each sample were aligned to the hg38 build of the human genome in Partek Flow (Partek; St. Louis, Missouri) with the SHRiMP2 aligner. Total miRNA counts within each sample were quantified with miRBase precursor and mature-microRNA v21. Poor quality reads (mean q score <30) were eliminated, and samples with total mature miRNA read counts less than 20,000 were excluded. Of the 2813 mature miRNAs aligned, we interrogated 527 miRNAs for differential expression among groups. The 527 miRNAs included: 1) those with robust expression (raw read counts greater than 10 in at least 10 % of samples; 375 miRNAs); and 2) those identified in previous ASD studies 17 and detectable in saliva (raw counts greater than 1 in 10% of samples; 152 miRNAs). Prior to statistical analysis, read counts were quantile-normalized, mean-centered, and divided by the standard deviation of each variable.

Statistical analyses
The primary outcome of this study was the identification of miRNAs that could differentiate children with ASD from children without ASD (including both TD and DD participants) on logistic regression analysis. Differences in medical and demographic characteristics between groups were compared using a two-tailed Student's t-test. In the training set (N=381) a nonparametric Kruskal-Wallis test and a partial least squared discriminant analysis (PLS-DA) were used to identify individual miRNA candidates for differentiating ASD children from TD and DD peers. The miRNAs with significant differences between groups (false discovery rate (FDR) <0.05), and/or PLS-DA weighted sum of absolute regression coefficients ≥ 2.0, were selected for biomarker testing. To control for confounding, medical and demographic characteristics were included in the logistic regression analysis as co-variates. In addition, we also explored the potential influence of RNA quality on any significant miRNA variables using analysis of covariance (ANCOVA) with Diagnosis and RNA Integrity Number (RIN) and their interaction used as main and interaction effects, respectively. Biomarker exploration was performed with Metaboanalyst R package (McGill University, Montreal, Canada, http://www.metaboanalyst.ca/faces/ModuleView.xhtml) using the biomarker workflow. 24 The training set was used to determine threshold (cut-off) concentrations for miRNAs, which were employed in ratios with selected medical/demographic covariates. To avoid "over-fitting" the model and to ensure that the miRNAs accurately differentiated participants with ASD, the algorithm was tested in a naïve replication set of 62 children. Performance was evaluated using AUC analysis from receiver operating characteristic (ROC) curves generated in the training and test sets.
Associations between salivary miRNA concentrations and ASD phenotypic characteristics were explored with Spearman's rank correlations (for dichotomous variables) or Pearson's correlations (for continuous variables), with FDR correction (FDR<0.05). The phenotypic characteristics of interest included: 1) adaptive behavior scores (VABS-II); 2) ASD traits (ADOS-II scores); and 3) medical co-morbidities (presence/absence of GI disturbance or ADHD). Relationships between salivary miRNA concentrations and confounding medical/ demographic characteristics (i.e. age, sex, ethnicity, body mass index, asthma, allergic rhinitis, time of collection, time of last meal, dietary restrictions) were also evaluated with Pearson's or Spearman's rank correlations. Any miRNA-variable association in which R>[0.25] and FDR<0.05 was reported as significant. Secondary analyses investigated the mRNA targets for two sets of miRNAs: 1) the miRNAs "altered" between ASD, TD, and DD groups based on initial Kruskal-Wallis testing; and 2) the miRNAs associated with ASD features on ADOS testing. For the latter, we also used multivariate regression to adjust the correlations by the RIN value and RNA sequencing quality (Q) scores. Functional analysis was performed for each miRNA set in DIANA mirPath v3 online software (http:// snf-515788.vm.okeanos.grnet.gr/). 25 The microT-CDS algorithm was employed to identify species-specific gene targets for each miRNA. DIANA mirPath identified KEGG pathways with significant (FDR < 0.05) target enrichment using a Fisher's Exact Test. A list of high confidence mRNA targets (experimentally validated miRNA-mRNA interaction with microT-CDS score ≥ 0.975) was interrogated for protein-protein interaction networks using moderate stringency settings (interaction score > 0.40) in String v10 software (http://stringdb.org). 26 Enrichment of mRNA target lists for the 961 autism-associated genes on the SFARI autism database (https://gene.sfari.org/database/human-gene/) 27 was explored using a chi-square test with Yates' correction. The number of over-lapping mRNAs was reported, along with enrichment relative to a random sampling of the ~20,000 coding mRNAs.
We also examined the potential differences in RNA quality metrics between the sample groups. Both ASD and Non-ASD groups had mean RIN values of approximately 4.4 in our samples with no significant difference between the ASD and Non-ASD groups (unpaired ttest p value = 0.7465) or between the 3 subgroups (ANOVA F = 0058, p = 0.943). This was also consistent with a lack of difference in the RNA sequencing quality Q scores between the ASD and Non-ASD groups (t-test p = 0.0611) or between all 3 groups (ANOVA F = 1.75, p = 0.173).

Expression of salivary miRNA
Concentrations of 527 mature miRNAs were explored in the saliva of children with ASD, TD, and DD in the training set. Among the 527 miRNAs, 80 were present in the saliva of every participant. The miRNA with the highest salivary concentrations across all participants was miR-203a-3p, accounting for 1.14×10 6 of the total 8.44×10 7 raw read counts in the experiment (1.4%). Kruskal-Wallis non-parametric testing identified 14 miRNAs with significant (FDR < 0.05) differences across ASD, TD, and DD groups ( Figure 1). The miRNA with the largest change was miR-28-3p (χ2=34.2, FDR=1.62E-5), which demonstrated down-regulation in children with ASD relative to both TD and DD groups. Four other miRNAs demonstrated relative down-regulation in the ASD group compared with both TD and DD groups (miR-148a-5p, miR-151a-3p, miR-125b-2-3p, and miR-7706). There were four miRNAs with relative up-regulation in the ASD group compared with TD and DD groups (miR-665, miR-4705, miR-620, and miR-1277-5p). One of the 14 miRNAs (miR-151a-3p) had been identified as "altered" in previous studies of miRNA expression in persons with ASD (21). The remaining 6/14 miRNAs identified on Kruskal-Wallis testing displayed intermediate concentrations in the ASD group (relative to TD and DD groups), or had nearly over-lapping expression patterns with either TD or DD groups.

Classification accuracy
A logistic regression analysis with 100-fold cross validation procedure was used to define a miRNA-based algorithm that differentiated the ASD group from the non-ASD group in the training set (N=381). Only the 28 miRNAs identified on PLS-DA/Kruskal-Wallis analyses were interrogated, and medical/demographic variables were included as covariates. An algorithm employing four miRNAs (miR-28-3p, miR-151-a-3p, miR-148a-5p, miR-125b-2-3p), while controlling for sex, family ASD history, disordered sleep, GI disturbance and presence/absence of chronic medical conditions, correctly identified 125/187 children with ASD and 129/194 children without ASD (Figure 3). This represented an AUC of 0.725 (95% CI: 0.650-0.785). Notably, the four miRNAs included in this algorithm were identified by both PLS-DA and Kruskal-Wallis analyses. Accuracy of the algorithm was prospectively assessed in the naïve test set (N=62). The same algorithm identified 33/37 ASD children and 8/25 non-ASD children in the test set, an AUC of 0.694. This represents a sensitivity of 89.2% and a specificity of 32.0%. Among non-ASD children in the test set, the algorithm was more accurate at differentiating those with TD (4/8) than those with DD (4/17).

Expression of salivary miRNA across ASD phenotypes
Salivary miRNA expression patterns were explored across ASD phenotypes for children with ASD in the training set (n=187; Table S1, available online). Significant correlations (R> [0.25], FDR<0.05) were identified (Table 2) between salivary miRNA levels and presence of GI disturbance (2 miRNAs), but not ADHD. Among all salivary miRNAs, five miRNAs correlated with standardized score on the socialization component of VABS-II testing, two of which (miR-379-5p and miR-221-3p) 14,31-32 had been previously identified in ASD studies. There were no miRNAs correlated with communication or activities of daily living scores on VABS-II testing. Eight miRNAs were correlated with social affect on the ADOS-II. Six of these miRNAs were previously identified in ASD studies (miR-223-3p, miR-142-3p, miR-182-5p, miR-142-5p, miR-181c-5p, miR-148b-3p), 17 and one displayed between-groups differences in the present study (miR-125b-2-3p) and was used in the logistic regression algorithm. Adjustment of these correlations based on either RIN scores or RNA sequencing Q scores did not change them substantially and all remained highly significant (not shown). Ten miRNAs correlated with restricted/repetitive behavior on the ADOS-II, and four of these had been identified in previous ASD studies (miR-136-3p, miR-106a-5p, miR-130a-3p, and miR-431-5p). 17 Notably, all 10 were positively correlated with restricted/repetitive behavior score. Finally, six miRNAs were correlated with total score on the ADOS-II, and all six had been identified in previous ASD miRNA studies. 17 As before, adjustment of these correlations based on RIN or Q scores did not change them substantially. All remained highly-significant (not shown). One of these miRNAs (miR-151a-3p) was down-regulated in children with ASD relative to children with TD and DD, and the miRNA was employed in the logistic regression algorithm.

Influences of clinical characteristics on miRNA expression
Associations of salivary miRNA expression and clinical/demographic characteristics were assessed in the training set (N=381) with Pearson's (continuous) or Spearman's Rank (dichotomous) correlation testing (Table S2, available online). There were no significant associations (R<[0.25], FDR<0.05) between expression of the 527 miRNAs and participant sex, ethnicity, body mass index, dietary restrictions, asthma status, or allergic rhinitis status. Time of saliva collection had the largest number of miRNA associations, compared with other medical/demographic variables tested (n=21). The strongest association was between miR-210-3p levels and time of saliva collection (R=−0.35; t-stat=−6.6; FDR=4.2E-8). One miRNA (miR-23b-3p) was associated with time since last meal (R=0.25; t-stat=4.2; FDR=0.012). Of the 22 miRNAs associated with time of collection or time since last meal, twelve had been identified as potential biomarkers in previous miRNA studies. 17 One was "altered" in the saliva of children with ASD in the current study (miR-151a-3p, R= −0.17, FDR=0.011). Given the importance of age in developing biomarker toolsets, it is worth noting that participant age was weakly (R<[0.25]), yet significantly (FDR<0.05) associated with 34 miRNAs. None of these miRNAs were utilized in the current biomarker panels, but 15 had been identified as potential targets in previous ASD miRNA studies. 17
Analysis of the eight miRNAs associated with ADOS-II total/socialization scores also revealed brain-related mRNA target pathways. The eight miRNAs had a total of 4147 mRNA targets, 3311 of which were unique (Table S3B, available online). There were two miRNAs (miR-182-5p, miR-142-5p) that accounted for 2064 (49.8%) of the total mRNA targets. The mRNA targets over-represented 47 KEGG pathways (Table S4B,

DISCUSSION
This prospective, case-control study of 443 children (ages 2-6 years) identified 28 salivary miRNAs with varying levels among children with ASD, TD, or DD. A panel employing four miRNAs distinguished ASD status in both training and naïve test sets. A subset of salivary miRNAs was associated with measures of adaptive and ASD behaviors. Together, these groups of miRNAs targeted genes strongly related to neurodevelopment and implicated in ASD pathogenesis (Table S5, available online).
There are a number of potential environmental factors which may disrupt levels of miRNAs in the oropharynx of ASD children. Certainly, dietary restrictions in children with ASD 18 may alter the salivary miRNA milieu. However, the current study found no associations between saliva miRNA levels and the presence of dietary restrictions, and only two miRNAs were strongly associated with GI disturbance. In addition, there was no difference in the rate of dietary restrictions between ASD, DD, and TD groups. A second potential mechanism for salivary miRNA disruption could be differences in dental hygiene, given the resistance of many children with ASD to teeth brushing. 33 For this reason, this study specifically excluded children with active dental infections or decay. There are alterations in the oral microbiome of children with ASD 34 which may drive a portion of salivary miRNA changes, but oral microbiome differences in children with ASD are largely unrelated to the bacteria implicated in dental carries. 35 Children with ASD experience difficulties with oral-motor (speech apraxia) and oral-sensory (food texture sensitivity) processing. 19,36 The cranial nerves which guide these processes may contribute to salivary miRNA patterns. Brain-relatedness of the salivary miRNAs identified in this study is supported by the functions of their mRNA targets, which include axonal guidance, neurotrophic signaling, GABAergic synapse, and addiction pathways (Tables S4A and S4B). For example, miR-148a-5p (utilized in the diagnostic panel of the current study) targets seven mRNAs involved in axon guidance (Table S3A, available online), and two of these (SLIT3 and SRGAP3) are autism candidate genes. 27 The SLIT3 protein product acts as a molecular guidance cue in axonal outgrowth by interacting with the protein product of another autism candidate gene, ROBO1. 37 Notably ROBO1 is a target of miR-944 (Table S3A, available online), a miRNA associated with ASD-status in the present study, and highly correlated with miR-148a-5p in both concentration ( Figure 1) and function ( Figure S1, available online). The parallel functions of miR-944 and miR148a-5p in axon guidance, coupled with their overlapping expression in children with ASD, highlight their potential significance in ASD pathophysiology.
The glymphatic system represents yet another potential route for salivary entry of brainrelated miRNAs. The anatomical proximity of the perivascular drainage spaces in the glymphatic system to the oropharynx creates a prospective avenue for gut-brain cross-talk and miRNA transfer (11). In light of the pronounced diurnal activity displayed by the glymphatic system, 38 indirect support for this transfer may lie in the surprising correlations between salivary miRNA levels and time of collection (Table S2, available online). In addition, the mRNA targets of ASD-associated miRNAs show enrichment for circadianrelated pathways (Table S4A, available online), which is notable since disordered sleep is a common medical condition among children with ASD. 39 The potential relevance of salivary miRNA levels to ASD behavior is underscored by the large number of salivary miRNAs associated with measures of ASD symptoms on the ADOS-II (Table S1, available online). Previous studies have described miRNAs as "altered" in persons with ASD relative to healthy control participants. 17 The increased power of the current investigation provides an opportunity to explore miRNA patterns among ASD phenotypes. Here we identify eight miRNAs associated with social affect and ten miRNAs associated with restricted/repetitive behavior. Such associations may be driven by robust miRNA "alterations" in a subset of children with a similar single nucleotide polymorphism or copy number variant. 40 In these children, phenotypic similarities may result from genetic mutations which produce a direct miRNA change, or lead to compensatory miRNA responses. One example is miR-106a-5p. 41 This miRNA has been previously identified in three separate ASD studies of post-mortem brain, 30 blood, 21 and lymphoblasts. 15 It targets 20 mRNAs involved in axon guidance ( Figure S1, available online), 25 including 4 autism candidate genes (SEMA5A, NTNG1, SRGAP3, and MAPK1). 27 Here, we show that miR-106a-5p levels are directly associated with restricted/repetitive behavior in children with ASD (Table 2). Thus, altered levels of miR-106a-5p could target key transcripts involved in brain development which underlie restricted/repetitive behaviors. Additional studies tracking expression patterns of such miRNAs alongside behavioral therapy interventions is warranted before strong conclusions can be drawn.
This study defines an algorithm using four miRNAs to differentiate children with ASD from peers with TD or DD (Figure 3). In a naïve test set the panel demonstrated 89% sensitivity and 32% specificity. This accuracy approaches that of subjective measures currently employed (e.g. MCHAT-R (7)) with the added benefit of being fast, objective, and noninvasive. Emerging biomarker work in eye-tracking, 3,42 imaging, 43 genetic, 44 and electrophysiologic markers 45 has also shown considerable promise for identifying ASD status. The future of ASD evaluation will likely involve a multi-factorial approach employing each of these components in concert. The results of this study suggest that salivary RNA biomarkers deserve strong consideration within this field. Indeed, bolstering the current algorithm with a poly-"omic" analysis of additional RNA families has led to an even more comprehensive and accurate approach. 46 Among the four miRNAs utilized in the diagnostic algorithm, two (miR-125b-2-3p and miR-151a-3p) were strongly associated with ASD traits on ADOS evaluation ( Table 2) and one (miR-151a-3p) was identified in previous studies. 17 Limited overlap with previous miRNA studies may have resulted because blood and lymphoblast miRNAs are not reliably transferred to (or expressed in) saliva. This finding may also reflect limited generalizability of small cohort studies to a large heterogeneous population of children with ASD. Levels of certain miRNAs can vary widely from child to child depending on many factors (e.g. time of collection, co-morbid medical conditions, age, sex). For this reason, "outlying" miRNA concentrations in just a few individuals could lead to the assumption that between-group differences exist, when the mean group expression is effectively biased by just a few samples. Small studies (i.e. nearly all the previous studies of miRNAs in persons with ASD) are particularly prone to this. Here, we use a large sample size as well as complementary Kruskal-Wallis and PLSDA approaches to select miRNAs which avoid this pitfall.
It is also notable that many previously identified miRNA biomarkers (11 miRNAs) demonstrated associations with time of collection (Table S2, available online). This factor that has not been routinely considered in previous ASD miRNA studies. Given recent findings that a significant proportion of serum-based miRNAs demonstrate diurnal variation, 47 these findings likely apply to blood-based biomarkers as well. Further studies examining the interaction between miRNA expression and circadian rhythm could be important in understanding the role of these molecules in sleep-wake cycles, and provide valuable information in the development of miRNA biomarkers for clinical application. Importantly, there were no differences in collection time between the ASD, TD, and DD groups in this study.
Surprisingly there was little overlap between the salivary miRNAs identified in our pilot investigation, and those identified in the present study. 20 This may have resulted from three important differences in study protocols; 1) the pilot study used expectorated saliva, while this investigation collected saliva with a swab technique. This change was made because children with ASD have difficulty producing expectorant on command. It may have led to differences in ratios of cell-derived and (vesicle) carrier-derived miRNA. 2) The pilot study involved children 5-14 years of age, while the current study enrolled children 2-6 years of age. This change was made to capture children at the age when ASD diagnosis is first made and screening/diagnostic testing is most needed (1). It may have influenced a subset of miRNAs with age-related expression.
3) The pilot study targeted children with "high functioning" ASD (average ADOS-II score = 10.6 ± 4.1), while this large follow-up study included all children with ASD regardless of severity (average ADOS-II score = 16 ± 6). Given that salivary miRNA expression is associated with levels of ASD symptoms (measured by ADOS-II) it is likely that expanding the current study to include a heterogeneous population of children with ASD led to changes in observable betweengroups differences.
There are numerous medical and demographic factors that must be considered when identifying and testing physiologic biomarkers. The prospective nature of the current study allows us to control for many of these factors by employing identical collection, storage, and sample processing techniques across groups. We have also attempted to match groups based on relevant factors such as age, gender, ethnicity, body mass index, and time of collection. Unfortunately, complete matching of all factors is nearly impossible. As a result, the training set displays between-group differences in age and sex. However, it is worth noting that the age range utilized in the present study (2-6 years) is extremely tight compared with many biomarker studies and the resulting age difference between ASD and TD groups (7 months) is unlikely to have significant bearing on miRNA expression. In addition, none of the miRNA biomarkers identified in this study demonstrated significant correlations with age or sex, and the multivariate regression algorithm controls for sex.
Another extremely important topic to consider when assessing the veracity of RNA research is nucleic acid integrity, and its potential influence on biomarker outcomes. Although we report RIN across the three groups of samples, it is important to note that this metric likely underestimates RNA quality in miRNA-enriched samples. Unlike longer messenger RNAs, small RNAs (e.g. miRNA, piwi-interacting RNA, small nucleolar RNA) are relatively resistant to salivary endonucleases. As a result, even samples with low RIN values (and presumably poor RNA quality) can demonstrate excellent miRNA yields on bio-analyzer output ( Figure S2, available online). Indeed, a study using human cell and tissue samples subjected to total RNA purification following longitudinal heat degradation has demonstrated that RIN values rapidly decrease with heat exposure and housekeeping messenger RNAs are lost to detection, while miRNAs remain remarkably stable over time. 48 In spite of the limits associated with RIN reporting, we note that: 1) average RIN for this data set exceeded RIN values reported in previous saliva RNA studies; 49 2) there was no difference in average RIN among ASD, TD, and DD groups; and 3) RNA-ADOS correlations were actually strengthened when RIN was added as a covariate. We encourage any future studies utilizing saliva RNA measures to employ stringent methods for RNA stabilization and extraction, and to carefully assess the influence of RNA integrity on biomarker findings.
This study provides large-scale evidence that salivary miRNA may be used to differentiate children with ASD from peers with typical development, or non-autism developmental delay. It shows that levels of salivary miRNAs are correlated with measures of adaptive and ASD behaviors and that these miRNAs target pathways implicated in ASD-pathogenesis. Improving specificity of the defined salivary miRNA algorithm is crucial for clinical utility. This has been achieved through a multi-modal approach, employing additional "-omic" measures. 46 Additional characterization of the factors that influence salivary miRNA expression will also be crucial.

Supplementary Material
Refer to Web version on PubMed Central for supplementary material.

Figure 1. Salivary microRNAs (miRNAs) are Differentially Expressed Across Groups
Note: The 14 miRNAs with differential expression (false discovery rate [FDR]<0.05) across autism spectrum disorder (red; n=187), developmental delay (DD) (green; n=69), and typically developing (TD) (blue; n=125) groups on Kruskal-Wallis testing are shown, along with χ2 statistics. Colored boxes represent relative group expression (measured by Pearson's distance metric) and miRNAs are clustered in the heatmap using a complete clustering algorithm. Hicks    b miRNAs with between-groups differences in the current study.