Your new experience awaits. Try the new design now and help us make it even better

ORIGINAL RESEARCH article

Front. Mol. Biosci., 24 November 2025

Sec. Molecular Diagnostics and Therapeutics

Volume 12 - 2025 | https://doi.org/10.3389/fmolb.2025.1696704

This article is part of the Research TopicUnraveling Molecular Mechanisms of Neurodegenerative Disease Development: Insights from Stem Cell-Derived 2D and 3D in vitro ModelsView all 7 articles

Decoding the peripheral transcriptomic and meta-genomic response to music in autism spectrum disorder via saliva-based RNA sequencing

  • 1Unidade de Xenética, Instituto de Ciencias Forenses, Facultade de Medicina, Universidade de Santiago de Compostela, and Genética de Poblaciones en Biomedicina (GenPoB) Research Group, Instituto de Investigación Sanitaria (IDIS), Hospital Clínico Universitario de Santiago (SERGAS), Santiago de Compostela, Spain
  • 2Genetics, Vaccines and Infections Research Group (GenViP), Instituto de Investigación Sanitaria de Santiago, Universidade de Santiago de Compostela, Santiago de Compostela, Spain
  • 3Centro de Investigación Biomédica en Red de Enfermedades Respiratorias (CIBER-ES), Madrid, Spain
  • 4Translational Pediatrics and Infectious Diseases, Department of Pediatrics, Hospital Clínico Universitario de Santiago de Compostela, Santiago de Compostela, Spain

Introduction: Behavioral interventions for autism spectrum disorder show variable outcomes, highlighting the need for complementary therapies. Music-based interventions are promising, yet their molecular mechanisms remain unclear. Saliva-based RNA sequencing (RNA-seq) provides a non-invasive framework to monitor neuroimmune and metabolic dynamics, but its application in autism remains underexplored.

Methods: We explored the buccal transcriptional effects of music exposure in five individuals with autism (8–37 years; 60% female). To overcome saliva-specific limitations, we combined Poly-A selection and Human-Enriched protocols preparation methods to enhance human transcript detection and reproducibility while capturing microbial signals.

Results: Individually, each dataset revealed a few differentially expressed genes, but integrated analysis improved biological resolution. Consistently modulated genes included HERC6, TSPAN5, and REM2, involved in neurodevelopmental and immune functions. Enrichment analyses highlighted pathways associated with immune regulation, oxidative phosphorylation, and epithelial differentiation, hallmarks of autism, such as immune dysregulation and mitochondrial dysfunction. Co-expression network analysis identified modules correlated with music exposure. The AKNA module, previously linked to autism, was downregulated and enriched for Ras-related GTPase and immune pathways, suggesting modulation of intracellular signaling and inflammation. Conversely, upregulation of the UBE2D3 module indicated activation of endoplasmic reticulum stress responses, a contributor to autism. Exploratory metagenomics identified 15 microbial species responsive to music exposure, including Acidipropionibacterium acidipropionici and Propionibacterium freudenreichii, producers of propionic acid, a metabolite associated with autism-like behaviors and neuroinflammation.

Conclusion: Saliva-based RNA-seq can stably capture transcriptomic and microbial responses to behavioral stimuli. Music exposure modulates neuroimmune pathways relevant to autism, supporting the biological plausibility of music therapy and demonstrating saliva-based RNA-seq as a viable, non-invasive tool for monitoring intervention outcomes.

Introduction

Autism Spectrum Disorder (ASD) is a neurodevelopmental condition involving deficits in communication, behavior, and motor skills (World Health Organization, 2021; Kamp-Becker, 2024; American Psychiatric Association, 2022). Symptoms typically emerge between 18 and 24 months, and diagnosis is based on persistent social interaction difficulties, communication challenges, and repetitive behaviors (Zeidan et al., 2022). ASD affects about 1.01% of U.S. children and 0.73% in Europe, with a male-to-female ratio of 4.3. This disparity may be partially due to underdiagnosis of females, who often exhibit subtler symptoms and camouflaging behaviors, suggesting that females may need a higher genetic or environmental threshold to express ASD traits (Volkmar et al., 2014). However, when diagnosed, females often present with more severe symptoms and intellectual disability (Delobel-Ayoub et al., 2020). ASD spans a wide spectrum, from severe impairments to mild social differences (Volkmar et al., 1989), and some individuals also show developmental regression, particularly in language and social skills (Ukkola-Vuoti et al., 2013). Educational, occupational, and social difficulties are frequent. However, others can exhibit exceptional abilities, such as enhanced pitch recognition or visual processing (Samson et al., 2012; Muth et al., 2014; Mazurek et al., 2017; Remington and Fairnie, 2017). Comorbid conditions are common, including intellectual disability (∼31%), sleep disturbances, seizures, anxiety, Attention-Deficit/Hyperactivity Disorder (ADHD), Obsessive-Compulsive Disorder (OCD), and mood disorders (Kielinen et al., 2004; Soke et al., 2018). Needs vary; some people live independently, while others require lifelong support (Seltzer et al., 2004).

ASD’s causes are multifactorial, involving genetic and environmental components. Brain abnormalities such as cerebellar malformations and early overgrowth are being studied as potential causal factors (Hazlett et al., 2017; Stoodley et al., 2017). Genetic research implicates genes involved in brain development and synaptic function (Zoghbi, 2003; McDougle et al., 2005; Voineagu et al., 2011). Environmental risk factors include prenatal exposure to valproic acid or thalidomide, older parental age, low birth weight, and cesarean delivery (Bjørk et al., 2018). Maternal inflammation during pregnancy is another risk factor for ASD. Cytokines like IL-6 and IL-17a, elevated during infections, may disrupt fetal brain development and promote ASD traits (Jash and Sharma, 2021). Postnatal immune profiles also differ in ASD individuals, implicating immune dysregulation. By contrast, folic acid supplementation during pregnancy appears to have a protective effect (Bjørk et al., 2018).

Behavioral interventions remain the cornerstone of ASD management and aim to improve communication, social interaction, and adaptive functioning. Common approaches include applied behavior analysis (ABA), social skills training, and speech or occupational therapy. While these interventions are effective for many individuals, their outcomes can vary and often require long-term, intensive engagement (Lord et al., 2020; Sandbank et al., 2020). Consequently, complementary approaches that enhance motivation, emotional engagement, and social responsiveness have gained attention. In this context, music therapy may represent a particularly promising avenue, as it engages widespread brain regions and influences emotion, cognition, and motor skills (Zatorre and McGill, 2005). It stimulates dopamine release, affects autonomic responses, and promotes neuroplasticity. Its effects vary by age, sex, culture, and musical training (Zatorre and McGill, 2005; Kanduri et al., 2015). Music therapy is used in diverse clinical contexts such as addiction, Post-Traumatic Stress Disorder (PTSD), dementia, and brain injury (The American Music Therapy Association, 2014). In ASD, music is promising due to the usually higher auditory sensitivity and pitch memory of the patients (Remington and Fairnie, 2017). While auditory sensitivity can lead to overload in ASD patients, it also allows deep engagement with music’s predictability and structure (Goris et al., 2020). Active and passive music therapy can support not only emotional regulation and social communication but also improve emotional processing, motivation, attention, and interaction (Cook et al., 2019; Navarro et al., 2025a; Navarro et al., 2025b; Salas et al., 2025). Structured programs, in turn, strengthen parent-child relationships while reducing symptom severity (Lense et al., 2020). Music stimulation has been shown to increase oxytocin and vasopressin, neuropeptides involved in empathy and bonding. Oxytocin enhances social cue recognition and has been genetically linked to ASD (Domes et al., 2007). Compared to speech therapy or medication, music therapy often yields greater improvements in life quality (Israel et al., 2008). Neuroimaging studies have shown that music therapy enhances brain connectivity, especially between auditory and motor regions (Sharda et al., 2018). However, most research studies focused on children and lacked standardization. Given individual variability in music sensitivity, personalized approaches are essential. In this context, the emerging field of sensogenomics investigates how sensory stimuli like music influence gene expression and brain function, offering new insights into ASD.

Many children with ASD exhibit significant fear and avoidance of medical procedures, particularly those involving needles (Sukhodolsky et al., 2008; Wolff and Symons, 2013) and, therefore, alternative less invasive approaches are under investigation. Though exposure therapy may help reduce stress related to blood draws, its application is limited by safety and logistical constraints (Meindl et al., 2019). In this context, saliva presents a non-invasive, easily collectable alternative for biomarker discovery, supplying DNA, RNA, and proteins for molecular analysis, including identification of polymorphisms linked to ASD. Saliva is increasingly recognized as a valuable non-invasive biofluid for biomedical research and biomarker discovery. In neurological conditions, saliva is being explored for monitoring brain health (Farah et al., 2018; Engeland et al., 2019). Parasympathetic and sympathetic stimulation of salivary glands triggers acetylcholine and noradrenaline release to the buccal cavity, respectively, promoting protein secretion; these neurotransmitters have been detected in rodent salivary glands (Murai et al., 1998; Nakamura et al., 2004; Proctor and Carpenter, 2007). The nervous regulation of the salivary glands supports the concept of a bidirectional oral-brain axis, where local inflammatory stages may influence each other (Sansores-Espana et al., 2021). In fact, different hormonal biomarkers, such as cortisol and oxytocin, have been studying in saliva from ASD patients (Majewska et al., 2014; Bakker-Huvenaars et al., 2020). Studies in different tissues, including saliva, have also shown that oxytocin levels are lower in ASD children (John and Jaeggi, 2021). Several studies have tried to “treat” ASD social impairments with oxytocin some of them showed improvements (Yatawara et al., 2016; Parker et al., 2017), while others showed no differences between the oxytocin treatment and the placebo treatment (Anagnostou et al., 2012; Munesue et al., 2016). In the context of ASD, small RNA molecules such as microRNAs (miRNAs), piRNAs, and snoRNAs are involved in neurodevelopment and synaptic regulation. In a 2018 study, 32 RNA markers in saliva distinguished ASD from controls with 85% accuracy. Changes in miRNAs like miR-146b-5p and piR-6463 were found diagnostically relevant (Hicks et al., 2018). miR-141-3p expression was found to be negatively correlated with oral bacteria of the genus Tannerella in ASD children’s saliva, suggesting that dysregulation of miRNAs and the oral microbiome dysbiosis may be linked to cognitive impairments (Frey et al., 2018; Ragusa et al., 2020; Kamble et al., 2024). Indeed, 70% of ASD children have gastrointestinal symptoms and are often associated with greater symptoms severity (Benach et al., 2012; Lefter et al., 2019). Numerous investigations have reported different microbial profiles in individuals with ASD compared to neurotypical controls, suggesting potential effects on ASD-related pathways (Iglesias-Vazquez et al., 2020; Peralta-Marzal et al., 2024). Although specific microbial differences associated with ASD remain inconsistent across studies, accumulating evidence supports a crucial role of the bidirectional microbiota gut–brain axis in ASD and other neurodevelopmental disorders (Cryan et al., 2019; Peralta-Marzal et al., 2021). Microbial metabolites can cross the intestinal barrier and enter systemic circulation, where they may exert effects on distant organs and tissues, such as the brain (Popoff, 2011). It has been hypothesized that some species of Clostridium can release neurotoxins that could influence the ASD pathology. For instance, Clostridium difficile is an emerging pathogen associated with antibiotic-related diarrhea and increased gut permeability. Individuals with ASD are frequently exposed to intensive antibiotic regimens, heightening their susceptibility to C. difficile infection (Kuhn et al., 2012). Studies in mouse models indicate that C. difficile can alter p-cresol levels and dopamine-β-hydroxylase activity, affecting dopaminergic signaling in the brain, suggesting a potential link between C. difficile and neurodevelopmental alterations in ASD (Kakabadze et al., 2022). However, other studies have reported contrasting findings, showing no significant association between ASD and C. difficile (Khalil et al., 2021).

Despite its advantages, saliva-based diagnostics face technical challenges. Only about 20% of RNA sequences align with the human genome, while 30% match known oral microbes, leaving 50% unclassified, likely representing fungi, viruses, or degraded fragments (Rylander-Rudqvist et al., 2006; Spielmann et al., 2012; Frey et al., 2018; Kamble et al., 2024). Human RNA is easily degraded, and bacteria produce abundant RNA. In datasets like RNAgE-WT, human RNA comprised just 0.1% of reads, while bacterial RNA exceeded 75% (Gosch et al., 2024). Low analyte concentrations, sample variability, and degradation during handling (e.g., freeze-thaw cycles) further complicate analysis (Gröschl et al., 2001; Christodoulides et al., 2005; Pfaffe et al., 2011).

Exploring music’s effect on gene expression in ASD through saliva offers a promising, non-invasive research avenue. However, it involves two major obstacles: the complexity of human vs. microbial RNA content, and the novelty of assessing transcriptional changes induced by music in a population with neurobiological heterogeneity. Addressing these challenges could advance personalized music-based therapies.

This pilot study evaluated the feasibility of using saliva to analyze human RNA changes before and after music exposure in individuals with ASD. We employed two complementary methods of RNA-seq library preparation, and the results revealed reproducible transcriptomic responses, supporting saliva’s potential as a non-invasive tool for exploring the neurobiological impact of music in ASD. These findings lay the groundwork for future research into personalized, music-based interventions for neurodevelopmental disorders.

Materials and methods

Participants and sample collection

