Dynamical features in fetal and postnatal zinc-copper metabolic cycles predict the emergence of autism spectrum disorder.

Metals are critical to neurodevelopment, and dysregulation in early life has been documented in autism spectrum disorder (ASD). However, underlying mechanisms and biochemical assays to distinguish ASD cases from controls remain elusive. In a nationwide study of twins in Sweden, we tested whether zinc-copper cycles, which regulate metal metabolism, are disrupted in ASD. Using novel tooth-matrix biomarkers that provide direct measures of fetal elemental uptake, we developed a predictive model to distinguish participants who would be diagnosed with ASD in childhood from those who did not develop the disorder. We replicated our findings in three independent studies in the United States and the UK. We show that three quantifiable characteristics of fetal and postnatal zinc-copper rhythmicity are altered in ASD: the average duration of zinc-copper cycles, regularity with which the cycles recur, and the number of complex features within a cycle. In all independent study sets and in the pooled analysis, zinc-copper rhythmicity was disrupted in ASD cases. In contrast to controls, in ASD cases, the cycle duration was shorter (F = 52.25, P < 0.001), regularity was reduced (F = 47.99, P < 0.001), and complexity diminished (F = 57.30, P < 0.001). With two distinct classification models that used metal rhythmicity data, we achieved 90% accuracy in classifying cases and controls, with sensitivity to ASD diagnosis ranging from 85 to 100% and specificity ranging from 90 to 100%. These findings suggest that altered zinc-copper rhythmicity precedes the emergence of ASD, and quantitative biochemical measures of metal rhythmicity distinguish ASD cases from controls.

INTRODUCTION

Autism spectrum disorder (ASD) affects 1 to 2% of the population in developed countries. Reports have shown that multiple nutrient elements and toxic metals are differentially absorbed and metabolized in children with ASD (–). However, the mechanisms underlying elemental dysregulation and ASD risk are not well understood, and it is not known whether this elemental dysregulation is present in fetal and early postnatal life before first clinical symptoms manifest.

Here, we tested whether zinc-copper pathways were disrupted in ASD and whether quantitative measures of zinc-copper dynamics could provide a biochemical basis to distinguish ASD cases from controls. We selected zinc-copper cycles as our primary target for testing the dysregulation hypothesis because these are highly conserved in evolution, are essential for maintenance of health, and also exert homeostatic control over known neurotoxicants (). Notable examples of disorders arising from perturbations of the metabolism of these metals include Wilson’s disease, a key component of which is copper dysregulation that leads to progressive neurological and cognitive dysfunction and psychosis-like symptoms also documented in adults with ASD (). In Wilson’s disease, realigning the zinc-copper interaction via zinc supplementation increases the expression of the metal-binding protein metallothionein, which helps regulate copper levels (). Recently, zinc-associated pathways have been implicated in ASD (, , ). For example, mutations in the gene coding for SHANK3, which is part of a zinc-mediated signaling system, are linked to increased ASD risk ().

Many physiologic processes follow cyclic patterns or rhythms that operate at a wide range of time intervals, from milliseconds (for example, neuron firing) to hours (for example, body temperature), days or circadian (for example, sleep cycles), weeks, or longer (for example, menstrual cycles) (, ). For any system in the body, multiple cycles operating at different timescales may interact, and the study of these processes requires methods that focus not solely on concentration but rather on the dynamic changes over time (, ). Zinc and copper levels are not stationary in the human body, exhibiting a cyclic pattern (–). We therefore focused on their cyclic properties rather than single-point measures of serum concentrations. We use the term cycles to indicate time-dependent rhythmic variation in the concentration of elements. Our methodology samples elemental concentrations in incremental zones of teeth (akin to growth rings in trees) at approximately 2- to 3-day intervals, which allows for the identification of cyclical processes varying at an approximately 7- to 10-day period. Contrasting these cycles between ASD cases and controls offers a way to uncover system-wide dysregulations in zinc-copper rhythmicity. Furthermore, because zinc and copper pathways are central regulators of multiple metals, their disruption would have downstream effects, including incomplete metabolism of other essential elements and toxic metals (, ). Therefore, we also studied the metabolic cycles of other metals during the fetal and early childhood periods (from the second trimester to approximately 1 year of age) (, ). To test our hypothesis, we applied the tooth biomarkers in four case-control samples: a discovery population of twins from Sweden and replication in two U.S. case-control samples and a birth cohort from the UK.

