Polygenic risk scores for schizophrenia are associated with oculomotor endophenotypes

Background Schizophrenia is a heterogeneous disorder with substantial heritability. The use of endophenotypes may help clarify its aetiology. Measures from the smooth pursuit and antisaccade eye movement tasks have been identified as endophenotypes for schizophrenia in twin and family studies. However, the genetic basis of the overlap between schizophrenia and these oculomotor markers is largely unknown. Here, we tested whether schizophrenia polygenic risk scores (PRS) were associated with oculomotor performance in the general population. Methods Analyses were based on the data of 2956 participants (aged 30–95) of the Rhineland Study, a community-based cohort study in Bonn, Germany. Genotyping was performed on Omni-2.5 exome arrays. Using summary statistics from a recent meta-analysis based on the two largest schizophrenia genome-wide association studies to date, we quantified genetic risk for schizophrenia by creating PRS at different p value thresholds for genetic markers. We examined associations between PRS and oculomotor performance using multivariable regression models. Results Higher PRS were associated with higher antisaccade error rate and latency, and lower antisaccade amplitude gain. PRS showed inconsistent patterns of association with smooth pursuit velocity gain and were not associated with saccade rate during smooth pursuit or performance on a prosaccade control task. Conclusions There is an overlap between genetic determinants of schizophrenia and oculomotor endophenotypes. Our findings suggest that the mechanisms that underlie schizophrenia also affect oculomotor function in the general population.


