Enhanced pupillary light reflex in infancy is associated with autism diagnosis in toddlerhood

Autism spectrum disorder (ASD) is a neurodevelopmental condition affecting around 1% of the population. We previously discovered that infant siblings of children with ASD had stronger pupillary light reflexes compared to low-risk infants, a result which contrasts sharply with the weak pupillary light reflex typically seen in both children and adults with ASD. Here, we show that on average the relative constriction of the pupillary light reflex is larger in 9–10-month-old high risk infant siblings who receive an ASD diagnosis at 36 months, compared both to those who do not and to low-risk controls. We also found that the magnitude of the pupillary light reflex in infancy is associated with symptom severity at follow-up. This study indicates an important role of sensory atypicalities in the etiology of ASD, and suggests that pupillometry, if further developed and refined, could facilitate risk assessment in infants.


Supplementary participant information 1.a. EASE sample
ASD diagnosis of the older sibling (hereafter: proband) in the HR group was confirmed by a psychologist led interview with parents and by reviewing medical records (in more than 70% of cases it was explicitly stated that the Autism Diagnostic Observation Schedule (ADOS) and/or the Autism Diagnostic Interview-Revised (ADI-R) were used as part of the diagnostic evaluation). Infants in the low-risk control group were recruited from a volunteer database.
Inclusion criteria included full-term birth and lack of any ASD within second-degree family members (as confirmed through parent interview regarding family medical history). All included low risk infants had at least one older sibling. All families received a compensation of € ~50 for their participation. The final EASE sample consisted of 15 infants in the TD group, 13 in the HR-ASD group, and 25 in the HR-no-ASD group.

1.b. BASIS sample
High-risk infants had at least one older sibling with a community clinical diagnosis of ASD.
Proband diagnosis was confirmed by an expert clinician (TC) based on information using the Development and Well-Being Assessment (DAWBA) 1

and the parent-report Social
Communication Questionnaire (SCQ) 2 . Parent-reported family medical histories were examined for significant medical conditions in the proband or extended family members, with no exclusions made on this basis. Infants in the low-risk control group were recruited from a volunteer database. Inclusion criteria included full-term birth, normal birth weight, and lack of any ASD within first-degree family members (as confirmed through parent interview regarding family medical history). All low-risk participants had at least one older sibling.
Screening for possible ASD in these older siblings was undertaken using the SCQ, with no child scoring above instrument cut-off for ASD. The final BASIS sample consisted of 25 infants in the TD group, 16 in the HR-ASD group, and 93 in the HR-no-ASD group.

Extraction of the pupillary light reflex parameters
Pupil samples outside the range 1-10 mm diameter were considered artefactual and removed.
Short gaps in the data series (<7 samples) with <0.2 mm pupil change were linearly interpolated, to account for occasional "flicker" in eye tracking data collection and create continuous pupil traces for artefact rejection. All samples with a change in pupil size of >0.3 mm/sample were removed, and all data segments shorter than 5 samples were removed, as such data islands typically deviated from the surrounding pupil readings. An additional linear interpolation of missing data (same settings as above) was performed. To achieve better temporal resolution 3 the data was resampled to 300 Hz. First and second order derivatives of the pupil size, yielding pupil velocity and acceleration respectively, were calculated using the MATLAB diff function, and 25-point moving average filters were applied to the pupil size, pupil velocity and pupil acceleration to reduce amplification of noise during derivation 4 .
The PLR latency was defined by the acceleration minima in the time interval 100 -400 ms relative stimulus onset 3,4 . The baseline pupil size for each trial was defined as the average pupil size in a 100 ms interval just before the latency time point. The relative constriction of the PLR was calculated as described in the main text, within the interval 500-1500 ms after the flash.
Trials were excluded if 1) the latency interval (100 -400 ms relative stimulus onset) had more than 75% interpolated data; 2) if the relative constriction interval (500-1500 ms) had more than 75% interpolated data; 3) if the relative constriction was outside a range of 5 to 80 percent, out of which the PLRs were too small to be visually distinguished from noise, or to high to be biologically plausible); 4) if the latency or max relative constriction time points were in or adjacent to interpolated data. This automatic exclusion procedure was validated by visual inspection blinded to outcome group (pupil traces from left and right eye for each trial was inspected separately), and trials that did not resemble a PLR were manually rejected.
Only one eye contributed with PLR data; data from the left or right eye that correlated best with the individual's mean trace (i.e. all valid trials averaged within subject) was selected automatically.
Infants with less than three valid trials at the 9-10 month time point were excluded from further analysis (total n=14; TD n=4; HR-no-ASD n=5; HR-ASD n=5). Q-Q plots were used to inspect data distributions for outliers and deviations from normality before statistical analysis and to inspect residuals after model fitting.
Supplementary note 1