This study was conducted as part of the Sensogenomics project (https://sensogenomics.com; Navarro et al., 2021; Gómez-Carballa et al., 2023; Navarro et al., 2023; Gómez-Carballa et al., 2025a; Gómez-Carballa et al., 2025b), which explores the biomolecular and physiological effects of music in diverse population groups. The experimental concert, Sensogenoma22, was performed by Real Filharmonía de Galicia, conducted by Baldur Brönnimann in the Auditorio de Galicia (Santiago de Compostela, Spain). A selected repertoire of 50 min, trying to capture attention and impulse emotions through the variety of timbres, tempos, styles, and tonalities: The Unsewered Question (C. Ives), The Merry Wives of Windsor (O. Nicolai), Slavonic Dances n° 2 and n° 3, op. 46 (A. Dvořák), Oblivion (A. Piazzolla), Hungarian Dance n° 5 (J. Brahms), The Barber of Seville: Overture (G. Rossini), Danzón n° 2 (A. Márquez). To minimize response bias, the musical program was not disclosed in advance. All participants remained seated throughout the performance to eliminate the influence of physical activity. The environment was carefully designed to be comfortable, calm, and non-invasive, thereby reducing stress or other external factors that could affect physiological responses. The audience consisted of individuals from distinct groups, including individuals diagnosed with ASD, healthy control participants, and general spectators. For the present pilot study, a subset of five individuals with ASD was selected. The mean age of the ASD subjects was 23.8 years (8–37 years), and 60% were females. To ensure the reliability and comparability of transcriptomic data, all participants were instructed to abstain from eating, drinking, smoking, chewing gum, or engaging in physical activity for at least 30 min before saliva collection.

Saliva samples were collected immediately before and after the concert, simultaneously for all participants. The collection was performed under the supervision of trained healthcare professionals to guarantee consistency and adherence to the protocol. Each participant used an Oragene RNA saliva collection device (ORE-100, DNA Genotek), which is optimized for the non-invasive collection, stabilization, and transport of high-quality RNA from saliva. Participants were instructed to rub the sponge tips of the device along both sides of their gums, avoiding contact with teeth, until sufficient saliva was absorbed. The sponge tips were then inserted into tubes containing 1 mL of stabilizing solution, which preserves RNA integrity. Samples were subsequently stored at room temperature until processing.

Written informed consent was obtained from the legal guardians or responsible parties of all participating individuals. The study was approved by the Ethics Committee of the Xunta de Galicia (registration code: 2020/021) and conducted in accordance with the principles outlined in the Declaration of Helsinki. Figure 1 shows a scheme of the experimental design used in the present study.

Figure 1
Panel A shows a scatter plot with Log₂ Fold Change values comparing the Human-Enriched and Poly-A datasets. A trend line with confidence intervals is included. Panel B contains three sets of scatter plots and histograms for Poly-A, Human-Enriched, and Merged datasets. Correlation coefficients are displayed for each comparison: 0.860, 0.901, and 0.988, indicating strong correlations.

Figure 1. Graphic scheme showing the design of the present work. The figure was built using Biorender resources (https://biorender.com/).

RNA isolation and RNA-sequencing

RNA was isolated from 500 μL of saliva using the RNeasy Micro Kit (Qiagen), using the protocol provided by the extraction kit as recommended by the Oragene tubes supplier. Samples were subjected to a DNase treatment to remove DNA. Further purification, concentration, and DNAse treatment were performed using the RNA Clean & Concentrator Kit (Zymo). Final RNA concentration was verified with Nanodrop, yielding >50 ng/μL, suitable for downstream applications.

RNA quantity and integrity (RIN) were assessed using the TapeStation 4,200 system (Agilent), and RNA concentration was quantified using the Qubit fluorometer (Thermo Fisher Scientific). Two different Illumina library preparation strategies were employed: i) Illumina RNA Prep with Enrichment (https://emea.illumina.com/products/by-type/sequencing-kits/library-prep-kits/rna-prep-enrichment.html), which uses hybrid capture to enrich for specific RNA targets, including non-coding RNAs, degraded RNA, or low-input samples (hereafter referred to as the Human-Enriched library), and ii) TruSeq Stranded mRNA, which selectively captures polyadenylated mRNA (poly-A selection), primarily targeting coding transcripts and thereby focusing on protein-coding gene expression (hereafter referred to as the Poly-A library). Human-Enriched libraries were prepared using the Illumina RNA Prep with Enrichment (ref. 20040537) workflow, which includes RNA fragmentation, reverse transcription into cDNA, adapter ligation, library amplification, purification, and target enrichment via hybridization capture using the Illumina Exome Panel (ref. 20020183). For the Poly-A library, strand-specific libraries were generated using the TruSeq Stranded mRNA technology (ref. 20020594). Indexing was performed using unique dual indexes with the IDT for Illumina TruSeq RNA UD Indexes kit (ref. 20022371).

Sequencing for both library types was performed on the NovaSeq 6000 System (Illumina) using a paired-end configuration (2 × 100 bp). Each sample generated approximately 30 million paired-end reads when using the Human-Enriched and 60 million paired-end reads with the Poly-A, with an average of 85% of bases achieving a quality score above Q30.

Transcriptomic data processing and statistical analysis

To assess the quality of the raw reads, FastQC (v0.11.9) and MultiQC were employed. Following the initial quality assessment, reads were processed using Trimmomatic (Bolger et al., 2014) to remove low-quality bases and residual adapter sequences. Reads were further filtered using the sliding window approach, where sequences are scanned from the 5′ to 3′ end and trimmed when the average quality within a defined window falls below a threshold. High-quality reads were aligned to the Homo sapiens reference genome (GRCh38/hg38) using STAR (Spliced Transcripts Alignment to a Reference) version 2.7.9a (Dobin et al., 2013). Gene-level quantification was performed using the HTSeq package (version 2.0.3) (Putri et al., 2022). Data normalization and differential expression (DE) analysis was conducted using the DESeq2 package (version 1.46.0) in RStudio (v4.4.2) (Love et al., 2014). To evaluate the robustness and reproducibility of the results, two datasets derived from different RNA isolation methods, Human-Enriched and Poly-A selected libraries, were compared. Exploratory data visualization was performed using histograms of read counts, Principal Component Analysis (PCA), hierarchical clustering heatmaps, and volcano plots. These analyses were executed with R packages ggplot2 (Wickham, 2016), ComplexHeatmap (Gu, 2022), and EnhancedVolcano (Blighe et al., 2020). DESeq2 package was used to detect differentially expressed genes (DEGs) after the musical stimuli when compared to the baseline. We include the patient ID in the model as a random effect to account for patient-to-patient differences in the baseline. Genes common to both datasets and displaying statistically significant differential expression (P-value <0.05) were examined via Pearson correlation to assess consistency across the two preparation methods. Upon confirming a strong correlation, a merged dataset was constructed by summing the read counts of shared genes. This integration strategy aimed to enhance the sensitivity of the analysis and mitigate false negatives associated with the limitations of each isolation protocol.

Functional enrichment analysis was performed on the merged dataset. First, a differential pathways analysis was performed using the Gene Set Variation Analysis (GSVA) R package (Hänzelmann et al., 2013). Gene sets corresponding to Gene Ontology (GO) biological processes were obtained from the Molecular Signatures Database (MSigDB) and used as the reference collection. Differentially regulated pathways (DRPs) were identified using the limma package (Ritchie et al., 2015), applying an adjusted P-value threshold of <0.05, with TP1 as the reference group. Additionally, we conducted an over-representation analysis (ORA) of GO terms with the ClusterProfiler R package (Wu et al., 2021; Yu, 2024), org.Hs.eg.db (Carlson, 2023), and AnnotationDbi (Pagès et al., 2024). These tools facilitated the identification of pathways significantly enriched among DEGs, providing insight into the biological processes modulated in response to the experimental conditions. For specific analyses, we generated a subset of neurobiological processes by systematically and manually curating the GO (biological processes) database, selecting entries that included neurobiologically related terms. This approach resulted in a total of 483 biological processes (Supplementary Table S1).

We next compared the DEGs found in the merged dataset with those indexed in the Simons Foundation Autism Research Initiative (SFARI) database (Abrahams et al., 2013). This database, developed for the ASD research community, compiles information on genetic variants associated with the etiology of ASD. Genes in the SFARI database are classified into four categories: i) syndromic (genes with mutations strongly linked to ASD and additional clinical features beyond core diagnostic criteria), ii) Level 1 (high-confidence ASD risk genes), iii) Level 2 (strong candidate genes), and iv) Level 3 (genes with moderate evidence based on previous studies). This comparison allowed us to assess the extent to which our findings align with the existing literature on ASD genetics.

Co-expression modules analysis

To explore gene co-expression patterns potentially associated with musical stimuli in the buccal mucosa of ASD patients, we constructed a signed weighted correlation network using the Weighted Gene Co-expression Network Analysis (WGCNA) R package (Langfelder and Horvath, 2008). Gene expression data of the merged dataset were first normalized and corrected for inter-individual variability. To focus on the most informative signals, we retained the top 75% of genes with the highest variance across samples. A soft-thresholding power of 26 was selected based on scale-free topology criteria, yielding a model fitting index above 0.80 for both datasets (Supplementary Figure S1A). We computed the Topological Overlap Matrix (TOM) and corresponding dissimilarities (1–TOM), followed by hierarchical clustering. Modules were defined using a minimum size of 30 genes and merged using a dendrogram cut height of 0.2 (Supplementary Figure S1B). Each resulting module was assigned a unique color identifier. Module eigengenes were correlated with the TP1-TP2 to identify modules significantly associated with musical stimulation, based on gene significance (GS). Within each co-expression module identified, Module Membership (MM) for each gene was calculated as the correlation between each gene’s expression profile and the module eigengene (the first principal component summarizing the expression pattern of all genes within that module). MM thus reflects the degree of intramodular connectivity, indicating how strongly a gene is associated with the overall expression trend of its module. Genes with high MM values were considered central, as they likely play key regulatory roles in the biological processes represented by the module. MM was used to identify hub genes (those with the highest connectivity). Each significant module was labeled according to its hub gene name.

To assess the biological relevance of the significant modules, a clustered over-representation was conducted using the compareCluster function from the clusterProfiler R package (Wu et al., 2021), which enabled the parallel evaluation of all genes belonging to each significant module, irrespective of their MM. GO biological process terms were used as a reference database. We reduce redundancy in GO terms with the simplifyEnrichment R package (Gu and Hubschmann, 2023), using the binary cut method for clustering the similarity matrix of biological terms.

Exploratory metagenomic analysis of salivary samples

To investigate potential shifts in bacterial composition associated with musical stimulation, we conducted an exploratory metagenomic analysis using reads that did not align to the human reference genome. These unmapped reads were extracted using the STAR aligner with the --outReadsUnmapped Fastx parameter. Read quality was assessed with FastQC and summarized using MultiQC, revealing consistently high base quality scores across all samples. Common metagenomic artifacts, such as per-base sequence content biases, sequence duplication, and adapter contamination, were detected. These issues were addressed using Trimmomatic, which efficiently removed adapter sequences and improved overall read quality. These unmapped reads were processed and classified against a bacterial genome database using the CLARK metagenomic classifier (v1.3.0.0) (Ounit et al., 2015; Ounit and Lonardi, 2016), which provides a taxonomic profile based on discriminative k-mers (short DNA sequences of fixed length k) and assigns bacterial labels based on discriminative k-mer comparisons against the NCBI/RefSeq bacterial genome database. CLARK constructs an index from reference genomes by identifying unique k-mers that are specific to each genome. Shared k-mers among multiple genomes are excluded to ensure that only genome-specific sequences are used for classification. During the read assignment, CLARK quantifies the number of shared discriminative k-mers between each read and reference genomes, classifying each read to the genome with the highest match. The analysis was conducted in full mode, which utilizes the complete set of discriminative k-mers and outputs detailed metrics, including read counts per taxon and confidence scores (ranging from 0.5 to 1.0). This mode was selected to retain maximum resolution and accuracy in taxonomic assignments.

DESeq2 package was used to normalize the data and detect differences in the abundance of microorganisms after the musical stimuli when compared to the baseline. We include the patient ID in the model as a random effect to account for patient-to-patient differences in the baseline using the limma package.

Graphical visualizations were generated in RStudio (v4.4.2) using the ggplot2 and ggwaffle packages (Rudis and Gandy, 2015). Visual outputs included bar plots of the most abundant bacterial genera and comparative plots of relative and absolute abundances before and after musical exposure.

Results

Quality assessment of raw data and comparative evaluation of the mapping from library preparation protocols

Initial quality control assessments revealed that the raw sequencing reads were of very high quality across all samples, with optimal values observed for most quality parameters (data not shown). However, two quality metrics were consistently flagged across all samples: “Per Base Sequence Content” and “Sequence Duplication Levels.” These warnings, though persistent, are not unusual in RNA-seq experiments and are not indicative of poor sample quality. In contrast, the presence of adapter contamination indicated a clear need for improvement. To address this, sequencing adapters and low-quality bases were removed using Trimmomatic. A comparative quality assessment between Human-Enriched and Poly-A RNA-seq libraries in both pre-processed and trimmed datasets revealed several differences across different sequencing metrics, reflecting intrinsic characteristics of the two library preparation strategies. The most remarkable difference between the two methods is related to the duplication rates. Human-Enriched library displayed higher levels of duplication, with a mean rate of 84.04% in the pre-processed data and 85.42% in post-trimming data (Supplementary Table S2). In contrast, Poly-A libraries exhibited markedly lower duplication, with an average of 61.88% pre-trimming and a moderate increase to 67.89% after trimming. GC content also differed between the two library types. Human-Enriched library maintained a consistently higher GC percentage, averaging 52%–53% across both pre-trimmed and trimmed data, but, in contrast, Poly-A library exhibited a lower GC content, ranging from 42.95% in the pre-processed data to 44.1% post-trimming (Supplementary Table S2). Finally, a substantial reduction in failed reads was observed following trimming in both library types.

To evaluate mapping efficiency and alignment quality, the output log files from STAR were aggregated and analyzed using MultiQC. Overall, poor alignment performance was observed in both the Human-Enriched and Poly-A datasets, characterized by a high proportion of reads that failed to map to the human reference genome (Figures 2A,B). This result was not unexpected, given the RNA source, saliva, which contains a substantial amount of microbial and other non-human genetic material. These unmapped reads were subsequently analyzed through a dedicated metagenomic pipeline to investigate the composition and modulation of the oral microbiome in response to musical stimulation. A comparative assessment of the two datasets revealed notable differences attributable to their respective library preparation methods. The Poly-A dataset generally yielded a higher total number of reads, both mapped and unmapped, as shown in Figure 2C (see also Supplementary Table S3). However, despite this greater sequencing depth, a smaller proportion of reads in the Poly-A dataset successfully aligned to human genes compared to the Human-Enriched dataset (Figure 2D).

Figure 2
Graphs comparing biological datasets. Panel A shows a scatter plot with a correlation of 0.93 between Human-Enriched and Poly-A datasets. Panel B features correlation and histogram comparisons among expr_polya, expr_HT, and expr_merged datasets, with correlations ranging from 0.927 to 0.981. Panel C presents a volcano plot of Log₂ Fold Change versus negative Log2 P-value, highlighting specific organisms. Panel D displays a PCA plot with PC1 and PC2 variances of 81 and 14 percent. Panel E is a heatmap showing hierarchical clustering and z-scores for various species across two time points, TP1 and TP2.