Background
Schizophrenia is a severe mental disorder with a lifetime prevalence of just under 1% (McGrath, Saha, Chant, & Welham, 2008). There is substantial evidence for a genetic basis of schizophrenia, with recurrence risk in families of about 8.6% (Lichtenstein et al., 2006) and heritability estimates of up to 81% (Sullivan, Kendler, & Neale, 2003). Genome-wide association studies (GWASs) thus far have identified 145 single nucleotide polymorphisms (SNPs) that are associated with schizophrenia (Pardiñas et al., 2018). However, the genetic variance of schizophrenia explained by these SNPs is low (Pardiñas et al., 2018). On the one hand, there are many SNPs that do not reach the genome-wide significant threshold in a GWAS (5 × 10 −8 ) but that could explain in sum a substantial proportion of genetic variance (International Schizophrenia Consortium, 2009). Another reason for this may be that schizophrenia is highly heterogeneous and encompasses a multitude of different syndromes that do not necessarily have a common biological basis (Braff, Freedman, Schork, & Gottesman, 2007).
Endophenotypes have been proposed as an approach to better understand the aetiology of schizophrenia . They are considered to link a disorder to its genetic basis and to be closer to the actions of genes than disease symptoms are (Gottesman & Gould, 2003). Furthermore, as current disease classification might not well reflect aetiology (The Brainstorm Consortium, 2018), endophenotypes may help identify more homogenous subgroups of patients with shared biological basis .
In the SPEM task, participants follow a slowly moving target with their eyes. Patients with schizophrenia have long been known to have SPEM impairments (Diefendorf & Dodge, 1908;Holzman, Proctor, & Hughes, 1973), characterised primarily by lower velocity gain (ratio of eye and target velocity) and higher saccade rates than controls (O'Driscoll & Callahan, 2008;Sereno & Holzman, 1995). Greater deficits in antisaccade error rate have also been observed in patients with schizophrenia who have impaired SPEM performance compared to patients with schizophrenia without SPEM impairments (Sereno & Holzman, 1995).
Here, we used SNPs previously associated with schizophrenia in GWASs to investigate whether genetic determinants of schizophrenia are associated with oculomotor endophenotypes in a large population-based cohort study. To this end, we calculated polygenic risk scores (PRS) based on the summary statistics of the largest schizophrenia GWAS to date (Pardiñas et al., 2018). The SNPs that have been identified so far in the schizophrenia GWAS have been associated with, inter alia, voltage-gated calcium channels, synaptic transmission, membrane depolarisation during action potentials, and fragile X mental retardation protein (FMRP) (Pardiñas et al., 2018). We hypothesised that higher PRS would be associated with worse antisaccade and SPEM performance (higher error rate, latency and saccade frequency during SPEM but lower amplitude gain and SPEM velocity gain) but unrelated to performance in the prosaccade control task (latency and amplitude gain).

Participants
We used data from participants of the Rhineland Study, a community-based cohort study in Bonn, Germany. All inhabitants of two geographically defined areas in Bonn who are 30 years or older are invited to participate in the Rhineland Study. Names and addresses were provided by the municipality. Study participation is possible upon invitation only and irrespective of health status. The only exclusion criterion is not having sufficient command of the German language to provide written informed consent. There are no financial incentives for study participation. The ethics committee of the Medical Faculty of the University of Bonn approved the study, which was carried out in accordance with the recommendations of the International Council for Harmonisation Good Clinical Practice standards (ICH-GCP). We restricted our sample to the first 4000 participants of the Rhineland Study. Since study recruitment is ongoing, we cannot provide information on general response rates, but 3523 participants (88.1%) of those first 4000 participants provided blood samples between March 2016 and July 2019. Of those, 3217 (91.3%) remained after quality control of genetic data (see section 'Genetic data and polygenic risk scores'). Of those, 250 participants (7.8%) had no SPEM and antisaccade data. Missing data were primarily due to technical issues during data acquisition and post-processing (74%), exclusion after visual inspection of data (9.6%), contraindications (8.8%), non-compliance (5.2%), refusal (0.8%), timeout (0.4%) or multiple of these reasons (1.2%). Finally, we excluded two individuals with a diagnosis of schizophrenia and nine individuals with a diagnosis of psychosis. Thus, we based our analysis on 2956 participants without schizophrenia or psychosis aged between 30 and 95 years.
frequencies (Price et al., 2006) [exclusion of N = 101 (2.9%) non-Caucasian participants]. Finally, missing SNPs were imputed using IMPUTE version 2 software (Howie, Donnelly, & Marchini, 2009) based on the 1000 Genomes reference panel (Auton et al., 2015). Imputation quality of the SNPs was checked using the info score metric [values of >0.3 are considered to indicate reliable imputation quality (Verma et al., 2014)]. PRS for schizophrenia were created using summary statistics from the largest schizophrenia GWAS to date, which included a discovery sample of 40 675 schizophrenia cases and 64 643 controls and an independent replication sample of 5762 cases and 154 224 controls (Pardiñas et al., 2018;Ripke et al., 2014). The results are publicly available (https://walters.psycm.cf.ac.uk; last retrieved at: 2021/05/31). We calculated PRS using PLINK (version 1.9) (Purcell et al., 2007) by first multiplying the number of risk alleles by the known effect size of each individual SNP locus and then aggregating the weighted effects of all SNPs under consideration (International Schizophrenia Consortium, 2009). We first created PRS based on the pre-specified SNPs from the GWAS, i.e. the 145 SNPs that reached genome-wide significance in the GWAS (Pardiñas et al., 2018). Then, we applied clumping to identify the most significant SNPs per linkage disequilibrium (LD) block (kilo base pair window: 250, LD r 2 < 0.1) (Chasioti, Yan, Nho, & Saykin, 2019) and created PRS at p value threshold ( p T ) for SNP inclusion of 0.01 and 0.05, since PRS at those thresholds were reported to have improved prediction accuracy (Jonas et al., 2019;Ripke et al., 2014;Toulopoulou et al., 2019;Zhang et al., 2019). For sensitivity analysis, we created two additional PRS. We created one PRS by first applying clumping and then using the genome-wide significant threshold for SNP inclusion ( p T = 5 × 10 −8 ), and another PRS using a more lenient threshold ( p T = 0.1), as in previous studies (Toulopoulou et al., 2019;Zhang et al., 2019).

Eye movement data
A detailed description of oculomotor data acquisition and processing has been published (Coors et al., 2021). In brief, eye movements were recorded using video-based infrared oculography (EyeLink 1000 and EyeLink 1000 Plus; SR Research Ltd, Ottawa, Canada) at 1000 Hz. After a horizontal-vertical five-point calibration task, participants performed fixation (not reported here), SPEM, prosaccade and antisaccade tasks in fixed order. SPEM outcomes were velocity gain (in %) and saccade rate (given in N/s, across the entire task duration). Prosaccade outcomes were latency (in ms) and amplitude gain (saccade amplitude divided by target step amplitude). Antisaccade outcomes were error rate (in %), latency and amplitude gain. Prosaccade and antisaccade outcomes were only calculated if there were at least seven valid trials. In case of more than four antisaccade errors, there also had to be at least one corrective saccade to ensure that participants understood the instructions. Additionally, latency and amplitude gain were only calculated if there were at least seven valid trials with directionally correct initial saccades. Before applying those criteria, we performed sensitivity analysis to rule out the possibility that they led to the exclusion of the participants with the highest PRS as this could have explained invalid or poor performance. Since we found no systematic pattern, we excluded those cases (51 participants for antisaccade error rate, 305 participants for antisaccade latency and amplitude gain and 14 participants for prosaccade latency and amplitude gain).

Statistical analyses
We hypothesised that high genetic risk for schizophrenia would be associated with worse antisaccade and SPEM performance but not with prosaccade outcomes (Calkins et al., 2008).
First, linear regression model assumptions were tested with diagnostic plots (scale-location plot and quantile-quantile plot) and by calculating the variance inflation factor [R package car (Fox & Weisberg, 2019), vif-function]. For prosaccade and antisaccade latency, the normality assumption was violated and therefore we log-transformed those variables.
Then, we assessed the associations between PRS and oculomotor outcomes with separate multivariable linear regression models for SPEM, antisaccade and prosaccade outcomes. Regression models included z-standardised PRS as a predictor and were adjusted for age, age 2 , sex and population stratification. For the latter, we calculated six principal components that we included as covariates in the model (Price et al., 2006). We used mean-centred age to reduce collinearity (Iacobucci, Schneider, Popovich, & Bakamitsos, 2016). Missing covariate data were imputed using predictive mean matching [Hmisc package, 10 bootstrap replicates (Harrell & Dupont, 2020)].
We did not correct for multiple testing as we had very specific a priori hypotheses regarding associations between the schizophrenia PRS and eye movement outcomes based on work that goes back decades (Diefendorf & Dodge, 1908;Fukushima et al., 1988;Holzman et al., 1973;Sereno & Holzman, 1995). Further, we included the prosaccade task as control condition and created additional PRS for sensitivity analysis. As argued elsewhere, correction for multiple testing is strongly context-dependent and can lead to misinterpretation of results if incorrectly applied (Rothman, 1990;Streiner & Norman, 2011). Multiple testing is considered inappropriate for a limited set of pre-specified hypotheses and becomes especially problematic if the statistical tests are not independent, which is clearly the case in our analyses where our predictors (i.e. PRS scores at different p value thresholds) are highly correlated (Streiner & Norman, 2011).
Given the large age range of our sample, we additionally tested whether the associations between PRS and eye movement outcomes varied with age using a likelihood ratio test. We also repeated the analyses in an age-truncated sample (participants aged 30-70 years).
Statistical analyses were performed in RStudio (version 1.1.447, R-base version 3.5.0), using an α level of 0.05.

Study sample
Sample characteristics are presented in Table 1.

Associations between PRS for schizophrenia and oculomotor performance
The main results from the multivariable regression models are listed in Table 2 and the results of the PRS that we created for sensitivity analysis are in online Supplementary Table A.1. PRS was positively associated with antisaccade error rate irrespective of inclusion criteria (pre-specified SNPs, p T = 0.01, p T = 0.05), but sensitivity analysis at p T = 5 × 10 −8 and p T = 0.1 was not significant. For antisaccade latency and amplitude gain, PRS at p T = 0.01 and p T = 0.05 were associated with those outcomes but PRS including the pre-specified SNPs and PRS created for sensitivity analysis were not. The associations were positive for antisaccade latency and negative for antisaccade amplitude gain. PRS at p T = 0.01 was positively associated with SPEM velocity gain but the other two PRS (pre-specified SNPs, p T = 0.05) were not significantly associated with it. In addition, sensitivity analysis revealed a positive association between PRS at p T = 0.1 and SPEM velocity gain. None of the PRS was associated with saccade rate. Regarding the prosaccade control task, none of the PRS was significantly associated with latency or amplitude gain.
Effects for all associations were small, with at most 0.22% of variance in oculomotor outcomes explained by the PRS.
We found no interaction effects between age and PRS. Effect estimates in the age-truncated analysis were highly comparable to those in the whole sample, but given the smaller sample size (N = 2636), some confidence intervals were wider. For the association between PRS at p T = 0.05 and antisaccade amplitude gain, this resulted in the inclusion of zero in the confidence interval, but the regression coefficient remained comparable (full sample: b = −1.279; 95% CI −2.491 to −0.067; age-truncated sample: b = −1.228; 95% CI −2.520 to 0.064).

Discussion
We investigated genetic determinants of schizophrenia in relation to oculomotor endophenotypes in a large, population-based cohort. We found that genetic variants that are associated with schizophrenia are also involved in the fine-regulation of particular aspects of oculomotor function. Schizophrenia-related genetic risk variants specifically affected antisaccade outcomes, but not saccade rate during SPEM or outcomes from the prosaccade control task. PRS showed inconsistent patterns of association with SPEM velocity gain. Whilst collectively these findings thus support the use of specific oculomotor endophenotypes as markers of those syndromes that are currently classified as schizophrenia, it should be noted that the effect sizes of the observed associations are small.
Our findings suggest that SNP inclusion thresholds of p T = 0.01 and p T = 0.05 were optimal for the detection of PRS correlates of eye movements, in line with the findings of previous studies (Jonas et al., 2019;Ripke et al., 2014;Toulopoulou et al., 2019;Zhang et al., 2019). Of the PRS that we calculated for sensitivity analysis, the p T = 0.1 cut-off might have been less optimal for this study as a more lenient threshold implies the inclusion of more uninformative SNPs and, therefore, an increase in noise (Chasioti et al., 2019). On the contrary, PRS at the genome-wide significant p T ( p T = 5 × 10 −8 ) might have excluded too many informative SNPs (Pardiñas et al., 2018). Our sensitivity analysis showed that applying clumping and then applying the genomewide significant threshold obscured the association between PRS and antisaccade error rate that we found when we used PRS based on only the pre-specified SNPs.
Supporting the endophenotype status of antisaccade latency and error rate, we found both to be positively associated with genetic risk for schizophrenia which is in line with the reports of deficits in these measures in patients with schizophrenia and their clinically unaffected relatives (Calkins et al., 2008;Curtis et al., 2001;Ettinger et al., 2004;Fukushima et al., 1988;Fukushima et al., 1990a;Fukushima et al., 1990b;Karoumi et al., 2001;Mazhari et al., 2011;Radant et al., 2007Radant et al., , 2010Reilly et al., 2014;Reuter et al., 2005;Sereno & Holzman, 1995). Antisaccade latency depends on cognitive processes such as attention, response-related decision-making and response execution (Hutton, 2008) and has been linked to activity in saccade neurons in the frontal eye fields and superior colliculus (Munoz & Everling, 2004). Reduced activation of frontal eye fields in individuals with schizophrenia compared to controls during eye movement tasks has been reported (Keedy, Ebens, Keshavan, & Sweeney, 2006) and may, therefore, partly account for the association between PRS and latency. Successful inhibition of antisaccade errors has been associated with activity in the dorsolateral prefrontal cortex (Munoz & Everling, 2004). Since individuals with schizophrenia have been found to have lower activity in prefrontal areas during antisaccade task performance than controls (McDowell et al., 2002), the association between genetic risk for schizophrenia and antisaccade error rate may be partly mediated by prefrontal areas, although the striatum may also play a role (Raemaekers, Ramsey, Vink, van den Heuvel, & Kahn, 2006. The negative association between genetic risk for schizophrenia and antisaccade amplitude gain is also in line with studies reporting lower antisaccade amplitude gain in schizophrenia patients (Ettinger et al., 2004;Karoumi et al., 2001;Radant et al., 2010) and their biological relatives (Ettinger et al., 2018,  N = number of participants, M = mean, S.D. = standard deviation. Education level was determined using the International Standard Classification of Education 2011 (ISCED) and was coded as low (lower secondary education or below), middle (upper secondary education to undergraduate university level) and high (postgraduate university study). Assessment of best-corrected visual acuity was based on visual scores from the right eye and was measured using an automated refractometer (Ark-1s, NIDEK CO., Tokyo, Japan). Categorisation of the visual acuity values was based on the guidelines of the International Council of Ophthalmology.

1614
Annabell Coors et al. 2006Coors et al. , 2004Karoumi et al., 2001), again supporting the endophenotype candidacy of this measure. However, in those studies, the mean antisaccade accuracy of controls typically ranged between 95% and 100% and, therefore, lower amplitude gain in patients or relatives indicated lower spatial accuracy. Instead, participants in our study on average tended to make hypermetric (overshooting) antisaccades (M = 111.7%, S.D. = 29.1; Table 1), implying that participants with higher PRS had in fact greater absolute spatial accuracy, as their scores were closer to 100%. Antisaccade amplitude gain values above 100% are not unusual and values comparable to ours have been reported by others (Sweeney, Rosanao, Berman, & Luna, 2001). Antisaccade spatial accuracy requires complex, non-standard sensorimotor transformations in the posterior parietal cortex (Herweg et al., 2014) and frontal eye fields (Moon et al., 2007). Together, these findings confirm a genetic overlap between schizophrenia and antisaccade amplitude gain and point to a general tendency of people with a higher genetic risk of schizophrenia to make antisaccades with lower amplitudes. However, since schizophrenia is a highly heterogeneous disease  and antisaccade performance has rather a low specificity for schizophrenia (Reilly et al., 2014), we cannot exclude the possibility that the associations between PRS for schizophrenia and antisaccade outcomes may also be partly accounted for by genetic risk variants of other psychiatric disorders.
The finding of higher SPEM velocity gain and, therefore, better performance in participants with higher PRS at p T = 0.01 was unexpected, given the highly consistent reports of lower velocity gain in patients with schizophrenia and their relatives compared to healthy controls (Calkins et al., 2008;O'Driscoll & Callahan, 2008). An explanation might be that higher genetic risk for schizophrenia but not having schizophrenia may be advantageous for performance in SPEM velocity gain. Genes associated with schizophrenia were found to be favoured by evolution which implies that they might be advantageous for (cognitive) functioning, at least to a certain degree (Banerjee et al., 2018;Srinivasan et al., 2016). However, it is unclear why those advantages should only exist in performance in one but not other oculomotor endophenotypes. In contrast to antisaccade outcomes, SPEM velocity gain was associated with PRS at p T = 0.01 and p T = 0.1 but unrelated to PRS at p T = 0.05. Since this is the only eye movement outcome for which the pattern was inconsistent for PRS at p T = 0.01 and p T = 0.05, and since the inclusion of more than the genome-wide significant SNPs in PRS creation may also increase the level of noise (Chasioti et al., 2019), an alternative explanation is that the relation of higher PRS with higher SPEM velocity gain may have been a false-positive finding.
Despite previous findings of higher saccade rate during SPEM in schizophrenia patients (O'Driscoll & Callahan, 2008) and their relatives (Calkins et al., 2008), we did not find an association of  The table displays the change in performance per one standard deviation increase in PRS for schizophrenia for different eye movement outcomes. b, unstandardised regression coefficient; 95% CI, 95% confidence interval. Unstandardised regression coefficients were obtained from the following multivariable linear regression model: Eye movement outcome ∼b 0 + PRS × b 1 + age + age 2 + sex + population stratification + residual error. R 2 refers to the variance explained in eye movement performance by PRS in per cent. In bold are those associations with a p value below 0.05.
PRS with that measure. This could be due to either a limited overlap in the genetic determinants between schizophrenia and saccade rate, or the current PRS capturing only a small proportion of those shared genetic factors. Our sensitivity analysis of the prosaccade control task showed that PRS were not associated with prosaccade latency and amplitude gain. Since performance in these outcomes is largely unaffected in patients with schizophrenia and their relatives (Calkins et al., 2008;Reilly et al., 2014), this finding corroborates our explicit a priori hypothesis that associations between PRS for schizophrenia and oculomotor outcomes are limited to those oculomotor outcomes that have been established as endophenotypes of schizophrenia.
From a genetic perspective, it is noteworthy that those SNPs that we included in the PRS were found to be associated with, inter alia, voltage-gated calcium channels and the FMRP (Pardiñas et al., 2018). Voltage-gated calcium channels play a key role in visual perception (Pangrsic, Singer, & Koschak, 2018) and mutations have been linked to visual deficits such as involuntary eye movements (Cain & Snutch, 2011). In mice, FMRP has been associated with prefrontal cortex dysfunction (Siegel et al., 2017), which is also a critical brain structure for successful antisaccade performance (Kaufman, Pratt, Levine, & Black, 2010). This suggests that shared mechanisms may underlie schizophrenia and oculomotor performance and that those shared mechanisms may become particularly evident in some specific oculomotor outcomes. However, we also know from genetics that there are systematic differences in allele frequencies between populations (Price et al., 2006). Thus, our findings may not be generalisable to non-Caucasian populations as we based our PRS on SNPs derived from a study including predominantly Caucasians (Wand et al., 2021).
Our sample included individuals aged between 30 and 95 years. Previous studies found that the majority of patients with schizophrenia develop the disease during adolescence and early adulthood (men: between age 10 and 25, women: between age 25 and 35) (Rajji, Ismail, & Mulsant, 2009). Approximately onequarter of patients with schizophrenia, and particularly women, experience their first episode after the age of 40 and very few patients are diagnosed after age 60 (Rajji et al., 2009). Thus, the probability that our population-based sample included individuals that are about to develop schizophrenia but have not yet been diagnosed is very low. Our large sample size and the young typical age of onset benefit our research aim as these factors lower the risk that the observed associations between schizophrenia PRS and eye movement performance were due to individuals about to develop schizophrenia. We found no evidence that the observed associations varied with age, yet we lacked the statistical power to run age group-specific analyses in more narrow age ranges. Further research, conducted in large samples with narrower age ranges, is needed to confirm the associations we found between genetic liability for schizophrenia and oculomotor measures.
The observation that only small amounts of variance in established oculomotor endophenotypes could be explained by PRS needs critical examination. One possibility is that current PRS do not fully capture the genetic basis of schizophrenia. The estimate for the common-variant SNP heritability of schizophrenia calculated in the largest GWAS to date is 24.4% if all SNPs are considered, and only 6% for PRS at p T = 0.05 (Pardiñas et al., 2018). These estimates are, therefore, well below family-based studies heritability estimates (Lichtenstein et al., 2009;Sullivan et al., 2003). Current PRS may capture only a fraction of genetic variance attributed to schizophrenia because even the largest GWASs to date were not sufficiently powered to detect all relevant common variants (Smeland, Frei, Dale, & Andreassen, 2020). In addition, part of the heritability results from copy number variants or rare variants that influence the boundaries of topologically associated domains (Halvorsen et al., 2020;Marshall et al., 2017), which are currently not tagged by conventional genotyping arrays and, therefore, not included in GWAS (Auer & Lettre, 2015).
Further, it should be remembered that heritability estimates for both schizophrenia (Lichtenstein et al., 2009;Sullivan et al., 2003) and eye movements (Bell et al., 1994;Greenwood et al., 2007;Hong et al., 2006;Katsanis et al., 2000;Litman et al., 1997;Macare et al., 2014;Malone & Iacono, 2002) are well below 100%. This, and the only modest sized oculomotor impairments in first-degree relatives of schizophrenia patients (Calkins et al., 2008), implies that oculomotor impairments in schizophrenia reflect not only genetic but also environmental factors as well as the interplay between genes and environment (Chakravarti & Little, 2003).
It should also be noted that our inclusion criterion of a minimum age of 30 years, combined with the typically rather early onset for schizophrenia (Rajji et al., 2009), may have led to the exclusion of some participants with very high genetic risk for schizophrenia, thereby reducing the variance in schizophrenia risk and oculomotor performance in our sample.
Taken together, the schizophrenia PRS alone is unlikely to fully account for differences in oculomotor performance. This also fits with our finding that current schizophrenia PRS do not have any predictive power for eye movement performance. Still, the molecular genetic confirmation implies that the role of those brain regions that are critically involved in antisaccade performance should be investigated more closely in the aetiology of schizophrenia. Thus, combining knowledge from eye movement and schizophrenia research could be beneficial to propel the field forward.

Conclusions
Using a molecular genetic approach, we confirm and extend previous findings from behavioural genetic studies, showing that antisaccade error rate, latency and amplitude gain have genetic overlap with schizophrenia. For SPEM outcomes, we found no association between PRS and saccade rate and inconsistent associations between PRS and velocity gain. As schizophrenia PRS based on currently available GWAS findings only accounted for <0.25% of variance in oculomotor endophenotypes, they currently have no predictive power. However, we expect that future studies using PRS that also include rare risk variants are likely to uncover a larger proportion of shared genetic determinants of schizophrenia and oculomotor performance.
Supplementary material. The supplementary material for this article can be found at https://doi.org/10.1017/S0033291721003251