RESULTS

Dysregulation of zinc-copper cycles in ASD

Our primary hypothesis predicts dysregulation of zinc-copper cycles in ASD, as quantified by cycle duration, complexity, and regularity. The mean diagonal length (MDL) measures the duration of a cycle and was used in our study to detect fragmentation of zinc-copper cycles (see Fig. 1 and Supplementary Materials and Methods). We found that MDL for zinc-copper cycles was significantly reduced across all cohorts (F = 52.25, P < 0.001), reflecting an average reduction of 15% in ASD twins (Fig. 2A). We also detected a significant ASD × study interaction (F = 3.83, P = 0.014), reflecting some across-study variability in the magnitude of ASD-related reductions, but post hoc tests confirmed significantly reduced MDL in ASD cases in all populations studied (Sweden, F = 38.05, P < 0.001; New York, USA, F = 8.07, P = 0.029; Texas, USA, F = 12.41, P = 0.005; UK, F = 8.57, P = 0.023). For additional details, see table S1.

Fig. 1

Overview of study design.

(A) Study participants recruited from Sweden (Swe) (twin–co-twin discovery set) and from the United States and UK (case-control replication sets). (B) Collected deciduous teeth were analyzed using laser ablation inductively coupled plasma mass spectrometry (LA-ICPMS) to generate temporal profiles of metal uptake during fetal and postnatal development. (C) Example exposure profile in one subject (top) ranging from −89 to 300 days since birth (TSB; time since birth in days). Dashed line indicates birth, while black arrows indicate a period of approximately 10 days. Bottom: Magnified region from −39 to 56 TSB to highlight cycles in elemental concentration varying on a roughly 10-day period. (D) Top: Example elemental exposure profiles in a single subject for two elements (Zn, blue line; Cu, green line) simultaneously sampled and overlaid from −126 to 290 TSB. Bottom: Magnified region from −126 to −64 TSB showing concentration of both elements rising and falling in synchrony. (E) Recurrence plot generated from single element trace in (C). This graphical analytical tool, analogous to a spectrogram, presents cyclical processes as diagonal lines to allow the timing and distribution of cycles to be analyzed and contrasted with singular moments that do not repeat (represented as white space), points that recur only once (singular black points), or periods of stability where concentrations are relatively constant over time (vertical or horizontal lines); these structures are emphasized in the circular inset. During recurrence quantification analysis (RQA), the duration of cycles is captured by measuring MDL, robustness (determinism), and complexity (entropy). (F) Cross-recurrence plot for the dual element cycles presented in (D).

Fig. 2

Disruption of zinc-copper cycles in ASD cases and controls.

(A to C) Mean diagonal length (A), entropy (B), and determinism (C) are reduced in autism cases (squares) compared to controls (circles) in all study populations, indicating that zinc-copper cycles are of shorter duration, lower complexity, and reduced regularity in cases. Pooled estimates generated by combining data from all studies. Data are means ± 95% confidence intervals.

Zinc-copper cycles in ASD cases also exhibited significantly reduced entropy (a measure of biochemical network complexity), in contrast to non-ASD controls who had more complex cycles (F = 57.3, P < 0.001; Fig. 2B). After correction for multiple comparisons, we found no evidence that this effect varied across populations; that is, the attenuation of entropy in ASD cases was consistent across studies.

Finally, we found that determinism, a measure of cycle regularity and resistance to perturbation, was significantly reduced in ASD patients (F = 47.99, P < 0.001; Fig. 2C). As with entropy, we found that the ASD-related attenuation of determinism was consistent across populations.