Figure 2. Alignment metrics and number or reads per sample across datasets. Multimapping, non-mapping and uniquely mapping reads illustrating the alignment performance of the Human-Enriched dataset (A) and the Poly-A dataset (B) against the human genome (hg38). (C) Boxplots comparing the total number of reads in the Human-Enriched dataset and the Poly-A dataset. (D) Boxplots comparing the total number of reads aligned to the human genome (hg38) in the Human-Enriched dataset and the Poly-A dataset.

Further analysis showed that the Poly-A dataset achieved broader transcriptomic coverage, detecting approximately 20,390 genes. However, this broader coverage came at the cost of lower average read counts per gene. In contrast, the Human-Enriched dataset, as expected, showed a more focused mapping pattern: while fewer genes were detected (approximately 17,730), the number of reads mapping to each gene was higher, indicating increased sequencing depth per gene.

Differential gene expression analysis

A differential expression analysis was carried out to detect DEGs after the musical stimuli in the two salivary datasets separately. To improve statistical power and reduce noise, genes with fewer than 10 counts in at least five samples were filtered out. After applying this filtering criterion, 9,455 genes were retained in the Human-Enriched dataset, while only 6,791 genes were retained in the Poly-A dataset. The lower number of retained genes in the Poly-A dataset likely reflects its lower average read count per gene and sample.

When considering genes with nominal P-values <0.05, the Human-Enriched and Poly-A datasets yielded 720 and 442 genes, respectively. Of these, 57 genes were shared between the two datasets. Among these common genes, 55 were protein-coding, one was a long non-coding RNA, and one was a transcribed unprocessed pseudogene. To assess the concordance in expression trends, we compared Log2FC values for all DEGs identified in both datasets. The Log2FC values of all DEGs identified in either dataset (n = 738), regardless of whether they were detected in one or both, showed strong agreement in the direction of gene expression changes (i.e., genes upregulated in Poly-A were also found to be upregulated in the Human-Enriched subset, and vice versa). This concordance was supported by high-performance metrics, including an accuracy of 0.85, precision of 0.87, recall of 0.90, an F1-score of 0.88, and a Cohen’s kappa coefficient of 0.68 (Supplementary Table S4). In addition, when examining only the shared DEGs detected in the two datasets (n = 57), these metrics improved substantially, with an accuracy of 0.98, precision of 0.98, recall of 1.00, an F1-score of 0.99, and a Cohen’s kappa coefficient of 0.92 (Supplementary Table S4). Finally, the concordance in Log2FC between the 57 shared DEGs revealed a significant correlation between the Human-Enriched and Poly-A datasets (R = 0.77; Figure 3A). These results suggest that, despite technical differences and inherent biases between the two library preparation methods, the integration strategy used to combine both datasets successfully captured overlapping transcriptional responses to the musical stimulus.

Figure 3
Panel A shows a volcano plot with Log2 fold change on the x-axis and negative Log₁₀ P-value on the y-axis, highlighting significant genes such as REM2, HERC6, and TSPAN5. Panel B presents box plots of normalized expression for HERC6, TSPAN5, and REM2 at two time points (TP1 and TP2). Panel C is a PCA plot of 50 common genes, with PC1 explaining 86% and PC2 7% of the variance. Panel D features a heatmap of 50 common genes, with clustering by time points TP1 and TP2, displaying varying expression levels.

Figure 3. Correlation analyses between significant expression changes (TP1 vs. TP2) from different datasets. (A) Plot showing the correlation of Log2FC values from 57 differentially expressed genes (P-value <0.05) common between the Human-Enriched and Poly-A datasets (correlation = 0.77; 95% confidence interval). (B) Correlation between Log2FC values from 50 differentially expressed genes (P-value <0.05) shared between the merged dataset and the 57 common genes previously identified from the analysis of Human-Enriched and Poly-A datasets.

Integration through a merged dataset and cross-dataset validation of expression trends

To further optimize the analysis, a merged dataset was created by summing the read counts of shared genes across samples that were present in both datasets. This integration approach aimed to combine the broader gene coverage offered by the Poly-A dataset with the higher per-gene sequencing depth of the Human-Enriched dataset. Each dataset was filtered individually using the same criteria (genes with ≥10 counts in at least five samples). Afterward, counts from both datasets were merged by summing their values, followed by normalization and differential expression analysis using DESeq2 (Supplementary Table S5). To further confirm the consistency of expression patterns in the two datasets, genes with P-values below 0.05 in the merged dataset (n = 433) were intersected with the 57 genes previously found to be shared between the Human-Enriched and Poly-A datasets. A total of 50 genes were common to all three analyses. The Log2FC values for these genes were again compared across datasets. The Pearson correlation coefficient was 0.86 when comparing the Human-Enriched and Poly-A libraries, and exceeded 0.90 when each was compared individually with the merged dataset (Human-Enriched vs. Merged: 0.988; Poly-A vs. Merged: 0.901). All three comparisons were highly statistically significant (Figure 3B). These correlation values indicate strong consistency in expression patterns and further support the validity of the merging approach.

The merged dataset revealed three DEGs after musical stimuli when compared to baseline, with adjusted P-values <0.1 (Figure 4A). Other genes previously associated to ASD were also significantly affected by the musical stimuli. Two genes, HERC6 and TSPAN5, were upregulated, while REM2 was downregulated (Figure 4B). Gene expression profiles of the 50 common genes could perfectly separate TP1 from TP2 in both the PCA (Figure 4C) and the heatmap analysis (Figure 4D).

Figure 4
Panel A shows a volcano plot indicating gene expression changes, with points highlighted in red and blue representing significant adjusted p-values. Panel B is a bar chart of neural pathways with log2 Fold Changes, where red and blue bars indicate pathways with an adjusted p-value below 0.05. Panel C displays a bar chart of the top ten GO terms, colored by p-adjust values, showing counts for immune-related terms.

Figure 4. Analysis of DEGs and their relevance to autism-related pathways. (A) Volcano plot of the merged dataset, highlighting genes listed in the SFARI database and encircling in red the three genes with adjusted P-value <0.1 that are not included in the SFARI database. The plot shows the Log2FC on the x-axis and the -log10 P-value on the y-axis. Genes with a P-value <0.05 and a |Log2FC| > 0.5 are colored in red, genes with a P-value <0.05 but a |Log2FC| < 0.5 in light blue, genes with an P-value >0.05 but a |Log2FC| > 0.5 in green and non-significant genes in dark grey. (B) Boxplot showing the expression distribution before (TP1) and after (TP2) the musical stimuli of the 3 DEGs listed in the SFARI database in the merged dataset (**P-value = 0.0079). (C) Principal Component Analysis (PCA) of the 50 differentially expressed genes (P-value <0.05) shared between the merged dataset and the 57 genes previously found in common between the Human-Enriched and Poly-A datasets (P-value <0.05), illustrating sample distribution based on gene expression profiles. (D) Heatmap displaying the expression patterns of the 50 genes from the intersection between the differentially expressed genes (P-value <0.05) in the merged dataset and the 57 genes previously found in common between the Human-Enriched and the Poly-A datasets (P-value <0.05).

Functional interpretation and pathway exploration

Although the number of statistically significant DEGs identified in this study was limited, the high consistency observed across multiple datasets and analytical approaches strongly supports the hypothesis that music exposure can elicit measurable transcriptional changes in the buccal mucosa of individuals with ASD. To further explore the potential biological relevance of these changes, we performed pathway-level analyses using both GSVA and ORA. These complementary methods allowed us to move beyond single-gene effects to identify broader biological pathways and molecular functions modulated by musical intervention.

The differential pathway analysis revealed 137 DRPs in ASD patients (adjusted P-value <0.05), of which 58% (n = 79) were upregulated and 42% (n = 58) were downregulated following music exposure (Figure 5A; Supplementary Table S6). Among the most significantly affected pathways (|Log2FC| > 0.5 and adjusted P-value <0.01), we identified biological processes such as aorta development, actin filament severing, endoplasmic reticulum to cytosol transport, and skeletal muscle tissue regeneration. Notably, 13 of the 137 DRPs were associated with neurobiological functions (Figure 5B; Supplementary Table S6).

Figure 5
Panel A and B show bar charts of read counts in millions across samples, categorized by read type: multimapped, non-mapping, and uniquely mapped. Panel C and D display box plots of total reads and total aligned reads per dataset for Human-Enriched and poly-A, respectively, highlighting variations between the two groups.

Figure 5. Comprehensive functional enrichment analysis of the merged dataset across multiple approaches. (A) Volcano plot of differential pathways analysis identified via Gene Set Variation Analysis (GSVA); The plot shows the Log2FC on the x-axis and the -log10 adjusted P-value on the y-axis. Genes with an adjusted P-value <0.05 and a |Log2FC| > 0.5 are colored in red, genes with an adjusted P-value <0.05 but a |Log2FC| < 0.5 in light blue, genes with an adjusted P-value >0.05 but a |Log2FC| > 0.5 in green and non-significant genes in dark grey. The numbers in the plot correspond to GSVA-enriched (adjusted P-value <0.05) neuronal biological functions (see Supplementary Tables S1, S5). (B) Barplot of significantly GSVA-enriched pathways related specifically to neuronal biological functions (adjusted P-value <0.05); see Supplementary Tables S1, S5 for details. (C) Barplot summarizing the top 10 significantly enriched pathways (adjusted P-value <0.05) identified through Over-Representation Analysis (ORA).

Interestingly, in contrast to the global trend toward upregulation, neurobiological pathways showed a predominance of downregulation: 12 out of 13 exhibited negative Log2FC values post-stimulation, suggesting reduced activity after music exposure. These included key processes such as chemical synaptic transmission, postsynaptic (Log2FC = −0.45, adjusted P-value = 0.015), regulation of neuroblast proliferation (Log2FC = −0.40, adjusted P-value = 0.015), and axo-dendritic transport (Log2FC = −0.31, adjusted P-value = 0.018), pointing toward a potential dampening of excitatory synaptic signaling and neurodevelopmental activity. The only upregulated neurobiological pathway (Log2FC = 0.31, adjusted P-value = 0.015) was the cellular response to corticosteroid stimulus. Additionally, regulation of GABAergic synaptic transmission was notably downregulated (Log2FC = −0.62, adjusted P-value = 0.030). Other affected pathways involved in synaptic plasticity, neurodevelopment, and neurodegeneration, such as regulation of amyloid-beta clearance, synaptic vesicle localization, and neuron migration, were also significantly suppressed.

Functional analysis using the ORA approach identified 19 differentially regulated pathways (Figure 5C; Supplementary Table S7), many of which were associated with immune-related processes. These pathways included those involved in immune T-cell proliferation and activation, migration and adhesion, as well as the maintenance of immune homeostasis.

Comparison between both unsupervised (GSVA) and ORA approaches revealed only a common DRP, namely “calcium ion import across plasma membrane”, which was downregulated after the musical stimuli (Log2FC = −0.74, adjusted P-value = 0.015). However, an additional pathway related to T-cell homeostasis also emerged from the GSVA results (“regulation of t cell differentiation in thymus”; Log2FC = −0.43, adjusted P-value = 0.04), reinforcing the involvement of T-cell-related pathways in the response to music.

Co-expression modules analysis in response to musical stimulation

WGCNA analysis was conducted on a total of 4,667 genes from the merged dataset, after filtering out less variable genes. The analysis detected 28 modules of co-expressed genes, and among them, eight modules were significantly correlated with the musical stimuli. Multiple test correction was applied (FDR < 0.05), resulting in three statistically significant correlated modules (|R| > 0.8); Supplementary Table S8. Two of these, bisque4 (hub gene: BCL10) and saddlebrown (UBE2D3), exhibited strong positive correlation, while one module, plum1 (AKNA), showed a strong negative correlation with the music stimuli (Figure 6A). Other modules, such as darkolivegreen (CDK2AP1) and steelblue (PIK3C2A), also demonstrated notable correlation (r = −0.75 and R = 0.76, respectively), although they did not retain statistical significance after correction (adjusted P-value = 0.065 for both). Equally, other additional modules, including ivory (C1GALT1), thistle1 (UXS1), and lightgreen (FYN) exhibited moderate correlations (r ranging from 0.64 to 0.67), but did not reach significance following adjustment (Supplementary Table S8). Correlation between module membership (MM) and gene-level correlation with the musical stimuli of the three significant modules demonstrated a robust and statistically significant association, supporting the idea that core genes within these modules are functionally aligned with the musical response (Figure 6B). The most substantial association was observed in the BCL10 module (r = 0.76; P-value = 6.3e-09). Heatmap analysis of gene expression patterns from genes included in each of the significantly associated modules evidenced their overall over-regulation (BCL10 and UBE2D3 modules) and under-regulation (AKNA module) after listening to music (Figure 6C). Eigengene values of individual samples also align with the module-trait correlation values and changes in expression patterns observed between TP1 and TP2 (Figure 6D).

Figure 6
Workflow diagram of ASD donors at two time points with musical stimuli. mRNA extraction transfers Poly-A and Human-Enriched mRNA for sequencing and mapping. Poly-A dataset undergoes DGE analysis of mapped reads, while both datasets undergo DAA analysis of unmapped reads. Poly-A and Human-Enriched datasets are tested for expression value correlation, resulting in a merged dataset.

Figure 6. Co-expression analysis in ASD donors. (A) Heatmap of correlation values between modules of co-expressed genes detected and the musical stimuli (TP1-TP2). The upper value corresponds to the individual correlation value of each of the modules. P-values of these correlations are displayed in brackets. (B) Comparison between MM (module membership) and musical stimuli correlation of genes from the most significant module detected (adjusted P-value <0.05). Only names from genes showing an MM > 0.9 and a correlation with the musical stimuli >0.9 are displayed. (C) Heatmap of gene expression profiles from the modules most significantly correlated with the musical stimuli (adjusted P-value <0.05). (D) Raincloud plots of differences in individual samples’ eigengenes between TP1 and TP2 from modules showing statistical significance (adjusted P-value <0.05). (E) Biological processes detected from most statistically significant modules (adjusted P-value <0.05). Modules with significant pathways are named as the hub gene names.

