Altered default mode network associated with pesticide exposure in Latinx children from rural farmworker families.

Pesticide exposure has been associated with adverse cognitive and neurological effects. However, neuroimaging studies aimed at examining the impacts of pesticide exposure on brain networks underlying abnormal neurodevelopment in children remain limited. It has been demonstrated that pesticide exposure in children is associated with disrupted brain anatomy in regions that make up the default mode network (DMN), a subnetwork engaged across a diverse set of cognitive processes, particularly higher-order cognitive tasks. This study tested the hypothesis that functional brain network connectivity/topology in Latinx children from rural farmworker families (FW children) would differ from urban Latinx children from non-farmworker families (NFW children). We also tested the hypothesis that probable historic childhood exposure to pesticides among FW children would be associated with network connectivity/topology in a manner that parallels differences between FW and NFW children. We used brain networks from functional magnetic resonance imaging (fMRI) data from 78 children and a mixed-effects regression framework to test our hypotheses. We found that network topology was differently associated with the connection probability between FW and NFW children in the DMN. Our results also indicated that, among 48 FW children, historic reports of exposure to pesticides from prenatal to 96 months old were significantly associated with DMN topology, as hypothesized. Although the cause of the differences in brain networks between FW and NFW children cannot be determined using a cross-sectional study design, the observed associations between network connectivity/topology and historic exposure reports in FW children provide compelling evidence for a contribution of pesticide exposure on altering the DMN network organization in this vulnerable population. Although longitudinal follow-up of the children is necessary to further elucidate the cause and reveal the ultimate neurological implications, these findings raise serious concerns about the potential adverse health consequences from developmental neurotoxicity associated with pesticide exposure in this vulnerable population.

Introduction

Farmworkers play a vital role in United States agriculture, performing much of the manual labor in this industry. These migrant and seasonal farmworkers experience chronic occupational exposure to pesticides, which places them at risk for acute and chronic health deficits (Hernandez et al., 2016; Kim et al., 2017; Kirkhorn and Schenker 2002). However, pesticide exposure is not limited to those who work in the field. Families of farmworkers experience significant exposure to pesticides as well due to living in close proximity to agricultural fields, through take-home pathways, and through prenatal exposure (Arcury and Quandt 2011; Fenske et al., 2002; Hyland and Laribi 2017; Schwartz et al., 2015). Children are particularly vulnerable to the adverse effects of pesticide exposure because they absorb higher doses of pesticides due to their larger surface to volume ratios, and also because the absorbed toxicants remain in their body longer, given their slower toxicant metabolism when compared to adults (Arcury et al., 2007). Higher doses of pesticides received by children, as well as reduced capacity for detoxifying xenobiotics, makes the neurotoxic effects potentially worse in this age group (London et al., 2012). Prenatal and postnatal exposure to pesticides in children have been associated with adverse effects on neurodevelopment in a growing body of studies (Gonzalez-Alzaga et al., 2014), including: adverse effects on cognitive development (Bouchard et al., 2011; Engel et al., 2011; Eskenazi et al., 2007; Kofman et al., 2006; Lizardi et al., 2008; Rauh et al., 2011; Ruckart et al., 2004), deficits in working memory (Rauh et al., 2011; Wang et al., 2016), and higher incidence of ADHD (Bouchard et al., 2010; Wagner-Schuman et al., 2015).

There is currently a dearth of studies examining neurobiological effects of pesticides on children’s brains. Although there are currently no published studies examining brain networks and pesticide exposure, there are two studies examining abnormal brain anatomy associated with prenatal pesticide exposure (Rauh et al., 2012; van den Dries et al. 2020). Using magnetic resonance imaging (MRI) data of 40 children, including 20 high- and 20 low-exposure children, Rauh et al. found altered brain gray matter anatomy along a large portion of the lateral temporal lobes as well as some smaller changes in lateral frontal and parietal regions. However, the findings presented in the first figure of that paper clearly demonstrate large scale changes on medial aspects of the brain. Most notably, the medial superior and orbital prefrontal cortex and the cuneus/precuneus regions exhibited pesticide exposure-related differences. The medial regions together with lateral parietal cortex make up the default mode network (DMN), an extensively studied brain subnetwork which is shown to play a critical role across a diverse set of cognitive processes (Broyd et al., 2009; Shen 2015), particularly higher-order cognitive processes (Smallwood et al., 2021). Although the DMN was not specifically discussed as a subnetwork in that manuscript, it is notable that the alterations found in the precuneus were significantly associated with cognition, as were some of the lateral brain changes. A recent study using Diffusion Tensor Imaging (DTI) data from 518 participants assessed the effects of prenatal exposure on brain morphology and white matter microstructure in preadolescents aged 9–12 years (van den Dries et al. 2020). This study demonstrated a strong association between prenatal exposure to organophosphate pesticides and altered white matter microstructure, indicating likely impacts on neurodevelopment in this population. Associations between maternal urinary organophosphate metabolites, particularly dimethyl alkyl phosphate, and white matter microstructure were concentrated in fiber tracts that interconnect regions of the DMN. It should be noted, however, that they did not identify gray matter associations as found by Rauh et al. In addition to the studies of brain anatomy, there is a single study using MRI to evaluate associations between brain activity and prenatal exposure to organophosphate pesticides (Binter et al., 2020). This study used fMRI data from 95 children aged 10–12 years, focusing on investigating the effects of prenatal organophosphate pesticide exposure on inhibitory cognitive control, and they found reduced activity in lateral frontal regions during a cognitive task. Given that the study primarily focused on activations and the motor inhibition regions, findings in the DMN could have been missed given that it was not a focus of the analysis and it is deactivated during this type of task.

Given the anatomic changes associated with pesticide exposure, the DMN is of particular concern in developing children exposed to pesticides. The DMN develops later in childhood and shows peak gray matter thickness between 8–10 years of age in typically developing children (Gogtay et al., 2004; Gogtay and Thompson 2010). These regions are among the latest cortical regions to develop along with the lateral frontal cortex. In addition, it has been shown that children’s cognitive performance is associated with the timing of cortical development in the medial frontal and lateral parietal regions (key components of the DMN). Those children with peak cortical thickness occurring between 9 and 12 years of age have the highest cognitive function (Shaw et al., 2006).

Analysis of brain networks obtained via resting-state functional MRI (rfMRI) data has shown great promise in characterizing abnormal changes in the brain that might eventually progress into cognitive abnormalities or neurological disorders across a range of populations (Bassett and Bullmore 2009). However, to our knowledge, no study has examined functional brain networks in children exposed to pesticides. We have previously shown the impacts of chronic exposure to pesticides on functional brain networks of adult male farmworkers (Bahrami et al., 2017), but studies in children are essential to assess the neurotoxic effects associated with being exposed to pesticides during active neurodevelopmental stages. Such work would complement the existing anatomical brain studies and could help lead to a better mechanistic understanding of the neurotoxic effects of pesticides in developing children.