Accuracy of elemental dynamics for classifying ASD

Given the significant differences we observed between ASD cases and neurotypically developing controls across multiple dynamical features of elemental cycles, we implemented an ensemble machine learning classification method generalized from weighted quantile sum (WQS) regression () to gauge how the full range of cyclical features classifies subjects as ASD cases or controls independently of their population of origin. Using optimal threshold criteria, this model was 90% accurate in predicting ASD cases, with 100% sensitivity for ASD diagnosis and 85% specificity to controls. Figure 3A shows the classification performance of this model, as applied to the holdout data set (15% of data). To confirm that robust classification could be achieved independently of the predictive algorithm applied, we also implemented an alternative approach based on penalized logistic regression [least absolute shrinkage and selection operator (LASSO)]. This approach achieved similar overall model accuracy (90%; model performance summarized in Fig. 3B), with slightly less sensitivity but improved specificity. Model performance relative to WQS is summarized in Fig. 3C; additional details on the construction of classification algorithms are presented in the Supplementary Materials.

Fig. 3

Performance of WQS-regression and penalized logistic regression algorithms in classifying ASD cases and controls.

(A) Receiver operating characteristic curve showing classification performance of the WQS algorithm with varying threshold values applied to the holdout data set (15% of data). AUC, area under the curve. (B) Classification performance of a penalized logistic regression algorithm applied to a holdout data set (15% of data) following 10-fold cross-validation in a training data set (85% of data). (C) Model performance characteristics.

DISCUSSION

We tested the hypothesis that disruption of zinc-copper dynamics underlies the elemental dysregulation that has been observed in ASD (–). We used tooth-matrix biomarkers to measure detailed temporal profiles of zinc, copper, and other elements from the second trimester to approximately 1 year postnatally. To characterize the coprogression of zinc and copper cycles throughout early development, we used cross-RQA (CRQA), a mathematical approach that originates in dynamical systems theory. We identified cycles varying on an approximately 10-day period and found strong evidence for abnormalities in zinc-copper cycles in ASD characterized by shorter duration, lower complexity, and less determinism. Further, we found that the integration of these dynamical features in an ensemble machine-learning algorithm allowed for a robust classification of ASD cases and controls, which we were able to replicate with similar accuracy with an alternative predictive algorithm based on penalized logistic regression. Together, these findings suggest that the normally well-coordinated and intricate metabolism of these essential elements is fragmented in ASD.

Metals have long been implicated in autism (–). Early life exposure to toxic metals and deficiencies of nutritional elements have been linked with several adverse developmental outcomes frequently associated with ASD, including intellectual disability and language, attentional, and behavioral problems (). Animal studies show that the effects of various metals on brain development could be mediated through dysregulation in neurotransmission and alterations in frontal and subcortical brain structures (), which have also been implicated in ASD (). Here, we found that measures of metal rhythmicity, but not metal concentration, were associated with ASD. These data emphasize that the dysfunctions apparent in cyclical properties are not evident in raw measures of concentration, which highlights the utility of our approach. This does not suggest that elemental concentrations are not relevant to the etiology of ASD, as we have previously published findings linking elevated levels of several elements (). However, that study dealt exclusively with a twin sample; without the inherent control of genetic and environmental confounds in nontwins, those differences in elemental concentrations may not be apparent.