Functional enrichment analysis identified different processes associated with the significantly correlated modules AKNA, FYN, UBE2D3 and UXS1 (Supplementary Table S9). Notably, the AKNA and UBE2D3 modules were among the three modules significantly correlated with the musical stimuli after false discovery rate (FDR) correction (FDR < 0.05); Figure 6E. The AKNA module, which was downregulated in response to the stimuli, was primarily enriched for processes involving GTPase-mediated signal transduction, including RAS protein signaling pathways. Additionally, this module was also associated with immune-related functions such as myeloid cell differentiation and activation, B-cell differentiation and signaling, and the regulation of inflammatory responses, including interleukin-1 (IL-1) production (Figure 6E). On the other hand, only two enriched pathways were identified among the genes included in the upregulated UBE2D3 module, namely “response to endoplasmic reticulum stress” and “cellular response to unfolded protein”. ORA analysis of the genes belonging BCL10 module did not reveal any significantly enriched pathways.

In addition, we further explored the functional roles of the hub genes from the significantly correlated modules (AKNA, UBE2D3 and BCL10) by retrieving the GO terms in which these genes are involved (Supplementary Table S10; Supplementary Figure S2). The AKNA gene is associated with 27 pathways, primarily related to nervous system processes (i.e. neurogenesis, forebrain development or stress response) as well as cell-cell adhesion. UBE2D3 gene participates in 50 pathways, mainly linked to ubiquitination, protein-related catabolic processes and response to BMP and TGF-beta signaling. Finally, BCL10 is involved in 157 pathways, predominantly associated with immune-related functions, cell–cell adhesion, and apoptotic processes.

ASD-related genes vs. genes altered by musical stimuli

To explore the potential relevance of musical stimulation to autism-related molecular mechanisms, we compared the gene expression data from our merged dataset with the SFARI gene database. Of the 6,305 genes identified, 519 matched entries in the SFARI database. Among these, ten genes showed significant differential expression (P-value <0.05 and |Log2FC| > 1): nine were upregulated and one was downregulated. Notably, genes such as EIF4E, DLGAP1, and CDKL5 were among the upregulated group, while PTBP2 was significantly downregulated. Several of these genes have strong or suggestive associations with autism, according to their SFARI scores (e.g., EIF4E–Score 2; CDKL5-Score 1S). The identified genes are implicated in various processes, including synaptic signaling, protein translation, and RNA processing, suggesting that music stimulation may influence molecular pathways commonly affected in ASD (Supplementary Table S11).

Exploratory analysis of salivary microbiota following musical stimulation

We also investigated whether musical stimulation induces detectable changes in the salivary microbiota by performing metagenomic inference on reads that failed to align to the human genome. As expected, a substantial proportion of these unmapped reads could not be matched to any bacterial reference genome (Supplementary Figure S3). This result is consistent with the complex and heterogeneous composition of saliva, which contains not only bacteria but also fungi, viruses, dietary RNA fragments, and environmental contaminants that may not be represented in the reference database used.

Applying a filtering criterion consistent with the human differential gene expression analysis (retaining only species with ≥10 reads in at least five samples), the number of species was reduced from 3,210 (Human-Enriched dataset) and 3,345 (Poly-A dataset) to 1,340 and 1,767 species, respectively. A total of 1,273 species were shared between both datasets.

Differential abundance analysis (DAA) was conducted using DESeq2 in the two datasets separately. No bacterial species passed the adjusted significance threshold (adjusted P-value <0.05) in either dataset. However, 52 species in the Human-Enriched dataset and 62 in the Poly-A dataset showed nominal significance (P-value <0.05), with 15 species overlapping between the two datasets (Supplementary Table S12). A Pearson correlation analysis of the Log2FC for these 15 shared species revealed a high concordance (r = 0.93), indicating strong agreement between the independent analyses (Figure 7A).

Figure 7
Composite image showing various data visualizations. Panel A: Heatmap of gene expression with color scale from blue (negative) to red (positive). Panel B: Scatter plots show module membership versus gene significance for three genes: AKNA, BCL10, UBE2D3. Panel C: Heatmaps for AKNA, BCL10, and UBE2D3 modules indicate phenotype differences with color-coded bars for TP1 and TP2. Panel D: Box plots compare module eigengenes between TP1 and TP2 for three genes, with p-values noted. Panel E: Chart lists biological processes, with a bubble plot indicating overrepresentation and gene ratios.

Figure 7. Differential abundance analysis of salivary microbial species pre- and post-musical stimulation in ASD subjects. (A) Scatter plot showing the Pearson correlation of Log2FC values for 15 microbial species with differential abundance (P-value <0.05) in both the Human-Enriched and Poly-A datasets. (B) Correlation between Log2FC values from the species from the intersection between differentially abundant species (P-value <0.05) from the merged dataset and the species with differential abundance previously found in common between the analysis of Human-Enriched and the Poly-A datasets, which were 15. (C) Volcano plot of the merged dataset displaying species with significant differential abundance. The plot shows the Log2FC on the x-axis and the -log10 P-value on the y-axis. Species with a P-value <0.05 and a |Log2FC| > 0.5 are colored in red, species with a P-value <0.05 but a |Log2FC| > 0.5 in light blue, species with an P-value >0.05 but a |Log2FC| > 0.5 in green and non-significant species in dark grey. (D) PCA plot based on abundance profiles from the merged dataset of the 15 commonly differentially abundant species. (E) Heatmap of the relative abundances from the merged dataset of the 15 commonly differentially abundant species.

Subsequently, a merged dataset was generated by aggregating read counts across datasets for the 1,273 shared species. This merged dataset showed also a high correlation in Log2FC values with the Poly-A dataset and Human-Enriched datasets for the 15 shared species between the two libraries (Figure 7B). DAA identified 68 species with significantly different abundances at TP2 compared to TP1 (P-value <0.05; Supplementary Table S13; Figure 7C). Notably, all 15 species previously found in common between the Human-Enriched and Poly-A datasets were also detected in the merged dataset. Changes in the abundance of these 15 species effectively distinguished the pre-stimulation and post-stimulation states in ASD subjects, as demonstrated by both PCA (Figure 7D) and heatmap analysis (Figure 7E). Unlike the separate analyses, this integrated approach identified one additional differentially abundant species, Acidipropionibacterium acidipropionici, with an adjusted P-value of 0.076 and a Log2FC of 2.5 (Supplementary Table S13).

Discussion

Although RNA-seq is a powerful technique for biomarker discovery, saliva poses significant challenges due to its high bacterial RNA content and low abundance of human RNA, requiring specialized extraction and enrichment protocols (Ostheim et al., 2020; Gosch et al., 2024). In this study, we investigated the transcriptional effects of music exposure in the buccal mucosa of individuals with ASD using two complementary RNA-seq datasets prepared via Poly-A selection and Human-Enriched library protocols. Human-Enriched and Poly-A RNA-seq libraries exhibited distinct sequencing characteristics that align with their differing biological targets and preparation protocols. The aforementioned heterogeneous composition of saliva directly impacts the alignment metrics, resulting in lower mapping rates compared to the expected 70%–90% typically observed for RNA-seq reads from other sample types, such as blood (Dobin et al., 2013). Human-Enriched library showed higher duplication and GC content compared to Poly-A library preparation. Although it detected fewer genes, it achieved higher read counts per gene, reflecting greater sequencing depth for targeted transcripts (Conesa et al., 2016). This is expected, as Poly-A generally captures polyadenylated transcripts for an unbiased, transcriptome-wide analysis of coding genes, while Human-Enriched relies on hybrid capture with a targeted panel (e.g., the Exome Panel), limiting analysis to predefined genomic regions. Although each dataset individually revealed a poor mapping performance against the human genome and a modest number of DEGs, the integration of both datasets enabled the identification of a more consistent transcriptional signature, with enhanced robustness and biological interpretability. Reproducibility among technical replicates is generally expected to be high (Mortazavi et al., 2008). Despite the inherent technical differences between the Poly-A and Human-Enriched library preparation methods, our comparative analysis demonstrated a strong concordance in gene expression patterns between the two datasets. The moderate-to-high concordance in gene expression trends between the Poly-A and Human-Enriched datasets, and the high consistency in the merged dataset, reinforce the validity of the transcriptional changes observed and the integration strategy applied. Although the number of DEGs detected at stringent thresholds was limited, likely due to sample size constraints and biological variability inherent in human studies, the robustness of these shared signals underscores the value of integrating different RNA-seq strategies.

Three genes, HERC6, TSPAN5, and REM2, emerged as significant from the merged Poly-A selection and Human-Enriched datasets, suggesting they may play important roles in the brain’s response to musical stimuli in individuals with ASD. The upregulation of HERC6 may point to a role for immune-related pathways in mediating the effects of music exposure. HERC6 encodes an E3 ubiquitin-protein ligase and is a member of the small HERC family, which, unlike its larger counterparts (HERC1, HERC2), is thought to have evolved independently. While in many species HERC6 participates in type I interferon responses by facilitating ISGylation, in humans ISGylation activity is mediated by HERC5 (Sanchez-Tena et al., 2016). However, emerging evidence suggests that HERC6 may still play roles in antiviral defense and immune regulation through alternative, non-ISGylation-based mechanisms (Oudshoorn et al., 2012). The modulation of such immune-related genes is particularly relevant given accumulating evidence that immune dysregulation is a common feature in ASD. Music, by influencing stress levels and immune signaling pathways, could indirectly affect the local expression of genes such as HERC6, potentially contributing to homeostatic recalibration in neuroimmune interactions. Perhaps more directly tied to neurodevelopmental processes, the upregulation of TSPAN5 presents a compelling finding. TSPAN5 encodes a glycoprotein of the tetraspanin family, which is involved in a broad array of cellular functions including adhesion, migration, and proliferation. Of particular importance is its role in the central nervous system, where it is highly expressed in the forebrain and cerebellum. Recent studies have demonstrated that TSPAN5 facilitates the maturation of dendritic spines by promoting the synaptic clustering of neuroligin-1 (NLG-1) (Moretto et al., 2019), a protein that has been closely linked to ASD-related social behavior deficits (Rawsthorne et al., 2021). The upregulation of TSPAN5 in the buccal mucosa of ASD patients following musical intervention may thus reflect a shift toward synaptic stabilization and maturation occurring in the brain. This aligns with previous reports indicating that music-based therapies can improve cognitive and social functioning in individuals with ASD, potentially through activity-dependent plasticity and enhanced synaptic integration. Conversely, REM2 was significantly downregulated after music exposure. REM2 belongs to the RGK subfamily of Ras-like GTPases and plays a critical role in neuronal differentiation and synaptic development. However, intriguingly, inhibition of REM2 has also been associated with increased dendritic branching and arborization, resulting in more structurally complex neurons (Ghiretti and Paradis, 2011). This dualistic effect complicates the interpretation of its downregulation, yet suggests that reduced REM2 expression may contribute to neuroplastic adaptations rather than dysfunction per se. In the context of this study, the simultaneous upregulation of TSPAN5 may offset potential drawbacks of REM2 suppression, together fostering a structural environment conducive to improved neural integration and function.

The concurrent modulation of HERC6, TSPAN5, and REM2 in the buccal cavity of ASD individuals after music exposure appear to span multiple functional domains, including immune regulation, synaptic maturation, and structural remodeling, all of which are relevant to the ASD phenotype. While the directionality and interplay between these pathways remain to be fully elucidated, our findings resonate with a growing body of literature suggesting that music can serve as a multisensory and emotionally salient stimulus capable of driving neuroplasticity.

Numerous studies have explored the relationship between ASD and alterations in the synaptic excitation-inhibition (E/I) ratio. Some have reported an increased E/I ratio in individuals with ASD (Rubenstein and Merzenich, 2003), while others, particularly in Down syndrome, have suggested that heightened inhibition may be the primary contributor to learning deficits (Belichenko et al., 2009). Collectively, these findings seem to indicate that homeostatic regulation of neural activity is often disrupted in ASD. Given the heterogeneity of the disorder, it is plausible that both hyperexcitation and hyperinhibition could give rise to distinct ASD subtypes or manifestations. GSVA revealed a coordinated downregulation of neurodevelopmental and synaptic pathways in the buccal mucosa of ASD patients after music exposure, contrasting with the global upregulation pattern observed across all significant pathways. We hypothesize that the observed downregulation of most neural DRPs following musical stimuli may represent a homeostatic response aimed at modulating excessive synaptic excitation, a characteristic frequently observed in ASD. This response could reflect a reduction in the E/I ratio, underscoring the potential of music as a modulator of neural excitability. The cellular response to corticosteroid stimulus was the only significantly upregulated pathway, which may represent a stress-related or neuroendocrine response to the stimulus. Corticosteroids are hormones involved in stress regulation and immune modulation, mainly through the hypothalamic–pituitary–adrenal axis (HPA) (Agorastos et al., 2019; Wong et al., 2022). HPA axis dysregulation is a well-documented feature of ASD (Spratt et al., 2012). HPA activity is correlated with cortisol levels, an indicator of stress. ASD often exhibits atypical cortisol secretion patterns and blunted or exaggerated stress reactivity (Spratt et al., 2012; Taylor and Corbett, 2014). Our data suggest that music exposure may transiently enhance corticosteroid-responsive gene expression, reflecting a temporary adjustment of stress-adaptive mechanisms rather than a stress response in these individuals. Further studies are required to clarify the functional significance of this pathway activation and its potential role in mediating the effects of music in ASD. In contrast, our ORA revealed significant enrichment in pathways associated with immune processes, particularly those involving the proliferation and differentiation of immune cells such as leukocytes, lymphocytes, and mononuclear cells. Immune dysregulation is a well-established co-occurring condition in ASD, with several studies reporting elevated levels of TNF-α, IL-6, granulocyte colony-stimulating factor (G-CSF), IFN-γ, and IL-8 in the brain, along with a heightened pro-inflammatory Th1 cytokine profile (Li et al., 2009; Wei et al., 2011; Robinson-Agramonte et al., 2022). Music exposure modulates immune-related pathways, suggesting an attenuation of pro-inflammatory states, supporting previous findings that music can exert immunoregulatory effects and stress-related responses.