The abnormal structural and functional brain findings coupled with the reports of deficits in cognitive development described above, necessitate further studies to identify mechanistic changes in brain function that are associated with pesticide exposure as noted in (Binter et al., 2020). Furthermore, the developmental trajectory of the DMN and the sensitivity of this circuit to prenatal pesticide exposure necessitates that this network be examined in pesticide-exposed children. Thus, the main objective of the current study was to compare functional brain networks between Latinx children from rural farmworker families (FW children) and urban Latinx children from non-farmworker families (NFW children) that are 8–9 years of age. The overarching premise was that the FW children will have a higher probability of being exposed to agricultural pesticides (as have been the focus of the vast majority of studies on neurodevelopment in children) and this would be associated with differences in brain network organization. Given that differences in brain networks can be locally restricted to particular subnetworks (Hallquist and Hillary 2019) or be global across the whole-brain connectome (Fornito et al., 2015), this study examined both local and global networks differences. The default mode network (DMN) was the focus of the local assessment, as it has been shown to have anatomical variations associated with prenatal pesticide exposure (Rauh et al., 2012; van den Dries et al. 2020). After comparing network topology between FW and NFW children, we evaluated associations between brain network topology and history of probable agricultural chemical exposure from prenatal to later childhood (96 months) periods. These secondary analyses were designed to see if higher exposure histories were associated with network differences comparable to those found in the main group comparison. Due to minimal historic exposure in NFW children, this analysis was restricted to FW children.

Methods

Study design and participants