Dynamic rhythms are central to human health and disease and have been reported in many systems (, ). Notably, much research has been carried out on imprints of rhythms in mammalian teeth (, ), making naturally shed deciduous human teeth an ideal biomatrix to study the role of metal rhythmicity and its disruption in childhood disorders, as we have carried out here for ASD. A key feature of our approach is the focus on dynamic rhythms in the internal metabolism of essential elements and toxic metals, rather than solely relying on an exposure paradigm that emphasizes external environmental exposures. By doing so, we could identify disrupted zinc-copper (and zinc-lead) cycles, even in populations where exposures were not excessive. In contrast to previous studies that have used blood or urinary biomarkers to investigate cyclical elemental dynamics in humans over periods less than a week, our methods allowed for a longitudinal analysis over a prolonged developmental period extending from the second trimester to approximately 1 year of age. An important consideration is that these differences in the developmental periods studied and temporal resolution between our study and those using repeated blood or urine measurements do not imply that the underlying processes being observed are dissimilar or unrelated; rather, we suggest that the cyclical processes we capture are superimposed on daily oscillations that have been previously described (, ). Independent of the duration of the cycles studied, we consistently found dysregulated cycle complexity and regularity in ASD cases across multiple replication samples, which bolsters our results against random measurement errors, and served as the basis of our classification analysis.

Critically, the focus on elemental dynamics represents a novel theoretical perspective in the consideration of factors predictive of ASD. While the role of elemental exposures in ASD is well studied, previous investigations have considered primarily the role of exposure intensities, that is, biomarker concentrations. Our results provide a new perspective, in that the temporal organization of elemental assimilation, as well as the interaction of multiple elemental pathways, appears critical to the emergence of ASD. This theoretical perspective, which we refer to as systemic elemental dysregulation (SED), opens new avenues for investigating elemental dysregulation in ASD and other disorders and yields testable hypotheses for future investigations to pursue. For example, our finding of dysregulated zinc-copper rhythmicity suggests underlying dysregulation in mechanisms involved in zinc-copper metabolism, which will be pursued in future studies with in vivo animal models and in vitro–induced pluripotent stem cells, where experimental manipulation of zinc-copper dynamics is possible. Consistent with this, a recent study () identified zinc deficiency and attenuated zinc transporter expression in humans with autism-related mutations in the SHANK3 gene, consistent with the SED etiological perspective of dysregulated metabolic processes (). In addition, our results emphasize the need to investigate other diseases where elemental assimilation is known to play a role, such as, attention deficit hyperactivity disorder (ADHD), which, to date, has been primarily studied through the lens of exposure intensity.

We also studied other elements because disruption of zinc and copper cycles may have downstream effects, including incomplete metabolism of other elements and toxic metals (, ). While we did not directly test biochemical pathways underlying interactions between zinc, copper, and other elements, we observed alteration in the joint rhythmicity of zinc and lead characterized by smaller and less complex cycles (fig. S1). Shorter and lower entropy interactions suggest that the regulatory effect of zinc on toxic metals, which is generally perceived as protective, is perturbed. We also found differences in other metals between cases and controls that were restricted to some, but not to all, study samples (see tables S1 and S2). This may be due to population-specific dynamics, such as differences in diet or exposure to pollution, which we have not analyzed in depth here.

Strengths of this study include the use of recently developed and validated methods to measure metals in teeth that provide detailed temporal profiles of elemental uptake during the fetal and early childhood periods. The use of population-derived discovery cohorts and the fact that our results could be replicated in ASD children from several geographical regions further bolster the validity of our findings. Limitations include small discovery and replication populations, but together, our sample size is approximately 200 participants, which was sufficient to detect and replicate significant elemental disruptions. Nonetheless, it is possible that other elements are also affected, which would only be detected with larger sample sizes. Our laser-based sampling method did not provide sufficient sampling points to separately study metal rhythmicity during the prenatal period from that observed after birth. Given that metal metabolism processes can incorporate toxic elements with similar chemical properties (cadmium and mercury, among others), the dysregulatory mechanism described herein may also play a role in their well-described neurotoxicity. Furthermore, other neuroactive metal regulatory pathways, such as iron and manganese, deserve to be studied as well.

In conclusion, in a discovery set of twins, including monozygotic twins discordant for ASD and three independent replication population samples, zinc-copper cycles were consistently altered in ASD. Furthermore, quantitative biochemical measures of the joint cyclic properties of zinc and copper accurately distinguished ASD cases from controls.