While further mechanistic studies are warranted, our GSVA and ORA results support the notion that sensory stimuli such as music may trigger compensatory transcriptional programs in individuals with ASD, influencing neural excitability and immune function. Cross-referencing music-responsive DEGs with the SFARI gene database revealed modulation of genes involved in synaptic structure, intracellular signaling, and RNA metabolism, pathways commonly disrupted in ASD. Notably, several upregulated genes, including DLGAP1 and CDKL5, play critical roles in synaptic scaffolding and neuronal morphogenesis (Chen et al., 2010; Shin et al., 2012). DLGAP1 encodes a postsynaptic protein localized at glutamatergic synapses, while CDKL5 is essential for the stability of excitatory synapses (Chen et al., 2010). The altered expression of these genes suggests a potential link between musical stimulation and enhanced synaptic plasticity in ASD-related contexts. In addition, genes such as EIF4E, STYK1, and CDKL5 are components of the PI3K/Akt/mTOR signaling pathway, a key regulator of synaptic protein synthesis and neural development. Importantly, we also observed the downregulation of PTBP2, a central regulator of neuronal exon splicing and RNA stability. Along with other RNA-binding genes such as UPF3B and FRG1, this finding indicates that music may exert broader effects on RNA processing and transcript diversity (Sun et al., 2011; Maracci et al., 2022). Taken together, these transcriptional changes suggest that music stimulation can reprogram gene expression patterns across interconnected molecular pathways, potentially contributing to neuroplastic and compensatory mechanisms relevant to ASD.

We identified different co-expression modules whose activity patterns were significantly correlated with exposure to musical stimuli in ASD donors. This suggests that musical stimuli may engage specific regulatory networks in the buccal mucosa, potentially reflecting broader systemic or brain-related transcriptional responses. From the three main significantly correlated modules, BCL10 and UBE2D3 showed a global upregulation pattern in response to music stimulation, while the AKNA module demonstrated an overall downregulation. Notably, hub gene AKNA has already been previously reported as a novel ASD risk gene, carrying de novo variants probably related to ASD (Kim et al., 2020). AKNA module was functionally enriched for pathways related to GTPase-mediated signal transduction, with a particular emphasis on Ras protein signaling, suggesting that musical stimulation affects specific intracellular signaling cascades related to these small GTPases. The Rho family of GTPases (members of the larger Ras superfamily) are small signaling proteins that play pivotal roles in neural development, including axon guidance, dendritic spine formation, and synaptic plasticity processes (Zhang et al., 2021). Importantly, dysregulation of Rho GTPase signaling has been increasingly recognized as a contributing factor in various neuropathological events, including ASD (Arrazola Sastre et al., 2020; Guo et al., 2020). Multiple lines of evidence have highlighted aberrant activity of Rho family members, such as Rac1, in ASD-related phenotypes both in human studies and animal models (Zeidan-Chulia et al., 2016). Interestingly, the hub gene of this module (AKNA) is directly involved in molecular pathways such as neurogenesis or forebrain development, indicating a potential regulatory role of this module in nervous system-related processes. In addition, this module was also enriched for genes involved in immune-related processes. The immunomodulatory effect of the musical stimuli appears to act through the suppression specific immune cells populations. The downregulation of inflammatory and IL-1 pathways could be of particular interest given the role of neuroinflammation in the pathophysiology of ASD (Ricci et al., 2013; Careaga et al., 2017; Saghazadeh et al., 2019). Module analysis also revealed a selective activation of stress-response mechanisms (via the UBE2D3 module). The accumulation of unfolded proteins in the endoplasmic reticulum (ER) triggers ER stress, specifically activating a protective unfolded protein response in the ER. Several studies have identified ER stress as a significant contributor to the pathophysiology of ASD (Fujita et al., 2010). ER stress has also been shown to impair neurite outgrowth, alter synaptic protein expression, and disrupt neuronal differentiation, all of which are relevant to ASD neurodevelopmental deficits (Kawada and Mimori, 2018). In addition, upregulation of ER stress-related genes has been reported in the brains of ASD individuals (Crider et al., 2017; Kawada and Mimori, 2018). Our results indicate that music exposure may enhance ER stress response mechanisms in buccal cells, potentially regulated by genes within UBE2D3 module, which could help to compensate for the impact of stress pathways implicated in ASD.

While the metagenomic analysis conducted in this study was exploratory and limited in scope, it yielded several notable observations. The large fraction of sequencing reads that failed to align with bacterial reference genomes highlights the inherent complexity of salivary metagenomic samples. Saliva is a biologically diverse medium that includes bacteria, fungi, viruses, dietary residues, and environmental RNA fragments, many of which are absent from current bacterial reference databases. Dysbiosis of the oral microbiome has previously been linked to ASD, potentially contributing to neurodevelopmental alterations via the gut–brain axis or immune modulation (Wong et al., 2021; D'Angelo et al., 2025). A subset of taxa showed consistent nominal differences between pre- and post-stimulation samples across both datasets. In particular, 15 species demonstrated shared differential abundance patterns with high correlation in Log2FC (r = 0.93), suggesting reproducibility across sequencing approaches (Human-Enriched vs. Poly-A datasets). Notably, the merged dataset identified Acidipropionibacterium acidipropionici as a potentially modulated species following musical stimulation (adjusted P-value = 0.076). This acid-producing bacterium is known for its role in dairy fermentation (Coronas et al., 2023), and closely related propionibacteria have been previously identified as part of the oral microbiota (Sato-Suzuki et al., 2020). We also identified Propionibacterium freudenreichii, a close relative involved in propionic acid production, among our top 15 candidate species. There is substantial evidence implicating propionibacteria in neurological effects via propionic acid production, an important short-chain fatty acid. Elevated propionic acid levels in the gut have been associated with ASD-like behaviors and neuroinflammation in animal models (Srikantha and Mohajeri, 2019). According to Li and Zhou (2016), 3-(3-hydroxyphenyl)-3-hydroxypropionic acid induces autism symptoms by emptying and depleting catecholamines in the brain. El-Ansary et al. (2012), El-Ansary et al. (2015) reported that the neurotoxicity of propionic acid could play a central role in the etiology of autistic biochemical features. This finding raises the intriguing possibility that auditory stimulation may exert an “echo effect” on the composition of the oral microbiome. These results invite further inquiry into the interplay between sensory stimuli and microbial communities. It remains unclear whether the observed microbial shifts are directly mediated by neural or physiological responses to music, or indirectly through changes in salivary secretion, flow rate, or composition. Future metagenomic studies with larger cohorts, increased sequencing depth, and broader reference annotations, including non-bacterial taxa, are needed to validate and expand on these preliminary findings. Further research should explore the potential of targeting the microbiota–gut–brain axis through musical stimuli as a therapeutic strategy in ASD, in parallel with approaches such as probiotic supplementation (Li and Zhou, 2016; Feng et al., 2023), or more broadly, interventions that manipulate the enteric microbiome to alleviate autism symptoms (Frye et al., 2015).

One of the most significant challenges encountered in the present study was the use of saliva as the biological matrix for transcriptomic analysis. While saliva sampling offers clear advantages, being non-invasive, stress-free, and highly accessible, particularly in sensitive populations such as individuals with ASD, it also presents substantial technical hurdles. Saliva contains a heterogeneous mix of host epithelial cells, immune cells, and a rich microbial community, which can result in low yields of high-quality human RNA and high proportions of unmapped or non-human reads. These factors complicate both the quality control and the interpretation of gene expression data. Despite these challenges, this study demonstrates that with careful preprocessing, rigorous quality filtering, and thoughtful experimental design, it is possible to extract biologically meaningful transcriptomic signals from saliva. This opens the door for future research into peripheral biomarkers of brain function and neuroplasticity in contexts where blood collection is difficult or ethically sensitive, such as pediatric or neurodiverse populations. To further advance the utility of saliva transcriptomics, future studies should aim to improve RNA stabilization protocols, integrate host-microbiome transcriptome analyses, and consider single-cell or spatial transcriptomics to better resolve cell-type-specific expression patterns. Combining saliva-based transcriptomics with behavioral, neuroimaging, or immune profiling could yield a more comprehensive understanding of how systemic biological pathways contribute to neurodevelopmental conditions like ASD.

This study is not without limitations. The sample size, while reasonable for a pilot analysis, limits the power to detect subtle effects and may also influence the stability of enrichment analyses, whose robustness depends on the reliability of the input gene list. Nonetheless, the enrichment of pathways consistently associated with ASD, such as immune regulation and mitochondrial function, supports the biological relevance of our findings. The exclusion of outlier samples was necessary to reduce noise but may have inadvertently excluded biologically informative variance. Additionally, the focus on peripheral tissue (saliva) restricts direct extrapolation to brain-specific mechanisms, although saliva has increasingly been recognized as a viable proxy for systemic and neuroimmune interactions. Future studies should build upon these findings by increasing cohort sizes, incorporating longitudinal designs, and using complementary modalities such as neuroimaging and behavioral assessments. Furthermore, investigating cell-type-specific transcriptomic responses or incorporating single-cell RNA-seq could help disentangle the contributions of diverse cell populations and improve the resolution of observed effects.

In conclusion, this study offers preliminary yet compelling evidence that music exposure can modulate the expression of genes involved in synaptic plasticity, neuronal architecture, and immune signaling in individuals with ASD. Notably, the upregulation of TSPAN5 and HERC6, alongside the downregulation of REM2, suggests a transcriptional profile that may favor enhanced neural adaptability and integration. These molecular findings enrich our understanding of the biological mechanisms through which music may exert therapeutic effects and point to promising avenues for developing targeted, non-invasive interventions in neurodevelopmental conditions. Furthermore, given the frequent reports of oral microbiota dysbiosis in ASD, the observed trends in microbial composition following musical stimulation, though exploratory, raise important questions about the broader systemic effects of sensory inputs. Understanding how music may influence the oral microbiome could open new perspectives on modulating host–microbe interactions through behavioral or environmental stimuli (Castelo-Martínez et al., 2025). Nevertheless, these initial findings should be interpreted with caution, and future research with larger cohorts and expanded metagenomic approaches will be essential to confirm and contextualize these effects.

Data availability statement

The raw count data from the human-merged and the metagenomic-merged datasets, along with the metadata, are publicly available in the Figshare repository under the accession: https://doi.org/10.6084/m9.figshare.29504588.

Ethics statement

The studies involving humans were approved by Ethics Committee of the Xunta de Galicia (registration code: 2020/021). The studies were conducted in accordance with the local legislation and institutional requirements. Written informed consent for participation in this study was provided by the participants’ legal guardians/next of kin. Written informed consent was obtained from the individual(s) for the publication of any potentially identifiable images or data included in this article.

Author contributions

AC: Investigation, Validation, Methodology, Visualization, Writing – review and editing, Supervision, Writing – original draft, Formal Analysis, Software, Data curation. NM: Writing – review and editing, Investigation, Data curation, Methodology, Software, Writing – original draft, Visualization, Formal Analysis, Validation. LN: Supervision, Conceptualization, Writing – review and editing. FM-T: Resources, Funding acquisition, Writing – review and editing, Conceptualization. AG-C: Supervision, Visualization, Investigation, Software, Conceptualization, Validation, Formal Analysis, Writing – review and editing, Methodology, Data curation, Writing – original draft. AS: Funding acquisition, Writing – review and editing, Project administration, Conceptualization, Supervision, Formal Analysis, Visualization, Writing – original draft, Methodology, Validation, Data curation, Investigation, Resources.

Funding

The author(s) declare that financial support was received for the research and/or publication of this article. This work was supported by: i) GAIN IN607B 2020/08 and IN607A 2023/02, and EUTERPE_adn (Programa de Cooperación Interreg-VI POCTEP; Ref. 0313_EUTERPE_ADN_1_E) (to AS), IIN607A2021/05 (to FM-T) and IN677D 2024/06 (to AG-C), and ii) Consorcio Centro de Investigación Biomédica en Red de Enfermedades Respiratorias (CB21/06/00103; to AS and FM-T). AG-C is supported by the Miguel Servet contract (CP23/00080), funded by the Instituto de Salud Carlos III (ISCIII) and co-funded by the European Union. The funders were not involved in the study design, collection, analysis, interpretation of data, the writing of this article, or the decision to submit it for publication.

Acknowledgments

Acknowledgements

The authors would like to express their appreciation to the patients and the study investigators of the Sensogenomics network [sensogenomics.com; Sensogenomics Working Group (see Supplementary File 1)], as well as the nursery and laboratory service at the Hospital Clínico Universitario de Santiago de Compostela, for their invaluable dedication and support. This research project was made possible through the access granted by the Galician Supercomputing Center (CESGA) to its supercomputing infrastructure. The supercomputer FinisTerrae III and its permanent data storage system have been funded by the Spanish Ministry of Science and Innovation, the Galician Government, and the European Regional Development Fund (ERDF).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Generative AI statement

The author(s) declare that no Generative AI was used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmolb.2025.1696704/full#supplementary-material

SUPPLEMENTARY FIGURE S1 | (A) Plots showing correlation between soft-thresholding powers and both the scale-free fit index and the mean connectivity (upper). (B) Dendrogram of genes and co-expression modules detected, represented by different colors.

SUPPLEMENTARY FIGURE S2 | Clustering of the semantic similarity matrix obtained from the Gene Ontology (GO; subcategory biological processes) terms in which hub genes from the significantly correlated modules are involved.

SUPPLEMENTARY FIGURE S3 | Pie charts representing the percentage of bacterial reads and unknown reads in the Human-Enriched dataset (A) and in the Poly-A dataset (B).

References

Abrahams, B. S., Arking, D. E., Campbell, D. B., Mefford, H. C., Morrow, E. M., Weiss, L. A., et al. (2013). SFARI Gene 2.0: a community-driven knowledgebase for the autism spectrum disorders (ASDs). Mol. Autism 4, 36. doi:10.1186/2040-2392-4-36

PubMed Abstract | CrossRef Full Text | Google Scholar

Agorastos, A., Heinig, A., Stiedl, O., Hager, T., Sommer, A., Muller, J. C., et al. (2019). Vagal effects of endocrine HPA axis challenges on resting autonomic activity assessed by heart rate variability measures in healthy humans. Psychoneuroendocrinology 102, 196–203. doi:10.1016/j.psyneuen.2018.12.017

PubMed Abstract | CrossRef Full Text | Google Scholar

American Psychiatric Association (2022). Diagnostic and statistical manual of mental disorders (5th ed., text rev.; DSM-5-TR). American Psychiatric Publishing.

Google Scholar

Anagnostou, E., Soorya, L., Chaplin, W., Bartz, J., Halpern, D., Wasserman, S., et al. (2012). Intranasal oxytocin versus placebo in the treatment of adults with autism spectrum disorders: a randomized controlled trial. Mol. Autism 3, 16. doi:10.1186/2040-2392-3-16

