The role of the microbiome on social interactions in mice

The original article by Stilling et al (2018), described the association between gut microbiome composition and social behaviour in mice. They compared the social interaction of germ-free mice (GF, devoid of a microbiome) against control mice (CON) and GF mice which had their microbiome reconstituted (exGF). All three groups had naive controls and a group of mice that were exposed to social interactions (CON-SI, GF-SI and exGF-SI). They observed an attenuation of the dynamic response to social interactions in GF and exGF mice, accompanied by an increased expression of splicing factors and alternative exon usage in the amygdala. They concluded that microbiome composition could affect neurological development and might be important in conditions such as autism.

For this article, I downloaded the RNAseq data generated for the experiment (GEO code: GSE114702) and injected it into the Omics Playground platform in order to replicate the main insights gained by the authors.

Clustering and Differential Expression Analysis

Partial clustering according to whether mice had been exposed to social interactions (social) or not (naive) could be observed using UMAP. However, no clear separation by microbiome and treatment state could be observed (Fig. 1).

**Figure 1**. The UMAP clustering plot shows a partial clustering of samples by exposure to social interactions (blue) or naive state (green). However, no clustering was apparent when considering the microbiome state. Symbol keys: ✫: control samples, ■: ex-GF samples, ●: GF samples.

Next, I performed a pairwise gene expression analysis, as shown in Table 1. Unlike the original article, I used the intersection of three differential expression analysis algorithms (Deseq2, EdgeR and limma) to identify significantly dysregulated genes (FDR<0.1, logFC>|±0.5|), which yielded a lower number of hits. The GF-SI vs GF comparison had the highest number of dysregulated genes, followed by the CON-SI vs CON comparison (Fig. 2), consistent with the findings obtained by the authors.

Samples Group	Reference Samples Group
GF	CON
exGF	CON
GF	exGF
CON-SI	CON
GF-SI	GF
exGF-SI	exGF
GF-SI	CON-SI
exGF-SI	CON-SI
GF-SI	exGF-SI

Table 1. Pairwise comparisons across experimental groups.

**Figure 2**. Volcano plots for all the available pairwise comparisons. The plots are based on the meta q values obtained from the three DE algorithms used in this analysis (DESEQ2, EdgeR, limma). The GF-SI vs GF and CON-SI vs CON comparisons yielded the greatest number of DEGs.

An intersection of of the comparisons involving GF and CON groups (GF-SI vs GF, CON-SI vs CON, GF vs CON and GF-SI vs CON-SI) revealed a number of genes that were activated or repressed upon exposure to social stimuli (51 UP and 17 DOWN), regardless of microbiome state. However, there were significant differences in dysregulated genes between control (101 UP and 103 DOWN) and germ-free comparisons (94 UP and 207 DOWN) (Fig. 3).

**Figure 3**. Venn diagram showing the intersection between the pairwise comparisons involving GF and CON groups. There is a significant number of dysregulated genes shared between the CON-SI vs CON and GF-SI vs GF comparisons. This indicates similar reactions in response to social stimulation. However, an even larger number of dysregulated genes is unique to either phenotype, with the GF-SI vs GF comparisons displaying the largest number of uniquely dysregulated genes. Labels: A) GF vs CON, B) CON-SI vs CON, C) GF-SI vs GF, D) GF-SI vs CON-SI. Upper numbers indicate upregulated genes, lower numbers indicate downregulated genes.

Gene Set Enrichment Analysis: KEGG and GO

Analysis of GO terms and of KEGG pathways confirmed the main findings of the authors (Fig. 4). There was an upregulation of pathways associated with the splicesome machinery, regulation of gene expression and protein folding in GF mice after exposure to social interactions. I also observed an upregulation of the MAPK pathway in control mice following social stimulation (Fig. 5).

**Figure 4**. KEGG pathway and GO terms geneset analysis indicates that pathways unique to GF mice are activated upon exposure to social stimuli. A) KEGG activation matrix, B) GO terms activation matrix.

**Figure 5**. Visualisation of the MAPK KEGG Pathway for the CON-SI vs CON pairwise comparison. The MAPK pathway was flagged as upregulated.

Furthermore, by using the geneset enrichment module in the platform, I could rapidly generate GSEA plots for the “Splicesome” KEGG pathway and the ”RNA splicing” GO term, replicating the results from the original manuscript (Fig. 6A-B).

**Figure 6**. GSEA plots for the spliceosome KEGG pathway (A) and the RNA splicing GO term (B).

Gene Set Enrichment Analysis: L1000 database

The gene perturbation L1000 database available on the platform can be used to explore how the expression changes in a single gene correlate with an experimental gene expression profile. I used it to identify if the expression or downregulation of an individual gene could be used to understand more about the MOA of microbiome removal on social interactions using the GF-SI vs CON-SI pairwise comparison.

Repression of transcription factor ETS2 (ERF) was found to mimic the differences in gene expression observed between GF-SI and CON-SI (Fig. 7A). ERF may be needed for correct neuronal function (Maroulakou and Bowe, 2001) and downregulation has been observed in subjects with autism spectrum disorder (Burraco et al, 2016). This appears consistent with the information provided by the platform.

Figure 7. GSEA plots obtained by comparing the GF-SI vs CON-SI signature against the L1000 individual gene perturbation database. (A) ERF gene silencing, (B) EDN3 expression inhibition.

Conversely, repression of EDN3 (endothelin 3) expression was found to generate an opposite profile to the one observed for the GF-SI vs CON-SI pairwise comparison (Fig. 7B). Endothelins act as neurotransmitters in the central nervous system and EDN3 has been observed to be upregulated in people suffering from autism spectrum disorder (Burraco et al, 2016). This result would thus suggest that upregulation of EDN3 may play a role in impaired social behaviour in GF-SI mice.

Finally, I explored drugs that could inhibit the gene expression profile generated by stimulated GF mice against stimulated CON mice. The main mode of action (MOA) among inhibitory drugs was mTOR pathway inhibition (Fig. 8). Interestingly, suppression of the AKT/mTOR pathway has been shown to rescue social behaviour in genetically modified mice (Xing et al, 2019), while increased mTOR pathway activation has been associated with autism (Rosina et al, 2019). The list also includes PI3K inhibitors, with the PI3K pathway strongly intertwined with the AKT/mTOR pathway. The second most common family of inhibitors, glycogen synthase kinase (GSK) inhibitors, is also interesting, since GSK3 overexpression has been shown to influence social preference and induce anxiety during social interactions (Mines et al, 2010) and GSK3 inhibition is recognised as a treatment for psychiatric disorders (Finkelman and Martinez, 2011).

**Figure 8**. Most common MOA for both matching and inhibitory compounds from the L1000 database. Highlighted in green are the mTOR and PI3K inhibitors MOAs, highlighted in orange is the glycogen synthase kinase inhibitor MOA.

Conclusion

I could recreate several of the original figures and findings with the Omics Playground platform and expanded on the original analysis to include comparisons with gene expression profiles from the L1000 database. This additional analysis indicated that alterations in expression of the ERF2 and EDN3 genes could mimic the profile obtained by comparing GF-SI vs. CON-SI samples. Both genes are involved in neuronal function and previous studies have observed alterations in the expression of these genes in autistic individuals. Furthermore, mTOR and GSK inhibitors were identified as potential treatments to restore normal social functions and their therapeutic potential was supported by scientific literature.