MATERIALS AND METHODS

Study participants

Our participants were recruited from four different studies being undertaken in three countries. Ethics clearances were obtained from the relevant institutional review boards at the coordinating site of each study. Key characteristics of the studies are outlined in Table 1, and details of participant recruitment and ASD case ascertainment are provided in Supplementary Materials and Methods. Briefly, a discovery analysis was undertaken in the Roots of Autism and ADHD Twin Study in Sweden (RATSS), which is designed to investigate the genetic and environmental determinants of ASD. Twins were recruited from nationwide registries and by advertisements in Sweden (). From the whole RATSS cohort, we collected and analyzed teeth from 75 participants (32 complete twin pairs and 11 individuals from twin pairs whose sibling did not donate a tooth). This subsample constituted 26% of the RATSS cohort and 50% of all participants in RATSS of tooth-shedding age.

Table 1

Characteristics of study participants.

N/A, not available.

Study	Location	Design	N (cases)	Male/female	Mean gestational days (SD)
RATSS	Stockholm, Sweden	National prospective twin cohort	75 (20)	46:29	247 (9)
ALSPAC	Bristol, UK	Case-control nested in prospective cohort	50 (25)	36:14	271 (21)
Seaver Autism Center	New York, USA	Hospital-based case-control	18 (10)	12:6	258 (21)
Autism Tooth Fairy Study	Texas, USA	Community-based case-control	50 (25)	25:25	N/A

For our primary replication analysis, we obtained teeth and supporting data from 25 ASD cases and 25 gender-matched controls enrolled in a long-term pregnancy cohort in the UK—the Avon Longitudinal Study of Parents and Children (ALSPAC) (). We undertook additional replication analysis by including ASD-diagnosed individuals (n = 10) at an autism clinic in New York, USA and their unaffected siblings (n = 8) (). Finally, we enrolled a random sample of 25 cases and 25 controls from a community-based study in Texas, USA—the Autism Tooth Fairy Project—that collects teeth and supporting data from parents through a community network ().

Laboratory assessments

Details of tooth collection protocols of each study population are provided in the Supplementary Materials. Our approach to measure metals in teeth using LA-ICPMS and assign developmental times was detailed elsewhere (, ). Briefly, teeth were sectioned, and the neonatal line (a histological feature formed in enamel and dentin at the time of birth) and incremental markings were used to assign temporal information to sampling points. A New Wave Research NWR-193 (ESI) laser ablation unit equipped with a 193-nm ArF excimer laser was connected to an Agilent Technologies 8800 triple-quadrupole ICPMS (Agilent Technologies). Helium was used as a carrier gas from the laser ablation cell and mixed with argon via Y-piece before introduction to the ICPMS. The system was tuned daily using NIST SRM 612 (trace elements in glass) to monitor sensitivity (maximum analyte ion counts), oxide formation (232Th16O+/232Th+, <0.3%), and fractionation (232Th+/238U+, 100 ± 5%). The laser was scanned in dentin parallel to the dentinoenamel junction from the dentin horn tip toward the tooth cervix. A preablation scan was run to remove any surface contamination. Data were analyzed as metal-to-calcium ratios (for example, 208Pb/43Ca) to control for any variations in the mineral content within a tooth and between samples. On average, each tooth was sampled at >150 locations. LA-ICPMS operating parameters are given in table S3.

Statistical analyses