PubMed Abstract | CrossRef Full Text | Google Scholar

Arrazola Sastre, A., Luque Montoro, M., Galvez-Martin, P., Lacerda, H. M., Lucia, A. M., Llavero, F., et al. (2020). Small GTPases of the ras and rho families switch on/off signaling pathways in neurodegenerative diseases. Int. J. Mol. Sci. 21, 6312. doi:10.3390/ijms21176312

PubMed Abstract | CrossRef Full Text | Google Scholar

Bakker-Huvenaars, M. J., Greven, C. U., Herpers, P., Wiegers, E., Jansen, A., Van Der Steen, R., et al. (2020). Saliva oxytocin, cortisol, and testosterone levels in adolescent boys with autism spectrum disorder, oppositional defiant disorder/conduct disorder and typically developing individuals. Eur. Neuropsychopharmacol. 30, 87–101. doi:10.1016/j.euroneuro.2018.07.097

PubMed Abstract | CrossRef Full Text | Google Scholar

Belichenko, P. V., Kleschevnikov, A. M., Masliah, E., Wu, C., Takimoto-Kimura, R., Salehi, A., et al. (2009). Excitatory-inhibitory relationship in the fascia dentata in the Ts65Dn mouse model of Down syndrome. J. Comp. Neurol. 512, 453–466. doi:10.1002/cne.21895

PubMed Abstract | CrossRef Full Text | Google Scholar

Benach, J. L., Li, E., and Mcgovern, M. M. (2012). A microbial association with autism. mBio 3, e00019-12. doi:10.1128/mBio.00019-12

PubMed Abstract | CrossRef Full Text | Google Scholar

Bjørk, M., Riedel, B., Spigset, O., Veiby, G., Kolstad, E., Daltveit, A. K., et al. (2018). Association of folic acid supplementation during pregnancy with the risk of autistic traits in children exposed to antiepileptic drugs in utero. JAMA Neurol. 75, 160–168. doi:10.1001/jamaneurol.2017.3897

PubMed Abstract | CrossRef Full Text | Google Scholar

Blighe, K., Rana, S., and Lewis, M. (2020). EnhancedVolcano: publication-ready volcano plots with enhanced colouring and labeling. version 180, R package.

Google Scholar

Bolger, A. M., Lohse, M., and Usadel, B. (2014). Trimmomatic: a flexible trimmer for illumina sequence data. Bioinformatics 30, 2114–2120. doi:10.1093/bioinformatics/btu170

PubMed Abstract | CrossRef Full Text | Google Scholar

Careaga, M., Rogers, S., Hansen, R. L., Amaral, D. G., Van De Water, J., and Ashwood, P. (2017). Immune endophenotypes in children with autism spectrum disorder. Biol. Psychiatry 81, 434–441. doi:10.1016/j.biopsych.2015.08.036

PubMed Abstract | CrossRef Full Text | Google Scholar

Carlson, M. (2023). org.Hs.eg.db: genome wide annotation for human. Available online at: https://bioconductor.org/packages/org.Hs.eg.db/.

Google Scholar

Castelo-Martínez, L., Mallah, N., Cavenaghi, A., Navarro, L., Martinón-Torres, F., Gómez-Carballa, A., et al. (2025). Harmonizing the oral-brain axis: music-induced microbiota shifts in age-related cognitive disorders and healthy aging. medRxiv. doi:10.1101/2025.1107.1129.25332408

CrossRef Full Text | Google Scholar

Chen, Q., Zhu, Y. C., Yu, J., Miao, S., Zheng, J., Xu, L., et al. (2010). CDKL5, a protein associated with Rett syndrome, regulates neuronal morphogenesis via Rac1 signaling. J. Neurosci. 30, 12777–12786. doi:10.1523/JNEUROSCI.1102-10.2010

PubMed Abstract | CrossRef Full Text | Google Scholar

Christodoulides, N., Mohanty, S., Miller, C. S., Langub, M. C., Floriano, P. N., Dharshan, P., et al. (2005). Application of microchip assay system for the measurement of C-reactive protein in human saliva. Lab. Chip 5, 261–269. doi:10.1039/b414194f

PubMed Abstract | CrossRef Full Text | Google Scholar

Conesa, A., Madrigal, P., Tarazona, S., Gómez-Cabrero, D., Cervera, A., Mcpherson, A., et al. (2016). A survey of best practices for RNA-seq data analysis. Genome Biol. 17, 13. doi:10.1186/s13059-016-0881-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Cook, A., Ogden, J., and Winstone, N. (2019). The impact of a school-based musical contact intervention on prosocial attitudes, emotions and behaviours: a pilot trial with autistic and neurotypical children. Autism 23, 933–942. doi:10.1177/1362361318787793

PubMed Abstract | CrossRef Full Text | Google Scholar

Coronas, R., Zara, G., Gallo, A., Rocchetti, G., Lapris, M., Petretto, G. L., et al. (2023). Propionibacteria as promising tools for the production of pro-bioactive scotta: a proof-of-concept study. Front. Microbiol. 14, 1223741. doi:10.3389/fmicb.2023.1223741

PubMed Abstract | CrossRef Full Text | Google Scholar

Crider, A., Ahmed, A. O., and Pillai, A. (2017). Altered expression of endoplasmic reticulum stress-related genes in the middle frontal cortex of subjects with autism spectrum disorder. Mol. Neuropsychiatry 3, 85–91. doi:10.1159/000477212

PubMed Abstract | CrossRef Full Text | Google Scholar

Cryan, J. F., O'riordan, K. J., Cowan, C. S. M., Sandhu, K. V., Bastiaanssen, T. F. S., Boehme, M., et al. (2019). The microbiota-gut-brain axis. Physiol. Rev. 99, 1877–2013. doi:10.1152/physrev.00018.2018

PubMed Abstract | CrossRef Full Text | Google Scholar

D'angelo, E., Fiori, F., Ferraro, G. A., Tessitore, A., Nazzaro, L., Serpico, R., et al. (2025). Autism spectrum disorder, oral implications, and oral microbiota. Child. (Basel) 12, 368. doi:10.3390/children12030368

CrossRef Full Text | Google Scholar

Delobel-Ayoub, M., Saemundsen, E., Gissler, M., Ego, A., Moilanen, I., Ebeling, H., et al. (2020). Prevalence of autism spectrum disorder in 7-9-year-old children in Denmark, Finland, France and Iceland: a population-based registries approach within the ASDEU project. J. Autism Dev. Disord. 50, 949–959. doi:10.1007/s10803-019-04328-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Dobin, A., Davis, C. A., Schlesinger, F., Drenkow, J., Zaleski, C., Jha, S., et al. (2013). STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29, 15–21. doi:10.1093/bioinformatics/bts635

PubMed Abstract | CrossRef Full Text | Google Scholar

Domes, G., Heinrichs, M., Michel, A., Berger, C., and Herpertz, S. C. (2007). Oxytocin improves “mind-reading” in humans. Biol. Psychiatry 61, 731–733. doi:10.1016/j.biopsych.2006.07.015

PubMed Abstract | CrossRef Full Text | Google Scholar

El-Ansary, A. K., Ben Bacha, A., and Kotb, M. (2012). Etiology of autistic features: the persisting neurotoxic effects of propionic acid. J. Neuroinflammation 9, 74. doi:10.1186/1742-2094-9-74

PubMed Abstract | CrossRef Full Text | Google Scholar

El-Ansary, A., Bhat, R. S., Al-Daihan, S., and Al Dbass, A. M. (2015). The neurotoxic effects of ampicillin-associated gut bacterial imbalances compared to those of orally administered propionic acid in the etiology of persistent autistic features in rat pups: effects of various dietary regimens. Gut Pathog. 7, 7. doi:10.1186/s13099-015-0054-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Engeland, C. G., Bosch, J. A., and Rohleder, N. (2019). Salivary biomarkers in psychoneuroimmunology. Curr. Opin. Behav. Sci. 28, 58–65. doi:10.1016/j.cobeha.2019.01.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Farah, R., Haraty, H., Salame, Z., Fares, Y., Ojcius, D. M., and Said Sadier, N. (2018). Salivary biomarkers for the diagnosis and monitoring of neurological diseases. Biomed. J. 41, 63–87. doi:10.1016/j.bj.2018.03.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Feng, P., Zhao, S., Zhang, Y., and Li, E. (2023). A review of probiotics in the treatment of autism spectrum disorders: perspectives from the gut-brain axis. Front. Microbiol. 14, 1123462. doi:10.3389/fmicb.2023.1123462

PubMed Abstract | CrossRef Full Text | Google Scholar

Frey, A. M., Ansbro, K., Kamble, N. S., Pham, T. K., and Stafford, G. P. (2018). Characterisation and pure culture of putative health-associated oral bacterium BU063 (Tannerella sp. HOT-286) reveals presence of a potentially novel glycosylated S-layer. FEMS Microbiol. Lett. 365. doi:10.1093/femsle/fny180

PubMed Abstract | CrossRef Full Text | Google Scholar

Frye, R. E., Slattery, J., Macfabe, D. F., Allen-Vercoe, E., Parker, W., Rodakis, J., et al. (2015). Approaches to studying and manipulating the enteric microbiome to improve autism symptoms. Microb. Ecol. Health Dis. 26, 26878. doi:10.3402/mehd.v26.26878

PubMed Abstract | CrossRef Full Text | Google Scholar

Fujita, E., Dai, H., Tanabe, Y., Zhiling, Y., Yamagata, T., Miyakawa, T., et al. (2010). Autism spectrum disorder is related to endoplasmic reticulum stress induced by mutations in the synaptic cell adhesion molecule, CADM1. Cell Death Dis. 1, e47. doi:10.1038/cddis.2010.23

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghiretti, A. E., and Paradis, S. (2011). The GTPase Rem2 regulates synapse development and dendritic morphology. Dev. Neurobiol. 71, 374–389. doi:10.1002/dneu.20868

PubMed Abstract | CrossRef Full Text | Google Scholar

Gómez-Carballa, A., Navarro, L., Pardo-Seco, J., Bello, X., Pischedda, S., Viz-Lasheras, S., et al. (2023). Music compensates for altered gene expression in age-related cognitive disorders. Sci. Rep. 13, 21259. doi:10.1038/s41598-023-48094-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Gómez-Carballa, A., Navarro, L., Mallah, N., Bello, X., Pischedda, S., Viz-Lasheras, S., et al. (2025a). Music elicits different gene expression responses in the buccal cavity of age-related cognitive disorders patients and healthy controls. Front. Aging Neurosci. 23, 1622816. doi:10.3389/fnagi.2025.1622816

PubMed Abstract | CrossRef Full Text | Google Scholar

Gómez-Carballa, A., Pischedda, S., Mallah, N., Nowak, W., Martinón-Torres, F., Navarro, L., et al. (2025b). Can music modulate gene expression involved in traumatic brain injury? An integrative transcriptomic and epigenomic proof of concept. BioRxiv. doi:10.1101/2025.1107.1102.662815

CrossRef Full Text | Google Scholar

Goris, J., Brass, M., Cambier, C., Delplanque, J., Wiersema, J. R., and Braem, S. (2020). The relation between preference for predictability and autistic traits. Autism Res. 13, 1144–1154. doi:10.1002/aur.2244

PubMed Abstract | CrossRef Full Text | Google Scholar

Gosch, A., Banemann, R., Dorum, G., Haas, C., Hadrys, T., Haenggi, N., et al. (2024). Spitting in the wind? the challenges of RNA sequencing for biomarker discovery from saliva. Int. J. Leg. Med. 138, 401–412. doi:10.1007/s00414-023-03100-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Gröschl, M., Wagner, R., Rauh, M., and Dorr, H. G. (2001). Stability of salivary steroids: the influences of storage, food and dental care. Steroids 66, 737–741. doi:10.1016/s0039-128x(01)00111-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Gu, Z. (2022). Complex heatmap visualization. Imeta 1, e43. doi:10.1002/imt2.43

PubMed Abstract | CrossRef Full Text | Google Scholar

Gu, Z., and Hubschmann, D. (2023). simplifyEnrichment: a bioconductor package for clustering and visualizing functional enrichment results. Genomics Proteomics Bioinforma. 21, 190–202. doi:10.1016/j.gpb.2022.04.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Guo, D., Yang, X., and Shi, L. (2020). Rho GTPase regulators and effectors in autism spectrum disorders: animal models and insights for therapeutics. Cells 9, 835. doi:10.3390/cells9040835

PubMed Abstract | CrossRef Full Text | Google Scholar

Hänzelmann, S., Castelo, R., and Guinney, J. (2013). GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinforma. 14, 7. doi:10.1186/1471-2105-14-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Hazlett, H. C., Gu, H., Munsell, B. C., Kim, S. H., Styner, M., Wolff, J. J., et al. (2017). Early brain development in infants at high risk for autism spectrum disorder. Nature 542, 348–351. doi:10.1038/nature21369

PubMed Abstract | CrossRef Full Text | Google Scholar

Hicks, S. D., Rajan, A. T., Wagner, K. E., Barns, S., Carpenter, R. L., and Middleton, F. A. (2018). Validation of a salivary RNA test for childhood autism spectrum disorder. Front. Genet. 9, 534. doi:10.3389/fgene.2018.00534

PubMed Abstract | CrossRef Full Text | Google Scholar

Iglesias-Vazquez, L., Van Ginkel Riba, G., Arija, V., and Canals, J. (2020). Composition of gut microbiota in children with autism spectrum disorder: a systematic review and meta-analysis. Nutrients 12, 792. doi:10.3390/nu12030792

PubMed Abstract | CrossRef Full Text | Google Scholar

Israel, S., Lerer, E., Shalev, I., Uzefovsky, F., Reibold, M., Bachner-Melman, R., et al. (2008). Molecular genetic studies of the arginine vasopressin 1a receptor (AVPR1a) and the oxytocin receptor (OXTR) in human behaviour: from autism to altruism with some notes in between. Prog. Brain Res. 170, 435–449. doi:10.1016/S0079-6123(08)00434-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Jash, S., and Sharma, S. (2021). In utero immune programming of autism spectrum disorder (ASD). Hum. Immunol. 82, 379–384. doi:10.1016/j.humimm.2021.02.002

PubMed Abstract | CrossRef Full Text | Google Scholar