The data for this study are part of a larger community-based participatory research project, called PACE5 (Preventing Agricultural Chemical Exposure – 5) that is being conducted in partnership between the North Carolina Farmworkers Project (Benson, NC; https://ncfwp.org/) and Wake Forest School of Medicine. This is a two-group, prospective study examining the health and cognitive effects of pesticide exposure for Latinx FW and NFW children. The study includes a sample of 76 children in rural, Latinx farmworker families (FW children) and a sample of 66 children in similar urban, Latinx families that do not include members employed as farmworkers (NFW children). The current analysis examines the baseline MRI data collected from a subsample of the entire cohort. Fifty-three of the FW children and 33 of the NFW children from the parent study agreed to participate; a total of 84 children, 52 FW children and 32 NFW children successfully completed the MRI protocol. The parents and children gave informed consent/assent to participate in the parent study and separately gave consent/assent to participate in the MRI study. The Wake Forest School of Medicine Institutional Review Board approved the PACE5 protocol and procedures. Participants received a $20 incentive for completing the baseline questionnaire and a $45 incentive for the MRI scan.

We put substantial effort into collecting data in both groups from children with similar socioeconomic status (SES) with respect to parental education, income, and occupation due to reported impact of SES on neurodevelopment (Hackman et al., 2010; Moriguchi and Shinohara 2019). Moreover, our multivariate modeling framework allowed for accounting for remaining differences in SES through incorporating them as confounding effects to further reduce the impact of SES on our results. Inclusion criteria for the children were similar for both the FW and NFW samples. All children had to be from families that self-identified as Latino or Hispanic and had household incomes below 200% of the US federal poverty line. Children had to be aged 8 years at baseline and had to have completed the first grade in the United States. For FW children, the mother or her partner living in the same home must have been employed in farm work on non-organic farms during the past three years (i.e., since the participating child was 5 y/o). For NFW children, adults living in the same home could not have been employed in an industry that involves routine exposure to pesticides (e.g., farm work, landscaping, and pest control) in the previous three years. Families of the NFW children could not have lived adjacent to agricultural fields in the previous three years. Children were excluded from the study if the primary language spoken in the home was other than Spanish or English; if the child had a life threatening illness, prior history of neurological conditions, physical condition or development disorder that would not allow them to complete or would interfere with the results of neurobehavioral tests; or MRIs; or refusal of mother/guardian to complete the questionnaires.

Participant recruitment

The rural FW families resided in eastern North Carolina, in the counties surrounding the town of Benson. The urban, NFW participants resided in central North Carolina, in the counties surrounding the city of Winston-Salem. For the FW children, North Carolina Farmworkers Project developed a list of families with an 8-year-old child, and the locations where they lived. This included families who had participated in a previous study of child growth and nutrition (Arcury et al., 2015) as well as families from the prior study of pesticide exposure in adult Latino farmworkers (Arcury et al., 2018). In addition, other community organizations that served FW families in the recruitment area were contacted for names of potential participants. For the NFW children, local recruiters and community members in Winston-Salem developed a list of families with an 8-year-old child and identified potential families at local community events. For both samples, mothers were contacted by a bilingual staff member who explained the overall study procedures and answered questions. For those parent/child dyads that agreed to participate, an in-person meeting was scheduled to obtain consent/assent and begin data collection. Because interviewers worked through community partners and attended community events to identify interested families, the number of potential participants or their parents who refused to participate is not known. All participants in the parent study were invited to participate in the MRI study, but not all agreed. We did not systematically assess reasons for refusing to participate in the MRI study, but time commitment and requirement to travel approximately two hours to the Wake Forest School of Medicine from eastern North Carolina were commonly noted.

Data collection

The primary measures used in the current study are functional brain images, but multiple other data have been used to compare the study populations and as covariates in the analyses. These additional data were collected using an interviewer-administered questionnaire assessing demographic and background data on the family and child completed by the child’s parent/guardian. Interviewers were native Spanish speakers; all spoke English, but with varying degrees of proficiency. They completed training before data collection began which included didactic instruction on the participant inclusion criteria, recruitment procedures, and questionnaire content. The interviewers completed audiorecorded practice interviews before the start of data collection. For most data collection the interviewers entered data in real time during the interviews using Research Electronic Data Capture (REDCap), a secure, web-based application for data collection (Harris et al., 2009).

For life history calendars used to assess historic pesticide exposure the data were collected in written form on long paper forms representing the child’s life span (for detailed methodology see (Quandt et al., 2020). For the main analyses presented here, months living near agricultural fields during pregnancy (0 prenatal months – 9 prenatal months), early childhood (birth – 35 months), and late childhood (36 months – 96 months) were used. The average months across the three exposure periods were used for each child. Since all three periods could be equally important, we first scaled the months of exposure within each period (centered and divided by standard deviation for that period), and then averaged across the three scaled values for each child.

We controlled for several demographic variables that were different between the FW and NFW children. These variables included: maternal education, maternal race, childhood preschool education status, and parents’ occupation. Childhood preschool education status indicated whether the child received early childhood education to aid in the cognitive and social development before they entered school, at places such as Migrant Head Start. Any history of a diagnosis of a learning disability in the child was assessed using a single question completed by their mothers. There were only 3 positive responses (1 rural child and 2 urban children) precluding useful incorporation into the statistical analyses. Maternal substance use during pregnancy (tobacco, alcohol, marijuana) was found in only 1 urban mother (alcohol); thus, they were not included in the statistical analyses.

Blood lead was measured through collecting finger stick samples. The samples were collected at the baseline visit using filter paper blood kits provided by Laboratory Corporation of America® Holdings (Labcorp.com, Test 79,128/CPT 83,655). After allowing a sample to dry completely, they were bagged and held at room temperature until submitted to Labcorp for analysis. The samples were assessed by Labcorp in micrograms per deciliter (μg/dl). Any sample with inadequate blood collection was recollected at the participant’s next study follow-up visit.

We should note that this study is specifically focused on using brain imaging as an independent tool to assess early brain changes associated with pesticide exposure in FW children. Such early changes, which predispose this vulnerable population at additional risk for developing cognitive disorders in later life, may not be necessarily reflected by conventional cognitive assessments. Although, as described above, we compared several other measures, such as child or maternal language, to control for their possible influence on our results, IQ and other cognitive assessments will be studied comprehensively in another independent part of this prospective project, and ultimately in combination with brain imaging data when the cognitive assessment data as well as longitudinal data from a 2–3 years follow up become available.

MRI acquisition

Study participants from FW families were scheduled for MRI scans in small groups (~4 children) and were transported to Wake Forest School of Medicine together with a parent in a large passenger van. These scans all took place on the weekends. The children from NFW families were from local urban areas and were scheduled individually for MRI scans based on their availability. Many of the local children were also scanned on weekends along with the groups of children from FW families. Some children also completed cognitive testing as part of the larger study on the same day as the MRI scan.

To reduce head motion and enhance data usability, each participant was familiarized with the scan procedure and environment prior to the MRI scan session. It has been shown that high success rates for MRI scans in children can be achieved using a home-built mock MRI scanner (Barnea-Goraly et al., 2014). We utilized a home-built mock scanner that also included a 3-D printed mock head coil to further replicate the real environment. As is done in the real MRI scanner, children were equipped with earbud headphones and were positioned in the head coil so that they could see a computer screen through a reverse-projection mirror. Children were instructed that they had to remain still and not move their head. Scanner noise was played through external speakers and the children were allowed to watch a cartoon on the computer screen during the simulation. Study staff monitored the participant and stopped the simulation when head motion or substantive body motion was visually apparent. The child was re-instructed about the importance of remaining still and the simulation was restarted. After the child was able to lie still for 5 consecutive minutes the simulation was ended, the child was praised for their performance, and instructed to do the same in the real MRI scanner. One child withdrew from the MRI scan portion of the study after the mock scanner experience. All other enrolled children completed the real scan session that occurred as close in time as possible to the simulation, but never more than about 60 min later.

The MRI scan was planned with the knowledge that children can have difficulty lying still for long periods in the MRI scanner. Scan sessions were planned for 30 min, with primary image sequence first (resting-state fMRI). During the functional scan the participants were instructed to focus on a fixation cross presented in the middle of the screen. The MRI technician evaluated the functional scan in real time to assess for excessive motion. In scans with visible motion, the functional scan was repeated. One child would not lie still during repeated attempts to collect the functional scan so the session was terminated. Once the resting-state scan was completed, a cartoon was played for the child to watch during the collection of structural images collection. Complete functional and high-resolution structural images were collected on 84 children.

All imaging was completed on a single Seimens 3T Magnetom Skyra scanner with a 32-channel head coil. The functional imaging used a whole-brain gradient echo-planar imaging (EPI) sequence to measure the blood-oxygenation level-dependent (BOLD) contrast (Ogawa et al., 1990). Images were acquired in the transverse plane with the following parameters: TR: 2000 ms, number of volumes: 157 (with 35 interleaved slices per volume), and spatial resolution: 4.0 mm × 4.0 mm × 5.0 mm. The scan duration was 5 min 14 s. Single-shot 3D MPRAGE GRAPPA2 sequence was used to acquire the T1-weighted structural data in the sagittal plane with the following parameters: TE: 2.99 ms, TR: 2.32 s, and number of slices: 192. The image was reformatted into the transverse plane at the scanner prior to exporting the data for processing.

MRI preprocessing

Image analyses were performed using Statistical Parametric Mapping 12 (SPM 12, Welcome Trust Center, London, UK: www.fil.ion.ucl.ac.uk/spm/), Advanced Normalization Tools (ANTs, (Avants et al., 2011)), and in-house Matlab scripts. The preprocessing included: i) segmenting the T1-wieghted structural images based on standard 6 tissue priors using SPM12; ii) combining gray matter (GM) and white matter (WM) tissue segments to create a brain tissue mask and the mask was thresholded at >50% and applied to the T1 image to create a masked T1; iii) manual cleaning of the masked T1 in MRIcron software (https://www.nitrc.org/projects/mricron) to correct any misclassified voxels; iv) normalizing the manually cleaned masked T1 to MNI standard space (www.mni.mcgill.ca) using ANTs. Several studies have introduced age-specific templates which could reduce the normalization error, specifically for studies involving children (Fonov et al., 2011; Sanchez et al., 2012). However, for consistency and to allow to better compare and contrast our results with other studies and the longitudinal studies in the future, we used a combination of a strong normalization tool (Avants et al., 2011) and the adult MNI template, which has been used at a substantially higher rate; v) transforming the Shen functional brain parcelation atlas (Shen et al., 2013) to each subject’s native space using the inverse transform to MNI space; vi) removing the first 10 vol for the functional scan; vii) realigning the functional scans to the first volume; viii) correcting for slice-time differences; ix) running AROMA motion correction (Pruim et al., 2015); x) filtering the functional data (voxel time series) with a standard band pass filter (0.009–0.08 Hz) to remove low frequency drift and physiological noise; xi) regressing out the motion parameters with 6 dof, within atlas mean signal, mean white mater (WM), and mean cerebral spinal fluid (CSF) signals; and xii) averaging time series within each of the 268 atlas regions. As the atlas was warped to each individual participant, all analyses were performed using fMRI data in the participants’ native brain space. Quality control of the preprocessing included: i) manually correcting misclassification of tissue segments for improved warping; ii) visually inspecting all warps and affine alignments; and iii) visually inspecting the results of the AROMA motion correction and excluding any scans that did not meet our rigorous standards. In total, 6 children, including: 4 FW and 2 NFW, were removed from subsequent analyses due to motion in resting state scans that could not be corrected using AROMA. Each scan with severe movement had at least three study team members check the output prior to exclusion. Supplemental Figure S1 shows the motion and time series data from an example child that was removed after failing to correct the excessive motion. The motion and time series data from a child who was retained after successful motion correction is shown in supplemental Figure S2.

A functional network for each participant was then generated by computing the Pearson (full) correlation between average time series from 268 brain regions of interest (ROIs) contained in the Shen atlas. The spatial distance between brain ROIs (used as an endogenous confounding effect (Friedman et al., 2014)) was computed using the functional atlas in each child’s native space. Negative correlation values were set to zero for the subsequent analyses due to the challenge of defining several network metrics, such as: clustering coefficient, in networks containing negative edges (Fraiman et al., 2009; Telesford et al., 2011). In addition, the neurobiological interpretations of positive and negative edges are very different, and the distributions of network variables (such as degree) are different for positive and negative edges (Parente et al., 2018; Schwarz and McGonigle 2011). Although negative correlations are not regularly used in network neuroscience, if they are used, positive and negative networks should be generated and assessed separately (Schwarz and McGonigle 2011).

We used a mask to perform region-specific assessments for DMN connectivity. The mask defining the DMN was made using 37 ROIs from the Shen atlas that were identified as located within DMN. These ROIs were identified by intersecting the Shen atlas with a functionally defined voxel-wise DMN mask generated from a prior study (Mayhugh et al., 2016). Any region in the Shen atlas with at least 50% of its voxels in the DMN mask was included in the final DMN mask. Supplemental Figure S3 shows the final DMN mask, and the 37 ROIs have been listed in a separately uploaded excel file named DMN_ROIs_Info.

Statistical analysis

Modeling framework

A two-part mixed-effects regression framework (Bahrami et al., 2019a; Simpson and Laurienti 2015) was used to simultaneously study group differences in the DMN and across all other regions of the brain, and to assess how history of childhood exposure is associated with brain network topology in the FW children. This regression framework was used to examine if and how desired covariates, such as: the grouping covariate separating FW and NFW children or the covariate representing historic childhood exposure to pesticides is related to the brain connectivity and topology in local (DMN) and global regions of the brain. This multivariate model also allowed for the inclusion of possible confounding covariates, including: maternal education, child preschool education status, and the spatial distance between brain regions as a potentially endogenous network-related confounding effect. The desired and confounding covariates were used as independent variables and were related to the probability (presence/absence) and strength of present brain connections, as dependent variables, in two separate models (part I: probability, part II: strength). Network metrics (clustering coefficient and global efficiency) were also included as independent covariates in both models, and their interactions with the main covariates of interest were used to find the topological differences between the two groups and to examine the effects of childhood exposure history on the brain network topology – i.e., the interactions showed if/how the relationship between brain connectivity (probability/strength) and network metrics (characterizing topology) was modified by the grouping covariate (FW vs NFW) or the childhood exposure history.

More specifically, let Rijk denote a binary variable which is one if the correlation value between node j and node k of the ith child’s network is positive, and zero otherwise. This variable specifies whether a connection exists between node j and node k of the ith child’s network. Also, let Yijk denote a continuous variable for positive correlation values (present connections) between node j and node k of the ith child’s network. Then, we can define the following conditional probabilities:

Where pijk is the probability of having a connection between node j and node k of the ith subject’s network, βr is the fixed effects (population) parameters vector (Note that the original paper which introduces this new extension uses θ for fixed effects parameters. We used β as it is a more familiar notation for fixed effect parameters.), and bri is the random-effects parameter vector modeling correlation (dependence) between repeated network features for subject i. If Sijk denotes a variable specifying present connections between node j and node k of the ith child’s network (i.e., Sijk = [Yijk ∣ Rijk = 1]), we can define the following two-part mixed modeling framework for the probability and strength of brain connections:

Where Tijk and Zijk are design matrices for the fixed- and random-effects, respectively, eijk captures the random noise (not captured by random effects) in the connection strength between node j and node k of the ith child’s network, and βs and bsi are analogous parameters to βr and bri, respectively, but for the connection strength. Eq. 3 is a logistic regression model that quantifies the relationship between the connection probability and sets of desired covariates. FZT is the Fisher’s Z-transform applied to ensure that normality assumption is met. By including network metrics (such as clustering coefficient and global efficiency) and regional covariates modeling brain subnetworks (DMN in this study), this framework allows testing hypotheses about the DMN and whole brain connectivity and network metrics and assessing how covariates of interest affect the connectivity and topology at local and global levels. This framework allows to simultaneously study changes in the DMN and across all other regions of the brain. For more detail, see the referenced paper. We used the WFU_MMNET toolbox (Bahrami et al., 2019b) and in-house Matlab scripts to generate the appropriate datasets for this modeling framework, and used SAS v.9.4 to fit the statistical models.

The differences in the DMN between the two groups and the effects of childhood exposure on this subnetwork were obtained by applying appropriate contrast statements (estimates of linear combinations of appropriate covariates) on already estimated residuals in post-hoc analyses (This will be explained more in the Planned post-hoc analysis (2.6.3) section). Two separate analyses were run with each analysis using the two-part model for the probability and strength of brain connections. The primary analysis compared the FW and NFW children to test for hypothesized group differences. This was followed by the secondary analysis examining associations between brain network topology and history of pesticide exposure quantified by months living next to agricultural fields (averaged across the three mentioned periods), based on prior work in these children (Quandt et al., 2020). The latter analysis was restricted to FW children as exposure values were nearly collinear with the FW vs. NFW classifications if all children were used.

Covariates

The utilized (independent) covariates are described in the following and are summarized in Table 1.

Table 1

Independent variables used for mixed-effects regression analyses.

Covariates	Parameters
	Probability Model	Strength Model
Covariate of Interest	β r,COI	β s,COI
Subnetwork covariate	β r,DMN	β s,DMN
Network metric covariates	β r,Clust	β s,Clust
	β r,Eglob	β s,Eglob
Interaction covariates	β r,COI×DMN	β s,COI×DMN
	β r,COI×Clust	β s,COI×Clust
	β r,COI×Eglob	β s,COI×Eglob
	β r,DMN×Clust	β s,DMN×Clust
	β r,DMN×Eglob	β s,DMN×Eglob
	β r,COI×DMN×Clust	β s,COI×DMN×Clust
	β r,COI×DMN×Eglob	β s,COI×DMN×Eglob
Confounding covariates	β r,Dist	β s,Dist
	β r,Dist2	β s,Dist2
	β r,MatRace	β s,MatRace
	β r,PreEdu	β s,PreEdu
	β r,MatEdu	β s,MatEdu
	β r,ParOccu	β s,ParOccu
Contrast Statements*	βr,COI×Clust + βr,COI×DMN×Clust	βs,COI×Clust + βs,COI×DMN×Clust
	βr,COI×Eglob + βr,COI×DMN×Eglob	βs,COI×Eglob + βs,COI×DMN×Eglob

Contrast statements were not used as additional independent variables. They were rather used in post-hoc analyses to test hypotheses (i.e., obtain inference) on combinations of estimated parameters using their already estimated residuals.

The r and s subscripts in this table and all subsequent tables are used for the probability and strength models, respectively. The following covariates were included: i) Covariate of interest: For the first analysis, this covariate was a binary variable separating FW (labeled as one) and NFW (labeled as zero) children. For the secondary analysis, this covariate was a continuous variable representing the average historic pesticide exposure from prenatal to 96 month old; ii) Subnetwork covariate: A binary variable separating the connections within DMN (labeled as one) from other connections in the brain as well as connections from DMN to other brain regions (all labeled as zero); iii) Network metric covariates: The average of clustering coefficient and global efficiency in each dyad (i.e., nodal pair); iv) Interaction covariates: The two-way interactions of the covariates described in (i) - covariate of interest and (ii) - subnetwork covariate, the two-way interactions of the covariates described in (i) and (iii) - network metrics covariates, the two-way interactions of the covariates described in (ii) and (iii), and the three-way interactions described in (i), (ii), and (iii). Incorporating all of these two-way and three-way interactions was required to test desired hypotheses directly or via contrast statements in post-hoc analyses.; and v) Confounding covariates: spatial distance and square of spatial distance between ROIs, childhood preschool education status – a binary variable specifying whether the child received (labeled as one) or did not receive (labeled as zero) education; maternal education – a categorical variable with three levels, including; 0–6 (level 1), 7–12 (level 2), and ≥13 (level 3) years of education with level 3 being used as the reference category, maternal race – a categorical variable with three categories, including: white, mixed, and other (African-American, Hispanic, etc.) with White being used as the reference category, and parents’ occupation – a categorical variable with 5 categories, including: farm work, construction, cleaning, two of these three occupations (e.g., father: construction, mother: cleaning), and other (none of these three) with two occupations being used as the reference category. These covariates and the parameters representing them in our models are further described in Table 1. Random effects were used for the connection probability, connection strength, clustering coefficient, global efficiency, spatial distance and square of spatial distance to capture nodal dependence for each individual. An unstructured variance covariance matrix for the random effects in terms of their Cholesky-root factors (Gentle 1998) was used.