Traditional statistical approaches to uncover the association of environmental stressors and disease outcomes, such as regression modeling, use concentrations at given time points as the fundamental units of analysis but do not resolve the dynamic cyclical nature of environmental exposures and their metabolism. To overcome this limitation, we used RQA and CRQA, which are nonlinear analytical methods for studying cyclical signal properties. These methods are well characterized with applications in diverse scientific fields, including physiology, molecular biology, geophysical sciences, and psychology [see reviews by Webber and colleagues (, ) and Marwan et al. ()]. The methods applied here, summarized in fig. S2, are described at length in Curtin et al. () and were implemented via the Cross-Recurrence Toolbox v5.16 () in MATLAB v2016b (MathWorks). Briefly, for RQA, time-series metal data from deciduous teeth were first used to construct RPs, with delay (τ) and embedding dimension (m) parameters calculated by the minimization of mutual information and false-nearest-neighbors algorithms, respectively (fig. S2). To facilitate cross-subject comparisons, threshold functions, ε, were constrained to fix signal recurrence rates to 10%. The temporal organization of features in the resulting RPs was then quantified via RQA to yield measures of MDL, recurrence time (RT), determinism, and entropy. These measures capture different features of diagonal line structures in RPs, which reflect cyclical signal components. MDL quantifies the average duration of these cycles, while RT measures the average interval between cycles. Determinism compares the ratio of cyclical to noncyclical recurrence points and is therefore a unit-less metric akin to a percentage, while entropy, a unit-less metric of predictability, measures the variability in cyclical lengths to capture the complexity of periodic signatures. For CRQA, these methods are extended to cross-recurrence plots, which capture the temporal organization of two-signal cyclical interactions over time (–).

Before inferential statistical tests, all variables were evaluated to confirm assumptions of normality in variable distributions. For analyses of recurrence features (RQA/CRQA), linear mixed models were used to test main effects of ASD diagnosis on features while also testing for potential interactions of ASD diagnosis and study population to evaluate population-specific effects. Familial relationships among subjects, including siblings and twins, were modeled as random variables to account for nonindependence of twins and siblings while also controlling for sex as a covariate. Additional covariates, including zygosity, gestational age, birth weight, and diagnosed comorbidities, including ADHD and intellectual disability, were also initially modeled, but these were ultimately excluded as they caused only negligible adjustments in our primary predictors and had no significant associations with outcomes. P values for study-specific differences reflect post hoc analyses of ASD × study interactions with Bonferroni adjustments for multiple comparisons. P values describing differences in the overall pool of subjects reflect the main effects of ASD diagnosis with false discovery rate P value adjustments.

Two classification models were undertaken: a generalization of WQS regression () to a binomial distribution to create a logistic classification algorithm or the application of penalized logistic regression (LASSO) (). Briefly, for WQS, data from all study populations were pooled and then divided into training (34% of sample size), validation (51%), and holdout (15%) subsets. The training data set was analyzed via WQS regression, yielding an empirically weighted index characterizing the mixture of recurrence variables, a β parameter estimating the association of this weighted index with ASD diagnosis, and an index of weights describing each individual variable’s contribution to the overall index. Weights derived from the training set were then used to construct and test a WQS index in the validation data set, and the resulting β estimates for the index were used to calculate a regression model for predicting outcomes in the holdout data set, which was entirely naïve to estimation of either weights or regression parameters. We randomly selected participants from all the populations for the classification analysis so that no single study population can drive the classification accuracy (or specificity or sensitivity) and significant features would be predictive independently of population.

For predictive models based on penalized logistic regression, data were initially randomly subset into training (85% of participants) and holdout (15% of participants) data sets. As in the implementation of WQS, all measures derived from RQA (MDL, entropy, determinism, and RT) for all elemental pathways and cross-recurrences studied were included in model building. Tenfold cross-validation was then applied in the training data set to identify a model with minimal deviance. The optimal LASSO model was then used to predict case/control status in the naïve holdout data set. Additional details on model formulation and implementation are provided in the Supplementary Materials.

In evaluating the performance of predictive classifiers, sensitivity was calculated as true-positive predictions/number of positive cases, specificity was calculated as true-negative predictions/number of negative cases, and accuracy was calculated as the number of correct predictions/number of total predictions. We confirmed the suitability of undertaking a pooled analysis by establishing the equivalence of the control subjects across the four study populations (see fig. S3 and the “Bioequivalence testing” section in the Supplementary Materials).