John, S., and Jaeggi, A. V. (2021). Oxytocin levels tend to be lower in autistic children: a meta-analysis of 31 studies. Autism 25, 2152–2161. doi:10.1177/13623613211034375

PubMed Abstract | CrossRef Full Text | Google Scholar

Kakabadze, E., Makalatia, K., Bakuradze, N., Grdzelishvili, N., Kereselidze, S., Ediberidze, T., et al. (2022). Gut health of children with autism spectrum disorder. World Acad. Sci. J. 4, 29. doi:10.3892/wasj.2022.164

CrossRef Full Text | Google Scholar

Kamble, N. S., Bera, S., Bhedase, S. A., Gaur, V., and Chowdhur, D. (2024). Review on applied applications of microbiome on human lives. Bacteria 3, 141–159. doi:10.3390/bacteria3030010

CrossRef Full Text | Google Scholar

Kamp-Becker, I. (2024). Autism spectrum disorder in ICD-11-a critical reflection of its possible impact on clinical practice and research. Mol. Psychiatry 29, 633–638. doi:10.1038/s41380-023-02354-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Kanduri, C., Raijas, P., Ahvenainen, M., Philips, A. K., Ukkola-Vuoti, L., Lahdesmaki, H., et al. (2015). The effect of listening to music on human transcriptome. PeerJ 3, e830. doi:10.7717/peerj.830

PubMed Abstract | CrossRef Full Text | Google Scholar

Kawada, K., and Mimori, S. (2018). Implication of endoplasmic reticulum stress in autism spectrum disorder. Neurochem. Res. 43, 147–152. doi:10.1007/s11064-017-2370-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Khalil, M., Azouz, H. G., Ahmed, S. A., Gad, H. A., and Omar, O. M. (2021). Sensory processing and gastrointestinal manifestations in autism spectrum disorders: no relation to clostridium difficile. J. Mol. Neurosci. 71, 153–161. doi:10.1007/s12031-020-01636-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Kielinen, M., Rantala, H., Timonen, E., Linna, S. L., and Moilanen, I. (2004). Associated medical disorders and disabilities in children with autistic disorder: a population-based study. Autism 8, 49–60. doi:10.1177/1362361304040638

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, N., Kim, K. H., Lim, W. J., Kim, J., Kim, S. A., and Yoo, H. J. (2020). Whole exome sequencing identifies novel de novo variants interacting with six gene networks in autism spectrum disorder. Genes (Basel) 12, 1. doi:10.3390/genes12010001

PubMed Abstract | CrossRef Full Text | Google Scholar

Kuhn, M., Grave, S., Bransfield, R., and Harris, S. (2012). Long term antibiotic therapy may be an effective treatment for children co-morbid with Lyme disease and autism spectrum disorder. Med. Hypotheses 78, 606–615. doi:10.1016/j.mehy.2012.01.037

PubMed Abstract | CrossRef Full Text | Google Scholar

Langfelder, P., and Horvath, S. (2008). WGCNA: an R package for weighted correlation network analysis. BMC Bioinforma. 9, 559. doi:10.1186/1471-2105-9-559

PubMed Abstract | CrossRef Full Text | Google Scholar

Lefter, R., Ciobica, A., Timofte, D., Stanciu, C., and Trifan, A. (2019). A descriptive review on the prevalence of gastrointestinal disturbances and their multiple associations in autism spectrum disorder. Med. Kaunas. 56, 11. doi:10.3390/medicina56010011

PubMed Abstract | CrossRef Full Text | Google Scholar

Lense, M. D., Beck, S., Liu, C., Pfeiffer, R., Diaz, N., Lynch, M., et al. (2020). Parents, peers, and musical play: integrated parent-child music class program supports community participation and well-being for families of children with and without autism spectrum disorder. Front. Psychol. 11, 555717. doi:10.3389/fpsyg.2020.555717

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, Q., and Zhou, J. M. (2016). The microbiota-gut-brain axis and its potential therapeutic role in autism spectrum disorder. Neuroscience 324, 131–139. doi:10.1016/j.neuroscience.2016.03.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, X., Chauhan, A., Sheikh, A. M., Patil, S., Chauhan, V., Li, X. M., et al. (2009). Elevated immune response in the brain of autistic patients. J. Neuroimmunol. 207, 111–116. doi:10.1016/j.jneuroim.2008.12.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Lord, C., Brugha, T. S., Charman, T., Cusack, J., Dumas, G., Frazier, T., et al. (2020). Autism spectrum disorder. Nat. Rev. Dis. Prim. 6, 5. doi:10.1038/s41572-019-0138-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Love, M. I., Huber, W., and Anders, S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550. doi:10.1186/s13059-014-0550-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Majewska, M. D., Hill, M., Urbanowicz, E., Rok-Bujko, P., Bienkowski, P., Namyslowska, I., et al. (2014). Marked elevation of adrenal steroids, especially androgens, in saliva of prepubertal autistic children. Eur. Child. Adolesc. Psychiatry 23, 485–498. doi:10.1007/s00787-013-0472-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Maracci, C., Motta, S., Romagnoli, A., Costantino, M., Perego, P., and Di Marino, D. (2022). The mTOR/4E-BP1/eIF4E signalling pathway as a source of cancer drug targets. Curr. Med. Chem. 29, 3501–3529. doi:10.2174/0929867329666220224112042

PubMed Abstract | CrossRef Full Text | Google Scholar

Mazurek, M. O., Lu, F., Symecko, H., Butter, E., Bing, N. M., Hundley, R. J., et al. (2017). A prospective study of the concordance of DSM-IV and DSM-5 diagnostic criteria for autism spectrum disorder. J. Autism Dev. Disord. 47, 2783–2794. doi:10.1007/s10803-017-3200-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Mcdougle, C. J., Erickson, C. A., Stigler, K. A., and Posey, D. J. (2005). Neurochemistry in the pathophysiology of autism. J. Clin. Psychiatry 66 (Suppl. 10), 9–18.

PubMed Abstract | Google Scholar

Meindl, J. N., Saba, S., Gray, M., Stuebing, L., and Jarvis, A. (2019). Reducing blood draw phobia in an adult with autism spectrum disorder using low-cost virtual reality exposure therapy. J. Appl. Res. Intellect. Disabil. 32, 1446–1452. doi:10.1111/jar.12637

PubMed Abstract | CrossRef Full Text | Google Scholar

Moretto, E., Longatti, A., Murru, L., Chamma, I., Sessa, A., Zapata, J., et al. (2019). TSPAN5 enriched microdomains provide a platform for dendritic spine maturation through neuroligin-1 clustering. Cell Rep. 29, 1130–1146.e8. doi:10.1016/j.celrep.2019.09.051

PubMed Abstract | CrossRef Full Text | Google Scholar

Mortazavi, A., Williams, B. A., Mccue, K., Schaeffer, L., and Wold, B. (2008). Mapping and quantifying mammalian transcriptomes by RNA-seq. Nat. Methods 5, 621–628. doi:10.1038/nmeth.1226

PubMed Abstract | CrossRef Full Text | Google Scholar

Munesue, T., Nakamura, H., Kikuchi, M., Miura, Y., Takeuchi, N., Anme, T., et al. (2016). Oxytocin for male subjects with autism spectrum disorder and comorbid intellectual disabilities: a randomized pilot study. Front. Psychiatry 7, 2. doi:10.3389/fpsyt.2016.00002

PubMed Abstract | CrossRef Full Text | Google Scholar

Murai, S., Saito, H., Masuda, Y., Itoh, T., and Kawaguchi, T. (1998). Sex-dependent differences in the concentrations of the principal neurotransmitters, noradrenaline and acetylcholine, in the three major salivary glands of mice. Arch. Oral Biol. 43, 9–14. doi:10.1016/s0003-9969(97)00088-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Muth, A., Honekopp, J., and Falter, C. M. (2014). Visuo-spatial performance in autism: a meta-analysis. J. Autism Dev. Disord. 44, 3245–3263. doi:10.1007/s10803-014-2188-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Nakamura, T., Matsui, M., Uchida, K., Futatsugi, A., Kusakawa, S., Matsumoto, N., et al. (2004). M(3) muscarinic acetylcholine receptor plays a critical role in parasympathetic control of salivation in mice. J. Physiol. 558, 561–575. doi:10.1113/jphysiol.2004.064626

PubMed Abstract | CrossRef Full Text | Google Scholar

Navarro, L., Martinón-Torres, F., and Salas, A. (2021). Sensogenomics and the biological background underlying musical stimuli: perspectives for a new era of musical research. Genes (Basel) 12, 1454. doi:10.3390/genes12091454

PubMed Abstract | CrossRef Full Text | Google Scholar

Navarro, L., Gómez-Carballa, A., Pischedda, S., Montoto-Louzao, J., Viz-Lasheras, S., Camino-Mera, A., et al. (2023). Sensogenomics of music and Alzheimer's disease: an interdisciplinary view from neuroscience, transcriptomics, and epigenomics. Front. Aging Neurosci. 15, 1063536. doi:10.3389/fnagi.2023.1063536

PubMed Abstract | CrossRef Full Text | Google Scholar

Navarro, L., Mallah, N., Nowak, W., Pardo-Seco, J., Gómez-Carballa, A., Pischedda, S., et al. (2025a). The effect of music interventions in autism spectrum disorder: a systematic review and meta-analysis. Front. Integr. Neurosci. 12. doi:10.3389/fnint.2025.1673618

PubMed Abstract | CrossRef Full Text | Google Scholar

Navarro, L., Mallah, N., Pardo-Seco, J., Gómez-Carballa, A., Pischedda, S., Nowak, W., et al. (2025b). Systematic review and meta-analysis reveal positive therapeutic effects of music in brain damage rehabilitation. medRxiv. doi:10.1101/2025.1107.1129.25332404

CrossRef Full Text | Google Scholar

Ostheim, P., Tichy, A., Sirak, I., Davidkova, M., Stastna, M. M., Kultova, G., et al. (2020). Overcoming challenges in human saliva gene expression measurements. Sci. Rep. 10, 11147. doi:10.1038/s41598-020-67825-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Oudshoorn, D., Van Boheemen, S., Sanchez-Aparicio, M. T., Rajsbaum, R., Garcia-Sastre, A., and Versteeg, G. A. (2012). HERC6 is the main E3 ligase for global ISG15 conjugation in mouse cells. PLoS One 7, e29870. doi:10.1371/journal.pone.0029870

PubMed Abstract | CrossRef Full Text | Google Scholar

Ounit, R., and Lonardi, S. (2016). Higher classification sensitivity of short metagenomic reads with CLARK-S. Bioinformatics 32, 3823–3825. doi:10.1093/bioinformatics/btw542

PubMed Abstract | CrossRef Full Text | Google Scholar

Ounit, R., Wanamaker, S., Close, T. J., and Lonardi, S. (2015). CLARK: fast and accurate classification of metagenomic and genomic sequences using discriminative k-mers. BMC Genomics 16, 236. doi:10.1186/s12864-015-1419-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Pagès, H., Carlson, M., Falcon, S., and Li, N. (2024). AnnotationDbi: manipulation of SQLite-based annotations in bioconductor. Available online at: https://bioconductor.org/packages/devel/bioc/manuals/AnnotationDbi/man/AnnotationDbi.pdf.

Google Scholar

Parker, K. J., Oztan, O., Libove, R. A., Sumiyoshi, R. D., Jackson, L. P., Karhson, D. S., et al. (2017). Intranasal oxytocin treatment for social deficits and biomarkers of response in children with autism. Proc. Natl. Acad. Sci. U. S. A. 114, 8119–8124. doi:10.1073/pnas.1705521114

PubMed Abstract | CrossRef Full Text | Google Scholar

Peralta-Marzal, L. N., Prince, N., Bajic, D., Roussin, L., Naudon, L., Rabot, S., et al. (2021). The impact of gut microbiota-derived metabolites in autism spectrum disorders. Int. J. Mol. Sci. 22, 10052. doi:10.3390/ijms221810052

PubMed Abstract | CrossRef Full Text | Google Scholar

Peralta-Marzal, L. N., Rojas-Velazquez, D., Rigters, D., Prince, N., Garssen, J., Kraneveld, A. D., et al. (2024). A robust microbiome signature for autism spectrum disorder across different studies using machine learning. Sci. Rep. 14, 814. doi:10.1038/s41598-023-50601-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Pfaffe, T., Cooper-White, J., Beyerlein, P., Kostner, K., and Punyadeera, C. (2011). Diagnostic potential of saliva: current state and future applications. Clin. Chem. 57, 675–687. doi:10.1373/clinchem.2010.153767

PubMed Abstract | CrossRef Full Text | Google Scholar

Popoff, M. R. (2011). Multifaceted interactions of bacterial toxins with the gastrointestinal mucosa. Future Microbiol. 6, 763–797. doi:10.2217/fmb.11.58

PubMed Abstract | CrossRef Full Text | Google Scholar

Proctor, G. B., and Carpenter, G. H. (2007). Regulation of salivary gland function by autonomic nerves. Auton. Neurosci. 133, 3–18. doi:10.1016/j.autneu.2006.10.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Putri, G. H., Anders, S., Pyl, P. T., Pimanda, J. E., and Zanini, F. (2022). Analysing high-throughput sequencing data in Python with HTSeq 2.0. Bioinformatics 38, 2943–2945. doi:10.1093/bioinformatics/btac166

PubMed Abstract | CrossRef Full Text | Google Scholar

Ragusa, M., Santagati, M., Mirabella, F., Lauretta, G., Cirnigliaro, M., Brex, D., et al. (2020). Potential associations among alteration of salivary miRNAs, saliva microbiome structure, and cognitive impairments in autistic children. Int. J. Mol. Sci. 21, 6203. doi:10.3390/ijms21176203

PubMed Abstract | CrossRef Full Text | Google Scholar

Rawsthorne, H., Calahorro, F., Feist, E., Holden-Dye, L., O'connor, V., and Dillon, J. (2021). Neuroligin dependence of social behaviour in Caenorhabditis elegans provides a model to investigate an autism-associated gene. Hum. Mol. Genet. 29, 3546–3553. doi:10.1093/hmg/ddaa232

PubMed Abstract | CrossRef Full Text | Google Scholar