Planned post-hoc analysis

To find significant differences or modification effects in the DMN, further post-hoc analyses were required. Specifically, we estimated the appropriate contrast statements of already estimated residuals for corresponding parameters, to obtain inference about if/how the relationship between network metrics and brain connectivity within the DMN is modified by the grouping covariate (i.e., if the relationship is different between FW and NFW children in the DMN) and if/how this relationship in the DMN is modified by the history of childhood exposure to pesticides. The contrast statements for continuous variables (as for our secondary analysis) are made through a unit change in that variable. The contrast statements are provided in Table 1, but, again note that these statements (i.e., combination of parameters) were not used as additional independent variables, and rather were used to test hypotheses (i.e., obtain inference) on combinations of estimated parameters by using their already estimated residuals. A summary of what each important interaction (i.e., interactions that include network metrics and the covariate-of-interest, or COI) shows, and how contrast statements for comparing the DMN are obtained via the interaction covariates are in the supplemental material. For more detail about comparing subnetworks through contrast statements within this framework please see the supplemental material and the referenced paper (Bahrami et al., 2019a).

Representative group network generation

Due to the non-Euclidean nature of complex brain networks, it is not possible to generate a group mean network without affecting the overall network topology (Hayasaka and Laurienti 2010). In other words, simply averaging connectivity matrices across study participants results in a single group network with topological properties that are not representative of the typical participant (Simpson et al., 2012). To generate a meaningful visualizations of the significant results from the current study, representative group networks were estimated with respect to the observed modification effects of the covariate of interest on the network variable – connectivity associations. Group networks were estimated from the fitted probability model since the strength model did not yield significant results (explained in more detail in the Results). For the first analysis where the COI was a binary variable separating FW and NFW children, a network was estimated for each group, and for the secondary analysis where the COI was a continuous variable (months of exposure), the representative group networks were estimated for the group of subjects with minimum and maximum exposure values. Using average network variables for each group and the estimated parameters along with their standard deviation, the probability of each edge existing (i.e., pijk) was derived which yielded a 35,778 × 1 (vectorized symmetric network with 268 nodes) vector of edge-wise probability values. We then simulated the existence of edges (presence/absence) for all 35,778 edges from a Bernouli distribution. To avoid potential bias from the simulation step for each edge and for each group, 20 probability values were derived (from a normal distribution using estimated parameters and their standard errors), and this whole process was repeated 20 further times and averaged to generate each final edge probability. The representative network was then thresholded to obtain the binary vector representing the presence/absence of edges (1 if >0.5 and 0 otherwise) and finally multiplied by the edge-wise probability vector (averaged across the 400 (20 × 20) realizations) to get a single probability-based weighted network. We only used covariates from the probability model that contributed to the observed modification effect of the COI for each analysis in estimating the representative group networks. A diagram is provided in the supplemental material (Figure S5) which demonstrates each step for estimating the group networks from the fitted models. All estimated networks were only used to illustrate the significant findings from the statistical analysis. The simulated group networks were not used to make inferences related to the hypotheses being tested.

