Dynamical Electrical Complexity Is Reduced during Neuronal Differentiation in Autism Spectrum Disorder

Abstract

Neuronal activity can be modeled as a nonlinear dynamical system to yield measures of neuronal state and dysfunction. The electrical recordings of stem cell-derived neurons from individuals with autism spectrum disorder (ASD) and controls were analyzed using minimum embedding dimension (MED) analysis to characterize their dynamical complexity. MED analysis revealed a significant reduction in dynamical complexity in ASD neurons during differentiation, which was correlated to bursting and spike interval measures. MED was associated with clinical endpoints, such as nonverbal intelligence, and was correlated with 53 differentially expressed genes, which were overrepresented with ASD risk genes related to neurodevelopment, cell morphology, and cell migration. Spatiotemporal analysis also showed a prenatal temporal enrichment in cortical and deep brain structures. Together, we present dynamical analysis as a paradigm that can be used to distinguish disease-associated cellular electrophysiological and transcriptional signatures, while taking into account patient variability in neuropsychiatric disorders.

Keywords: autism spectrum disorder; dynamical complexity; minimum embedding dimension; multielectrode array; neurodevelopmental disorder models.

Graphical abstract

Figure 1

Study Design and Visualizing Dynamical Complexity (A) Conceptual overview of the study design. iPSCs from seven neurotypical controls and eight ASD cases were differentiated into neurons. Neurons were characterized through RNA sequencing and MEA recordings to identify disease signatures that integrate across both levels of analysis. (B) Dynamical complexity may vary even when standard spiking variables are constant. For these pedagogical examples, the firing rate is constant at 20 spiking events in the interval. Nevertheless, differences in the organization of the spiking events may result in large differences in complexity, as measured by dynamical analysis. These differences may distinguish groups of interest, such as ASD, and correlate to gene expression changes related to neuronal dysfunction. See also Figure S1.

Figure 2

Electrical Analysis of Neuronal Lines, Minimum Embedding Dimension, and Clinical Correlations (A) Raw MEA data. Raw electrode waveforms (top) are used for spike detection (bottom) and downstream analyses. (B) Spike interval and bursting variables highlight electrophysiological differences. Group average values with 95% confidence intervals are given for several time points of two variables. Variation in the interspike interval and number of bursting electrodes does show significant group-based differences for day 15 and 17 between control (blue) and ASD (red). For each relevant panel, significance was as tested with Welch's two-sided t test and indicated by asterisks. ^⁎p ≤ 5 × 10⁻², ^⁎⁎p ≤ 5 × 10⁻³, ^⁎⁎⁎p ≤ 5 × 10⁻⁵, and ^⁎⁎⁎⁎p ≤ 5 × 10⁻⁷. (C) Overview of group-wise differences across spiking measures and MED. This binary matrix indicates at which time points a measure was able to detect a significant difference between cases and controls (black squares). The MED outperforms all measures in distinguishing groups, but it does overlap with some bursting and spike interval variables, as shown in (B). Relatively common measures, such as the mean firing rate, fail to distinguish control and ASD activity over any day of the recording period. (D) MED offers a sustained and statistically significant separation of groups. Average MED values and 95% confidence intervals for both groups are depicted for the first eight time points. The control subjects (blue) are associated with higher MED complexity score in comparison with the ASD subjects (red). The dotted rectangle indicates the values of the MED during day 15 of the MEA recordings, when RNA sequencing was performed. (E) MED is correlated to clinical endpoints of interest. Low MED scores, which are associated with ASD diagnosis, are correlated to lower Vineland adaptive behavior score (left) and lower early nonverbal IQ (right). The gray shading represents a 95% confidence interval around the fitted curve, as estimated by a linear model. See also Figure S2.

Figure 3

Gene Expression Signatures of MED and ASD (A) MED is associated with differential expression in 1,423 genes before multiple comparison correction and 53 genes after correction, as depicted in the heatmaps. Hierarchical clustering of the samples based on gene expression separates low- and high-complexity samples, as indicated by the MED color bar. (B) Examination of ASD genes in a combined MED and ASD interaction model identifies 761 differentially expressed genes before correction and 7 differentially expressed genes after multiple comparison correction. (C) The log₂ fold changes of the 53 differentially expressed genes after multiple comparison correction for MED are shown. The red bars correspond to genes that are negatively correlated to MED, and the blue bars correspond to genes that are positively correlated to MED. Underlined genes are also implicated in ASD risk or brain function in previous studies. Most differentially expressed genes are negatively correlated to MED, and positively correlated genes have less variability in their log fold changes. See also Figure S4.

Figure 4

Enrichment of MED Gene Expression Signature with ASD-Related Genes, Biological Pathways, and Spatiotemporal Analysis To interpret the broad biological relevance of MED and ASD expression signatures, differentially expressed genes before FDR correction were further analyzed. (A) The MED (black) and ASD (gray) differentially expressed genes were tested for enrichment in five lists: (1) highly brain-expressed genes, (2) putative ASD risk genes, (3) gold standard ASD risk genes, and (4) highly liver-expressed genes (negative control). Significance was tested using the binomial test with a cutoff of p = 0.05 (red line). Both ASD and MED genes were enriched in brain, Princeton, and SFARI lists, although the effect size of the MED signature was notably stronger. Neither set overlapped with liver-expressed genes, a negative control. (B) Biological processes implicated with MED. The gene expression signature of MED affects neurodevelopmental, cell migration, and growth-related ontologies. Significance of enrichment was calculated using Fisher's exact test with false discovery rate control and a cutoff of p_FDR = 0.05 (red line). (C) Gene lists from 16 neuroanatomical regions and 13 developmental stages were tested for enrichment with the MED genes. The neuroanatomical regions include the inferior temporal cortex (ITC), thalamus (THA), dorsal frontal cortex (DFC), medial frontal cortex (MFC), cerebellar cortex (CB), ventral frontal cortex (VFC), hippocampus (HIP), superior temporal cortex (STC), primary visual cortex (V1C), amygdala (AMY), inferior parietal cortex (IPC), primary auditory cortex (A1C), striatum (STR), primary motor cortex (M1C), olfactory cortex (OFC), and primary somatosensory cortex (S1C). Fisher's exact test was performed for each region-stage pair and plotted on a heatmap after false discovery rate correction. A variety of cortical and deeper structures are implicated in MED, primarily at prenatal time points. Strong enrichment of the DFC and MFC is also indicated during late childhood. The grayed-out regions represent structures that are not present in the early prenatal brain. (D) Visualization of cortical and interior MED-associated regions during late prenatal development (25–28 post conceptual weeks). Raw enrichment p values (p_raw) were plotted to show the range of structural involvement at this time point. Cortical regions such as the A1C, S1C, and ITC are enriched for the MED signature; however, deeper lying structures such as the thalamus and striatum are also revealed.