1.a. MSEL
To investigate whether the main finding could be confounded with a general intelligence factor we included the MSEL as a covariate in the main model (i.e. relative constriction at 9-10 months as dependent variable, group as fixed factor, and MSEL at 9-10 months as covariate. The results showed no effect of the MSEL: F(1, 183)=0.024, P=0.495, ηp²=0.003, and the other effects remained unchanged. Although the MSEL did not distinguish between groups, the scores might share some meaningful relationship with the PLR. To assess this at the most general level we performed Pearson correlations of all participants pooled together (n=187). The results did not show any significant association between MSEL and the relative constriction: r=0.012, P=0.872 (see Supplementary Fig. 1). Similarly, we found no relationship between relative constriction (9-10 months) and MSEL at 36 months ( Supplementary Fig. 2).

1.b. Missing data
Another possible confound is that the data quality from participants varied in systematic ways between groups. To investigate this issue, we calculated the amount of missing data for each trial and used the average within participants as a dependent variable in a GLM with group as fixed factor. Because the dependent variable did not show equal error variances between groups (Levene's test) we ran the analysis again on the log transformed dependent variable.
This model did not show any significant differences between groups: F(2, 191) = 0.360, P=0.698, ηp²=0.004, and neither did any pairwise comparison between groups. Further, adding the amount of missing data as a covariate in the main relative constriction analyses did not show any interaction effect with group: F(2, 181)=0.088, P=0.916, ηp²=0.001, and the main effect of group and pairwise comparisons between groups remained unchanged. This suggests that differences in data quality between groups did not confound our main findings.

1.d. Sensitivity analysis (reference group)
As noted in the main text, we used the TD group as our first choice for site normalization, as we assumed that this group was most likely to be similar across our two sites given the similar context (recruited from volunteer database for "babylabs", both being situated in major cities in Europe). However, as the group comparison for the combined dataset is influenced by this choice, we repeated the main analysis using the other groups as normalization as well. With

a. PLR latency
The same GLM model structure as in the main text was used to assess differences between groups, with group and fixed factor and TD normalized latency as dependent variable. As

2.b. Baseline pupil size
Because there could be group differences in arousal that affect the pupil responses we tested whether the baseline pupil size differed between the groups. In this analysis the TD normalized baseline was used as the dependent variable and group was used as fixed factor. Also, adding the TD normalized baseline values as a covariate in the main analysis of relative constriction did not change the main effect of group or any of the pairwise comparisons, which suggests that our main finding was not biased by baseline values.

Supplementary note 3 3. Site effects 9-10 month
The two sites presented different types of stimuli to elicit PLR responses, and we accounted for site related variability by normalizing all values by each site's average TD value. To investigate group differences within each site, we performed the same GLM as in the main analysis for each site separately. In sum, these analyses show that descriptively, the same pattern of group differences (TD < HR-no-ASD < HR-ASD) replicates in both sites, but that the pattern of significant results in specific tests differ between sites. In the following analyses the GLM models used TD normalized relative constriction as dependent variable and group as fixed factor.

3.a. EASE relative constriction
The results showed a significant main effect of group, F(2, 50) = 11.958, P<0.001, ηp²=0.324, and planned comparisons showed a difference between the HR-ASD and the TD group Because the PLR was measured following a black gap between two white stimuli in the BASIS sample, we investigated whether the gap length, which varied slightly between subjects as stated in the main text, had any impact on the relative constriction results in the BASIS sample. We added the gap length as a covariate to the GLM model and investigated both main and interaction effects with group. The results did not show any significant main effect of gap length: F(1, 128)=1.731, P=0.191, ηp²=0.013, nor a significant interaction: F(2, 128)=1.531, P=0.220, ηp²=0.023), suggesting that this possible error source did not influence the other results.
In the BASIS sample, after the first time point 54 of the high risk participants was enrolled in a randomized controlled trial (RCT) of parent-mediated intervention 5 and another four participants in a similar non-RCT intervention. To rule out any confounding contributions from the intervention programmes we performed the following analyses using two binary factors as such: 1) To account for differences due to different sampling of the groups recruited or not recruited for intervention, we tested for a main effect of a recruitment factor (1=participating in the RCT, 0 =not part of the RCT); 2) To account for effects of intervention treatment on the relative constriction, we tested for a main effect of the treatment factor (1=treated, 0=not treated); 3) To account for any moderating effect treatment might have had on the relationship between relative constriction and diagnostic outcome, we tested for an interaction effect between treatment and diagnostic outcome. Because there were no significant effects in any these analyses, neither the recruitment factor nor the treatment factor were used in any further analysis.

Supplementary note 4 4. Longitudinal analysis of relative constriction
As stated in the main text, our finding of increased PLR in HR-ASD at 9-10 months of age are opposite to previous findings of decreased PLR in children with ASD 6 . We therefore asked whether developmental changes occur sometime between infancy and middle childhood which would lead to a crossover effect between HR-ASD and TD groups. Our data allowed for a longitudinal analysis of the relative constriction at 9-10 months and 14-15 months.
Despite the short time span we found a significant difference in the change scores between the