Results

Study cohort

A total of 84 children (52 FW and 32 NFW) completed functional and structural scans. However, as detailed above, 6 children were removed due to excessive head motion. The final data included 78 children, 48 FW and 30 NFW children. Table 2 shows the demographic information of the 78 children used in this paper with comparisons between the two groups. As Table 2 (last column) shows, the two groups are different in childhood preschool education status, maternal education, maternal race, and parents’ occupation. We thus controlled for these variables by including them as confounding covariates in all subsequent analyses to reduce potential bias in our results. Also, as supplemental Table S5 shows, the MRI sample used here represents the parent study population in every demographic variable.

Table 2

Study cohort.

Variable	FW Children (N = 48)	NFW Children (N = 30)	p-value
Child Age (yrs)			0.1673
mean±std	8.34±0.29	8.44±0.35
Range	8.01 – 8.99	8.01 – 9.00
Child Birth Weight (lbs)	6.94±1.36	6.82±1.66	0.7455
Child Gender			0.9572
Male	24	14
Female	24	16
Childhood Preschool Education			0.0123
Yes	11	16
No	37	14
Childhood Learning Disability			0.6753
Yes	1	2
No	47	28
* Childhood Blood Lead Levels			0.5018
< 1 ug/dl[#]	16	11
= 1 & < 2 ug/dl	29	11
≥ 2 ug/dl	3	2
Child Birth Country			0.8441
United States	42	28
Mexico	4	1
Other	2	1
Child Language			0.6918
English	48	30
English & other[#]	46	30
Maternal Education (yrs)			0.0451
0 – 6	22	6
7 – 12	23	19
≥ 13	3	5
Maternal Birth Country			0.7947
United States	3	3
Mexico	41	24
Other	4	3
Maternal Race			0.0004
White	31	12
Mixed	17	1
Other	0	17
Maternal Language			0.3294
Spanish	48	30
Spanish & other[#]	27	21
Maternal Alcohol Consumption			0.8113
Yes	0	1
No	48	29
Maternal Tobacco Consumption			1.0000
Yes	0	0
No	48	30
Maternal Marijuana Consumption			1.0000
Yes	0	0
No	48	30
Parents’ Occupation			0.0005
Farm work	22	0
Construction	5	16
Cleaning	4	4
Two of above	14	0
Other	3	10

FW children: Children from rural farmworker families, NFW children: Children from urban non-farmworker families. Bold values show significant differences.

the measurements for blood lead levels were missing for 22 children.

ug/dl: micrograms per deciliter

English & other: English and Spanish for almost everyone in this language group.

Spanish & other: e.g., Spanish and English for almost everyone in this language group.

Comparing FW and NFW children

Table 3 presents the key results for our first analysis: comparing FW and NFW children. Parameter estimates and (adaptive) False Discovery Rate (FDR)-corrected p-values of the important covariates for the probability (Eq. 3) and strength (Eq. 4) models are shown in this table (Full results are shown in Table S6). Below we briefly explain a summary of the essential findings both within the DMN (local) and across all other brain regions (global). The results for other parameters are fully explained in the supplemental material (see Comparing FW and NFW Children – Full Results in supplement).

Table 3

FW vs NFW Children.

Probability Model Outputs			Strength Model Outputs
Parameter	Estimate	*p-value	Parameter	Estimate	*p-value
β r,COI×Clust	0.0435	0.1600	β s,COI×Clust	−0.0017	0.7501
β r,COI×Eglob	−0.0265	0.3588	β s,COI×Eglob	0.0042	0.4375
β r,COI×DMN×Clust	0.4790	<0.0001	β s,COI×DMN×Clust	0.0003	0.9422
β r, COI×DMN×Eglob	−0.4793	<0.0001	β s, COI×DMN×Eglob	0.0001	0.9945
βr,COI×Clust + βr,COI×DMN×Clust	0.5226	<0.0001	βs,COI×Clust + βs,COI×DMN×Clust	−0.0014	0.7752
βr,COI×Eglob + βr, COI×DMN×Eglob	−0.5058	<0.0001	βs,COI×Eglob + βs,COI×DMN×Eglob	0.0042	0.6167

COI: A binary variable separating FW and NFW children.

Adjusted using the adaptive FDR procedure described in (Benjamini et al. 2000). Bold values show significant COI – related inferential results.

βr,COI×Clust, βs,COI×Clust: The relationship between the clustering coefficient and connection probability/strength across the entire brain excluding the DMN is not significantly different between FW and NFW children.

βr,COI×Eglob, βs, COI×Eglob : The relationship between the global efficiency and connection probability/strength across the entire brain excluding the DMN is not significantly different between FW and NFW children.

(βr,COI×Clust + βr,COI×DMN×Clust), (βs,COI×Clust + βs,COI×DMN×Clust): The relationship between the clustering coefficient and connection probability within the DMN is significantly different between FW and NFW children. However, the same relationship between the clustering coefficient and connection strength is not significantly different. Within the DMN, FW children are more likely (have higher probability) to have connections between regions with higher clustering coefficient when compared to NFW children (positive estimate/slope).

(βr,COI×Eglob + βr,COI×DMN×Eglob), (βs,COI×Eglob + βs,COI×DMN×Eglob): The relationship between the global efficiency and connection probability within the DMN is significantly different between FW and NFW children. However, the same relationship between the global efficiency and connection strength is not significantly different. Within the DMN, FW children are less likely (have lower probability) to have connections between regions with higher global efficiency when compared to NFW children (negative estimate/slope).

In Fig. 1, the findings from the statistical model detailed above have been transformed back into brain space to illustrate differences in the topology of the DMN. Fig. 1 also shows the line plots for the connection probability – clustering/efficiency for the two groups along with their 95% confidence intervals. Nodes in Fig. 1.A are sized by their actual clustering values and colored by the sum of their connection probability – clustering slopes, and in Fig. 1.B, they are sized by their actual global efficiency values and colored by the sum of their connection probability – global efficiency slopes. Fig. 1.A shows that in FW children, the DMN nodes have an overall greater probability of being connected as clustering increases (red nodes). This is consistent with the line plots shown on the left and the overall model results where the probability/clustering slope difference between the two groups (FW - NFW) was positive. In NFW children, the probability of the DMN nodes to be connected to other nodes is essentially independent of their clustering values (green nodes) indicating no relationship between connection probability and clustering. This is also clearly demonstrated by the line plot. Fig. 1.B shows that in FW children, the probability of DMN nodes to be connected to other nodes is nearly independent of their global efficiency values (teal nodes), indicating no relationship between connection probability and global efficiency, as clearly demonstrated by the line plot on the left. In NFW children however, the DMN nodes have an overall greater probability of being connected as global efficiency increases (red nodes) as found in the statistical model results and the line plot in the left. The same color scale has been used for all networks in this Figure (and Fig. 2) to better contrast the similarities and differences.

Fig. 1.

Visualization of differences in DMN organization between FW and NFW children. This figure (including line plots and brain networks) was created using the coefficients from the probability model in Table 3. Line plots along with their confidence intervals (magenta shading) and representative group networks mapped back into the brain space are shown. Plots of the relationship between the network metrics and connection probability demonstrate the significant group differences in the slopes for clustering (A) and global efficiency (B) within the DMN identified using the statistical models. For connection probability/clustering, FW children have a significantly steeper slope compared with NFW children, while for the connection probability/global efficiency, a reverse pattern is observed. This is further illustrated at the level of individual nodes in the DMN using representative group networks. Nodes are colored by the sum of their connection probability/clustering slopes (A) and connection probability/global efficiency slopes (B) with the same color scale. Also, nodes are sized by their actual clustering (A) and global efficiency (B) values (scaled independently for each network) to aid in interpreting the slope differences between the two groups. * The y-axis in this figure and Fig. 2 is the log odds of connection probability but the axis is labeled as connection probability for simplicity. Also, the min and max values in the y-axis are the same for both line plots to better contrast the differences.

Fig. 2.

Visualization of the associations between historic pesticide exposure and DMN organization in FW children. The surface plots show how the relationships of connection probability/clustering coefficient (A) and connection probability/global efficiency (B) are associated with history of childhood exposure to pesticides. Increasing levels of pesticide exposure is associated with a more positive relationship between connection probability and clustering coefficient, but with a more negative relationship between connection probability and global efficiency in the DMN. Group representative networks are shown for the subjects with minimum (dashed line on the surface) and maximum (solid line on the surface) exposure values for each group (For node color and size, see Fig. 1 caption.). As this figure clearly demonstrates, with respect to both clustering and efficiency in the DMN, the brain network of FW children with low exposure is similar to the brain network of NFW children (in Fig. 1), and the brain network of FW children with high exposure is similar to the brain network of FW children (in Fig. 1). We have used the same color scale for all networks shown in this Figure and Fig. 1 as the color bars show to better contrast the differences and similarities.

History of pesticide exposure in FW children

In our secondary analysis, using data from FW children only, we examined the association of a history of childhood exposure to pesticides (averaged across pregnancy through 96 months of age) with brain network topology. Table 4 presents the key results (estimates and (adaptive) False Discovery Rate (FDR)-corrected p-values) for the probability and strength models (Full results are shown in supplemental Table S7). A brief description of important inferential results identifying network topological properties associated with childhood exposure to pesticides in the DMN (local) and across other brain regions (global) are presented below. The estimates for other parameters are fully explained in the supplemental material (see History of Pesticide Exposure in Children – Full Results in supplement).

Table 4

Childhood exposure to pesticides among FW.

Probability Model Outputs			Strength Model Outputs
Parameter	Estimate	*p-value	Parameter	Estimate	*p-value
β r,COI×Clust	−0.0160	0.3956	β s,COI×Clust	0.0019	0.4330
β r,COI×Eglob	−0.0009	0.9598	β s,COI×Eglob	−0.0005	0.8009
β r,COI×DMN×Clust	0.1444	<0.0001	β s,COI×DMN×Clust	−0.0030	0.3566
β r,COI×DMN×Eglob	−0.0705	0.0304	β s,COI×DMN×Eglob	0.0009	0.7818
βr,COI×Clust+ βr,COI×DMN×Clust	0.1284	0.0003	βs,COI×Clust+ βs,COI×DMN×Clust	−0.0010	0.7504
βr,COI×Eglob+ βr, COI×DMN×Eglob	−0.0759	0.0466	βs,COIβEglob+ βs, COI×DMN×Eglob	0.0003	0.9269

COI: A continuous variable representing months of exposure.

Adjusted using the adaptive FDR procedure described in (Benjamini et al. 2000). Bold values show significant COI – related inferential results.

βr,COI×Clust, βs,COI×Clust: The relationship between the clustering coefficient and connection probability/strength across the entire brain excluding the DMN is not affected by a history of childhood exposure to pesticides in FW children.

βr,COI×Eglob, βs,COI×Eglob: The relationship between the global efficiency and connection probability/strength across the entire brain excluding the DMN is not affected by a history of childhood exposure to pesticides in FW children.

(βr,COI×Clust + βr,COI×DMN×Clust), (βs,COI×Clust + βs,COI×DMN×Clust): The relationship between the clustering coefficient and connection probability within the DMN is significantly affected by the history of childhood exposure to pesticides. However, the same relationship between clustering coefficient and connection strength is not affected by the history of childhood exposure to pesticides. Within the DMN, higher childhood exposure to pesticides is associated with a greater likelihood (higher probability) of having connections between regions with higher clustering coefficient.

(βr,COI×Eglob + βr,COI×DMN×Eglob), (βs,COI×Eglob + βs,COI×DMN×Eglob): The relationship between the global efficiency and connection probability within the DMN is significantly affected by the history of childhood exposure to pesticides. However, the same relationship between global efficiency and connection strength is not affected by the history of childhood exposure to pesticides. Within the DMN, higher childhood exposure to pesticides is associated with a lower likelihood (lower probability) of having connections between regions with higher global efficiency.

The significant findings described above are illustrated by the surface plots and representative networks in Fig. 2. The surfaces are colored by the slope of connection probability/clustering relationship (Fig. 2.A), and connection probability/global efficiency relationship (Fig. 2.B) at each exposure value. As the surface plots in Fig. 2 clearly illustrate, the relationship between connection probability and clustering increases as the childhood exposure increases. For the relationship between connection probability and global efficiency, a reversed pattern is observed for increasing the childhood exposure to pesticides. The brain space figures illustrate relationships between history of pesticide exposure and topology of the DMN. Since the exposure was modeled as a continuous variable, for illustrative purposes, the relationship between the connection probability in the DMN and network metrics (clustering coefficient (A) and global efficiency (B)), are shown for minimum (dashed black line on the surface) and maximum (solid black line on the surface) exposure values. As in Fig. 1, nodes in Fig. 2 are colored by the sum of their connection probability – clustering (A) and connection probability – global efficiency (B) slopes.

The findings in FW children were qualitatively compared to the differences observed between FW and NFW children (i.e., the differences found from our previous analysis). Comparisons of the brain maps between this analysis using historic pesticide exposure measures in FW children with the differences found between these children and those from NFW children (summarized in Table 5) are highly consistent. This indicates that the topology of the DMN in FW children with low historic pesticide exposure is more similar to the topology of the DMN in NFW children than FW children with high historic exposure. Additional supporting analyses are presented in the supplemental materials using historic exposure assessments from three different time periods: pregnancy (0 prenatal months – 9 prenatal months), early childhood (birth – 35 months), and late childhood (36 months – 96 months).

Table 5

Summary of Results from Tables 3 and 4.

	FW vs NFW	Exposure effect in FW
Connection Probability - Clustering relation in DMN	* more positive for FW than NFW	* more positive for FW with higher exposure
Connection Probability - Efficiency relation in DMN	* more positive for NFW than FW	* more positive for FW with lower exposure
Connection Probability - Clustering relation in other regions	No difference	No effect
Connection Probability - Efficiency relation in other regions	No difference	No effect
Connection Strength - Clustering relation in DMN	No difference	No effect
Connection Strength - Efficiency relation in DMN	No difference	No effect
Connection Strength - Clustering relation in other regions	No difference	No effect
Connection Strength - Efficiency relation in other regions	No difference	No effect

Significant relationship (adjPval < 0.05).

Discussion

Children who live in agricultural communities are exposed to pesticides through multiple pathways. While low-level acute exposure to pesticides may only cause minor immediate implications, such as rash, nausea, and vomiting, this exposure can be associated with serious long-term health consequences within this vulnerable population (Bouchard et al., 2011; Nichols et al., 2014; Quandt and Arnold 2020). A major consequence for concern with childhood pesticide exposure is the potential impact on neurodevelopment. Prenatal and childhood exposure to pesticides have been associated with adverse neurological consequences (de Joode et al. 2016; Munoz-Quezada et al., 2016; Zhang et al., 2019). However, the data needed to provide a more profound understanding of possible effects on brain physiology are still limited, making further investigation critical.

In this study, we compared connective/topology of functional brain networks between FW and NFW children. After making the group comparisons, we studied associations between childhood historic exposure and brain network organization. We examined both global changes across the brain as well as local changes in the DMN, a critical subnetwork in the brain (Buckner et al., 2008; Smallwood et al., 2021). We used a state-of-the-art modeling framework that allowed modeling the DMN within the context of the whole brain rather than isolating it and treating it as an independent network (Bahrami et al., 2019a). This model also allowed controlling for important sources of confounding effects such as maternal education and race, childhood preschool education, etc. The reliability for studies using mixed-effects regression models, which already incorporate the within-subject correlations, is often assessed with respect to goodness of fit measures. Since the most appropriate way to assess the goodness of fit of statistical models in the network context is through simulation analyses (Hunter et al., 2008), we did simulation analyses to evaluate the performance of our model fits. For more detail, see supplemental section Reliability Analysis and Table S11 in the supplemental file.

The overarching goal of this study was to test the hypothesis that brain connectivity/topology is different between FW and NFW children, and more importantly, to test the hypothesis that pesticide exposure was associated with any observed difference. Results from our mixed regression framework, summarized in Tables 3 and 4 (as well as supplemental tables S6 and S7 for full results), raise serious concerns about the impact of pesticides on the studied children. To aid in understanding and interpreting our results, we used Figs. 1 and 2 to better illustrate them. While some care is needed in interpreting the results as this methodology is relatively new and the complex network topology cannot always be easily summarized by a single value, our statistical modeling framework allowed for identifying more complex changes in the brain among children that might ultimately progress into more serious abnormalities.

Our results demonstrated significant topological differences between FW and NFW children, which were identified in the DMN but not the remainder of the brain. In FW children, connectivity of the DMN was primarily associated with clustering coefficient with connection probability increasing as clustering increased. However, the connectivity of the DMN in this group of children had essentially no relationship to global efficiency. Thus, in children from farmworker families, the DMN is made up of a core of highly clustered nodes that are highly interconnected. For NFW children, the opposite pattern was observed. Connections in the DMN of these children were primarily associated with global efficiency, with increased connection probability being associated with increased global efficiency. The probability that nodes within the DMN were connected had no relationship with clustering in NFW children. This pattern of connectivity indicates that the DMN in NFW children is composed of a core of highly interconnected nodes that also have high global efficiency.

Clustering coefficient and global efficiency are key measures of segregation (presence of highly interconnected regions supporting regional specialization) and integration (wide-spread connectivity interconnecting specialized regions) in the brain, respectively (Rubinov and Sporns 2010). Studies of structural and functional brain networks have shown that the balance between segregation and integration of neural information processing in the brain is essential for normal cognitive processes (Bullmore and Sporns 2009; Sporns 2013). Numerous studies over the past two decades have shown disruptions of this architecture, affecting segregation, integration, or both, in neurological disorders such as Alzheimer’s disease (Dai et al., 2019; Sanz-Arigita et al., 2010; Zhao et al., 2012), schizophrenia (He et al., 2012; Liu et al., 2008), and Parkinson’s disease (Herrington et al., 2017; Wu et al., 2018). Alterations of clustering coefficient and global efficiency within the DMN, which disrupts the optimal communication balance, has been shown to be associated with adverse neurological effects (Akiki et al., 2018; Garrity et al., 2007; Hafkemeijer et al., 2012). Our study indicates that the topology of the DMN is different between FW and NFW children. The dominance of clustering on connectivity in FW children and the dominance of global efficiency on connectivity in the NFW children suggests differential segregation and integration in these groups. Given the complex association between network topology and information processing, further work is clearly warranted to understand how these group differences may be associated with differences in cognitive processing.

The cross-sectional nature of this two-group comparisons precludes us from concluding why the groups differ, and from interpreting the findings in FW children as detrimental compared to the NFW children. However, a critical additional finding of this study concerns the association between historic childhood pesticide exposure and brain network topology in FW children. We found that histories of low exposure was associated with brain networks that resemble the NFW children. Whereas, the children with high exposure more closely resembled the FW children as a group. This pattern of low exposure resembling NFW children and high exposure resembling FW children lends support to the hypothesis that the FW children have altered DMN organization due to adverse effects of exposure to neurotoxic pesticides on neurodevelopment in these children who experience chronic exposure to pesticides (Sapbamrer and Hongsibsong 2019). The ultimate implications of the observed DMN connectivity, which is engaged across a diverse set of cognitive processes (Broyd et al., 2009; Buckner et al., 2008), particularly higher-level cognitive processes (Smallwood et al., 2021), will require further study. Much of the maturation of this circuit occurs later and parallels maturation of high-order cognitive and emotional processing (Fan et al., 2021).

The current study is not without weaknesses. The data presented here are from the baseline assessments from an ongoing longitudinal study. Thus, the main outcome of comparing groups cannot establish causality. While associations between historic pesticide exposure and brain network topology suggest that pesticide exposure may be a contributing factor, further work is clearly needed. Future analyses of the longitudinal data currently being collected for this study may help clarify the role of pesticide exposure in the observed DMN connectivity differences. Life history calendars have an important and valued role in assessing environmental exposures that occurred prior to data collection (Hoppin et al., 2001; Rodriguez et al., 2012; Zahm et al., 2001). Our use of these historic exposure data allowed for the examination of relationships between probably exposure and brain network organization. However, there is a risk of fabrication or faulty memory in collecting historic early life data from the children’s parents. Also, in this study, we examined whether any history of pesticide exposure regardless of the developmental period (as averaged across prenatal, early childhood and late childhood) was associated with the altered network topology. Although we present supporting supplemental results for different developmental time periods, future studies are needed for a more thorough assessment of developmental age at time of pesticide exposure. In addition, even though we put considerable effort into collecting data from populations matched with respect to their socioeconomic status and critical demographics, there are always potential effects from such variables even after controlling for them as confounding effects. Finally, the study sample was restricted to North Carolina and exposure profiles from the children included may not be representative of children living in different geographic areas.

In summary, this study provides compelling evidence for the likely adverse effects of pesticide exposure on one of the most critical subnetworks of the brain, the DMN, in children from rural farmworker families, and provides evidence implicating historic childhood pesticide exposure as a potential cause. Although the ultimate implications of the observed topological differences require further studies, these findings suggest that pesticide exposure in Latinx children from rural farmworker families could place this vulnerable population at a higher risk for future neurobiological, neurocognitive, or neurobehavioral problems. Continued efforts to reduce pesticide exposure in all children must be continued given the growing evidence of adverse health effects.