Remington, A., and Fairnie, J. (2017). A sound advantage: increased auditory capacity in autism. Cognition 166, 459–465. doi:10.1016/j.cognition.2017.04.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Ricci, S., Businaro, R., Ippoliti, F., Lo Vasco, V. R., Massoni, F., Onofri, E., et al. (2013). Altered cytokine and BDNF levels in autism spectrum disorder. Neurotox. Res. 24, 491–501. doi:10.1007/s12640-013-9393-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Ritchie, M. E., Phipson, B., Wu, D., Hu, Y., Law, C. W., Shi, W., et al. (2015). Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43, e47. doi:10.1093/nar/gkv007

PubMed Abstract | CrossRef Full Text | Google Scholar

Robinson-Agramonte, M. L. A., Noris Garcia, E., Fraga Guerra, J., Vega Hurtado, Y., Antonucci, N., Semprun-Hernandez, N., et al. (2022). Immune dysregulation in autism spectrum disorder: what do we know about it? Int. J. Mol. Sci. 23, 3033. doi:10.3390/ijms23063033

PubMed Abstract | CrossRef Full Text | Google Scholar

Rubenstein, J. L., and Merzenich, M. M. (2003). Model of autism: increased ratio of excitation/inhibition in key neural systems. Genes Brain Behav. 2, 255–267. doi:10.1034/j.1601-183x.2003.00037.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Rudis, B., and Gandy, D. (2015). Waffle: create waffle chart visualizations. In CRAN: contributed packages. doi:10.32614/CRAN.package.waffle

CrossRef Full Text | Google Scholar

Rylander-Rudqvist, T., Hakansson, N., Tybring, G., and Wolk, A. (2006). Quality and quantity of saliva DNA obtained from the self-administrated oragene method--a pilot study on the cohort of Swedish men. Cancer Epidemiol. Biomarkers Prev. 15, 1742–1745. doi:10.1158/1055-9965.EPI-05-0706

PubMed Abstract | CrossRef Full Text | Google Scholar

Saghazadeh, A., Ataeinia, B., Keynejad, K., Abdolalizadeh, A., Hirbod-Mobarakeh, A., and Rezaei, N. (2019). A meta-analysis of pro-inflammatory cytokines in autism spectrum disorders: effects of age, gender, and latitude. J. Psychiatr. Res. 115, 90–102. doi:10.1016/j.jpsychires.2019.05.019

PubMed Abstract | CrossRef Full Text | Google Scholar

Salas, A., Navarro, L., Martinón-Torres, F., and Gómez-Carballa, A. (2025). Beyond behavioral studies: exploring the multi-omics impact of music on human health. Phys. Life Rev. doi:10.1016/j.plrev.2025.1009.1002

CrossRef Full Text | Google Scholar

Samson, F., Mottron, L., Soulieres, I., and Zeffiro, T. A. (2012). Enhanced visual functioning in autism: an ALE meta-analysis. Hum. Brain Mapp. 33, 1553–1581. doi:10.1002/hbm.21307

PubMed Abstract | CrossRef Full Text | Google Scholar

Sanchez-Tena, S., Cubillos-Rojas, M., Schneider, T., and Rosa, J. L. (2016). Functional and pathological relevance of HERC family proteins: a decade later. Cell Mol. Life Sci. 73, 1955–1968. doi:10.1007/s00018-016-2139-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Sandbank, M., Bottema-Beutel, K., Crowley, S., Cassidy, M., Dunham, K., Feldman, J. I., et al. (2020). Project AIM: autism intervention meta-analysis for studies of young children. Psychol. Bull. 146, 1–29. doi:10.1037/bul0000215

PubMed Abstract | CrossRef Full Text | Google Scholar

Sansores-Espana, L. D., Melgar-Rodriguez, S., Olivares-Sagredo, K., Cafferata, E. A., Martinez-Aguilar, V. M., Vernal, R., et al. (2021). Oral-gut-brain axis in experimental models of periodontitis: associating gut dysbiosis with neurodegenerative diseases. Front. Aging 2, 781582. doi:10.3389/fragi.2021.781582

PubMed Abstract | CrossRef Full Text | Google Scholar

Sato-Suzuki, Y., Washio, J., Wicaksono, D. P., Sato, T., Fukumoto, S., and Takahashi, N. (2020). Nitrite-producing oral microbiome in adults and children. Sci. Rep. 10, 16652. doi:10.1038/s41598-020-73479-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Seltzer, M. M., Shattuck, P., Abbeduto, L., and Greenberg, J. S. (2004). Trajectory of development in adolescents and adults with autism. Ment. Retard. Dev. Disabil. Res. Rev. 10, 234–247. doi:10.1002/mrdd.20038

PubMed Abstract | CrossRef Full Text | Google Scholar

Sharda, M., Tuerk, C., Chowdhury, R., Jamey, K., Foster, N., Custo-Blanch, M., et al. (2018). Music improves social communication and auditory-motor connectivity in children with autism. Transl. Psychiatry 8, 231. doi:10.1038/s41398-018-0287-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Shin, S. M., Zhang, N., Hansen, J., Gerges, N. Z., Pak, D. T., Sheng, M., et al. (2012). GKAP orchestrates activity-dependent postsynaptic protein remodeling and homeostatic scaling. Nat. Neurosci. 15, 1655–1666. doi:10.1038/nn.3259

PubMed Abstract | CrossRef Full Text | Google Scholar

Soke, G. N., Maenner, M. J., Christensen, D., Kurzius-Spencer, M., and Schieve, L. A. (2018). Prevalence of Co-occurring medical and behavioral conditions/symptoms among 4- and 8-Year-Old children with autism spectrum disorder in selected areas of the United States in 2010. J. Autism Dev. Disord. 48, 2663–2676. doi:10.1007/s10803-018-3521-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Spielmann, N., Ilsley, D., Gu, J., Lea, K., Brockman, J., Heater, S., et al. (2012). The human salivary RNA transcriptome revealed by massively parallel sequencing. Clin. Chem. 58, 1314–1321. doi:10.1373/clinchem.2011.176941

PubMed Abstract | CrossRef Full Text | Google Scholar

Spratt, E. G., Nicholas, J. S., Brady, K. T., Carpenter, L. A., Hatcher, C. R., Meekins, K. A., et al. (2012). Enhanced cortisol response to stress in children in autism. J. Autism Dev. Disord. 42, 75–81. doi:10.1007/s10803-011-1214-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Srikantha, P., and Mohajeri, M. H. (2019). The possible role of the microbiota-gut-brain-axis in autism spectrum disorder. Int. J. Mol. Sci. 20, 2115. doi:10.3390/ijms20092115

PubMed Abstract | CrossRef Full Text | Google Scholar

Stoodley, C. J., D'mello, A. M., Ellegood, J., Jakkamsetti, V., Liu, P., Nebel, M. B., et al. (2017). Altered cerebellar connectivity in autism and cerebellar-mediated rescue of autism-related behaviors in mice. Nat. Neurosci. 20, 1744–1751. doi:10.1038/s41593-017-0004-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Sukhodolsky, D. G., Scahill, L., Gadow, K. D., Arnold, L. E., Aman, M. G., Mcdougle, C. J., et al. (2008). Parent-rated anxiety symptoms in children with pervasive developmental disorders: frequency and association with core autism symptoms and cognitive functioning. J. Abnorm Child. Psychol. 36, 117–128. doi:10.1007/s10802-007-9165-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Sun, C. Y., Van Koningsbruggen, S., Long, S. W., Straasheijm, K., Klooster, R., Jones, T. I., et al. (2011). Facioscapulohumeral muscular dystrophy region gene 1 is a dynamic RNA-associated and actin-bundling protein. J. Mol. Biol. 411, 397–416. doi:10.1016/j.jmb.2011.06.014

PubMed Abstract | CrossRef Full Text | Google Scholar

Talantseva, O. I., Romanova, R. S., Shurdova, E. M., Dolgorukova, T. A., Sologub, P. S., Titova, O. S., et al. (2023). The global prevalence of autism spectrum disorder: a three-level meta-analysis. Front. Psychiatry 14, 1071181. doi:10.3389/fpsyt.2023.1071181

PubMed Abstract | CrossRef Full Text | Google Scholar

Taylor, J. L., and Corbett, B. A. (2014). A review of rhythm and responsiveness of cortisol in individuals with autism spectrum disorders. Psychoneuroendocrinology 49, 207–228. doi:10.1016/j.psyneuen.2014.07.015

PubMed Abstract | CrossRef Full Text | Google Scholar

The American Music Therapy Association (2014). Music therapy fact sheets and bibliographies. Available online at: https://www.musictherapy.org/research/factsheets/.

Google Scholar

Ukkola-Vuoti, L., Kanduri, C., Oikkonen, J., Buck, G., Blancher, C., Raijas, P., et al. (2013). Genome-wide copy number variation analysis in extended families and unrelated individuals characterized for musical aptitude and creativity in music. PLoS One 8, e56356. doi:10.1371/journal.pone.0056356

PubMed Abstract | CrossRef Full Text | Google Scholar

Voineagu, I., Wang, X., Johnston, P., Lowe, J. K., Tian, Y., Horvath, S., et al. (2011). Transcriptomic analysis of autistic brain reveals convergent molecular pathology. Nature 474, 380–384. doi:10.1038/nature10110

PubMed Abstract | CrossRef Full Text | Google Scholar

Volkmar, F. R., Cohen, D. J., Bregman, J. D., Hooks, M. Y., and Stevenson, J. M. (1989). An examination of social typologies in autism. J. Am. Acad. Child. Adolesc. Psychiatry 28, 82–86. doi:10.1097/00004583-198901000-00015

PubMed Abstract | CrossRef Full Text | Google Scholar

Volkmar, F., Siegel, M., Woodbury-Smith, M., King, B., Mccracken, J., State, M., et al. (2014). Practice parameter for the assessment and treatment of children and adolescents with autism spectrum disorder. J. Am. Acad. Child. Adolesc. Psychiatry 53, 237–257. doi:10.1016/j.jaac.2013.10.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Wei, H., Zou, H., Sheikh, A. M., Malik, M., Dobkin, C., Brown, W. T., et al. (2011). IL-6 is increased in the cerebellum of autistic brain and alters neural cell adhesion, migration and synaptic formation. J. Neuroinflammation 8, 52. doi:10.1186/1742-2094-8-52

PubMed Abstract | CrossRef Full Text | Google Scholar

Wickham, H. (2016). ggplot2: elegant graphics for data analysis. New York: Springer-Verlag.

Google Scholar

Wolff, J. J., and Symons, F. J. (2013). An evaluation of multi-component exposure treatment of needle phobia in an adult with autism and intellectual disability. J. Appl. Res. Intellect. Disabil. 26, 344–348. doi:10.1111/jar.12002

PubMed Abstract | CrossRef Full Text | Google Scholar

Wong, G. C., Montgomery, J. M., and Taylor, M. W. (2021). The gut-microbiota-brain axis in autism spectrum disorder. Brisbane (AU): Exon Publications.

CrossRef Full Text | Google Scholar

Wong, A. S. K., Burns, S., and Woodruff, E. (2022). Examining the impact of social stressor stimuli in eliciting physiological reactivity in children and adolescents with autism spectrum disorder: a systematic review and meta-analysis protocol. BMJ Open 12, e060048. doi:10.1136/bmjopen-2021-060048

PubMed Abstract | CrossRef Full Text | Google Scholar

World Health Organization (2021). International classification of diseases, eleventh revision (ICD-11).

Google Scholar

Wu, T., Hu, E., Xu, S., Chen, M., Guo, P., Dai, Z., et al. (2021). clusterProfiler 4.0: a universal enrichment tool for interpreting omics data. Innov. (N Y) 2, 100141. doi:10.1016/j.xinn.2021.100141

PubMed Abstract | CrossRef Full Text | Google Scholar

Yatawara, C. J., Einfeld, S. L., Hickie, I. B., Davenport, T. A., and Guastella, A. J. (2016). The effect of oxytocin nasal spray on social interaction deficits observed in young children with autism: a randomized clinical crossover trial. Mol. Psychiatry 21, 1225–1231. doi:10.1038/mp.2015.162

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, G. (2024). Thirteen years of clusterProfiler. Innov. (Camb) 5, 100722. doi:10.1016/j.xinn.2024.100722

PubMed Abstract | CrossRef Full Text | Google Scholar

Zatorre, R., and Mcgill, J. (2005). Music, the food of neuroscience? Nature 434, 312–315. doi:10.1038/434312a

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeidan, J., Fombonne, E., Scorah, J., Ibrahim, A., Durkin, M. S., Saxena, S., et al. (2022). Global prevalence of autism: a systematic review update. Autism Res. 15, 778–790. doi:10.1002/aur.2696

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeidan-Chulia, F., De Oliveira, B. N., Casanova, M. F., Casanova, E. L., Noda, M., Salmina, A. B., et al. (2016). Up-regulation of oligodendrocyte lineage markers in the cerebellum of autistic patients: evidence from network analysis of gene expression. Mol. Neurobiol. 53, 4019–4025. doi:10.1007/s12035-015-9351-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, H., Ben Zablah, Y., Zhang, H., and Jia, Z. (2021). Rho signaling in synaptic plasticity, memory, and brain disorders. Front. Cell Dev. Biol. 9, 729076. doi:10.3389/fcell.2021.729076

PubMed Abstract | CrossRef Full Text | Google Scholar

Zoghbi, H. Y. (2003). Postnatal neurodevelopmental disorders: meeting at the synapse? Science 302, 826–830. doi:10.1126/science.1089071

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: ASD, microbiome, meta-genomics, transcriptomics, gene expression, RNA-seq, sensogenomics

Citation: Cavenaghi A, Mallah NEZ, Navarro L, Martinón-Torres F, Gómez-Carballa A and Salas A (2025) Decoding the peripheral transcriptomic and meta-genomic response to music in autism spectrum disorder via saliva-based RNA sequencing. Front. Mol. Biosci. 12:1696704. doi: 10.3389/fmolb.2025.1696704

Received: 04 September 2025; Accepted: 21 October 2025;
Published: 24 November 2025.

Edited by:

Hem Chandra Jha, Indian Institute of Technology Indore, India

Reviewed by:

Nitin Kamble, University of Cincinnati Medical Center, United States
Wu Yuqi, Microbiota I-Center Limited, Hong Kong SAR, China

Copyright © 2025 Cavenaghi, Mallah, Navarro, Martinón-Torres, Gómez-Carballa and Salas. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Antonio Salas, YW50b25pby5zYWxhc0B1c2MuZXM=

These authors have contributed equally to this work

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.