Causal effects of cingulate morphology on executive functions in healthy young adults.

In this study, we want to explore evidence for the causal relationship between the anatomical descriptors of the cingulate cortex (surface area, mean curvature-corrected thickness, and volume) and the performance of cognitive tasks such as Card Sort, Flanker, List Sort used as instruments to measure the executive functions of flexibility, inhibitory control, and working memory. We have performed this analysis in a cross-sectional sample of 899 healthy young subjects of the Human Connectome Project. To the best of our knowledge, this is the first study using causal inference to explain the relationship between cingulate morphology and the performance of executive tasks in healthy subjects. We have tested the causal model under a counterfactual framework using stabilized inverse probability of treatment weighting and marginal structural models. The results showed that the posterior cingulate surface area has a positive causal effect on inhibition (Flanker task) and cognitive flexibility (Card Sort). A unit increase (+1 mm2 ) in the posterior cingulate surface area will cause a 0.008% and 0.009% increase from the National Institute of Health (NIH) normative mean in Flankers (p-value <0.001), and Card Sort (p-value 0.005), respectively. Furthermore, a unit increase (+1 mm2 ) in the anterior cingulate surface area will cause a 0.004% (p-value <0.001) and 0.005% (p-value 0.001) increase from the NIH normative mean in Flankers and Card Sort. In contrast, the curvature-corrected-mean thickness only showed an association for anterior cingulate with List Sort (p = 0.034) but no causal effect.

INTRODUCTION

The neural basis of executive functions (EFs) is a key topic in neuroscience research. Up to this date, we do not have a mechanistic explanation of the anatomic structures that cause individual differences. One of the most relevant structures associated with EF is the cingulate cortex (CC), as evidenced by studies based on functional magnetic resonance imaging (MRI) (Amanzio et al., 2020; Palermo et al., 2018) for group average and clinical populations (Braden et al., 2017; Zhan et al., 2018).

Instead of focusing on functional brain research (fMRI), here we will focus on structural correlates. The literature shows that structure follows functions (Wenger et al., 2017), and the structural differences may underpin the intra‐ and inter‐individual differences in cognitive functions. An extensive recent study found some of the strongest associations between overall intelligence and the volume of the insula and posterior cingulate (Deary et al., 2021). However, they did not use the components of the volume: thickness and cortical surface area (Deary et al., 2021). Other studies relate the EFs with structural correlates of the CC in the clinical (Kempton, 2011; McClintock et al., 2010; Saleh et al., 2017; Shinde et al., 2021) and healthy population (Amiez et al., 2018; Bento‐Torres et al., 2019).

These studies are associational and based on the clinical populations and very few controls or small sample sizes. Association is much less than causal analysis, specifically under the modern framework of potential outcomes and counterfactuals. We have implemented the causal modeling under the counterfactual framework using stabilized inverse probability of treatment weighting (IPTW) and marginal structural models (MSMs). We want to assess how the structural properties of the CC causally affect the performance of EFs. Individual differences in the performance of cognitive tasks related to the EFs like flexibility, inhibitory control, and working memory can be measured in well‐known tasks such as the Dimensional Change Card Sort (DCCS) (Zelazo, 2006), the Ericksen's Flanker task (Eriksen & Eriksen, 1974), and List Sorting (Tulsky et al., 2013). The current study is based on a cross‐sectional sample of healthy young subjects of the Human Connectome Project (HCP).

To the best of our knowledge, no studies demonstrate the causal relationship between cingulate structural anatomy and individual performance of EFs under the modern statistical framework of counterfactual analysis.

MATERIAL AND METHODS

This section describes the data and the details of statistical analysis. All the analysis was done in R (R Core Team, 2020). Figure 1 shows an overview of the methodology.

FIGURE 1

Overview of methodology

Participants

Our data sample consisted initially of 1206 subjects from Human Connectome Data (HCP) (WU – Minn Consortium Human Connectome Project, 2017), from which we excluded 307 subjects (one subject from each pair of mono and dizygotic twins) to avoid inter‐sibling correlations. We have included 899 subjects in this study, aged 22–36 years. Of all subjects, 47.6% were male; see Table 1 for descriptive data statistics. All subjects are healthy and do not have significant neurological or psychiatric disorders. More details can be found at http://humanconnectome.org/storage/app/media/documentation/s1200/HCP_S1200_Release_Reference_Manual.pdf).

TABLE 1

Summary of data

Variables	Mean (SD)	Variables	Mean (SD)
Number of subjects	899	List Sort	111.02 (11.30)
Gender = M (%)	428 (47.6)	ACC area	5700.01 (781.91)
Age	28.68 (3.76)	PCC area	3218.10 (448.81)
Handedness	66.34 (43.79)	ACC volume	18,434.99 (2467.55)
Gray matter volume	686,095.96 (67,213.24)	PCC volume	10,239.83 (1469.63)
Card Sort	114.66 (10.56)	ACC thickness	2.99 (0.11)
Flanker	111.28 (10.09)	PCC thickness	2.68 (0.10)

Abbreviations: ACC, anterior cingulate cortex; PCC, posterior cingulate cortex.

Variables

Demographics

We include age, gender, handedness, race, and father, mother ID of each participant in the analysis.

Age in years.

Gender: female or male.

Father/mother ID: a unique numeric ID to identify siblings and family.

Handedness has a value [−100, 100], where negative values were for left‐handed people and positive values for right‐handed people. We dichotomized handedness for ease of analysis.

Race was a six‐factor variable with possible values of

White

Black/African American

Asian/Hawaiian/Pacific Islanders

American Indian/Alaskan

More than one

Unknown

Cognitive tests for executive functions

The National Institute of Health (NIH) designed the cognitive tasks employed here. The brain's EF depends, among other skills, on cognitive flexibility, inhibition, and working memory (Diamond, 2013). We selected three tasks to assess EF using the neuropsychological battery in the HCP for our study. The tasks were: DCCS (cognitive flexibility), Flanker (inhibition), and List Sorting (working memory). For each task, we have used unadjusted scores. A score of 100 indicates the mean score for NIH normative sample. A higher value indicates better performance (Cole et al., 2021). Further explanation about tasks can be found at https://wiki.humanconnectome.org/display/PublicData/HCP+Data+Dictionary+Public‐+Updated+for+the+1200+Subject+Release#HCPDataDictionaryPublicUpdatedforthe1200SubjectRelease‐Instrument:ExecutiveFunction/Inhibition(FlankerTask).

Cognitive flexibility (Dimensional Change Card Sort)

DCCS measures cognitive flexibility where the subject matches the initial pictures with other subsequent picture cards in terms of either shape or color. Participants switch between the shape and colors as per instruction in different trials. These changing instructions (sorting based on color or shape) measure cognitive flexibility. The test scores are acquired based on the subject's reaction time and accuracy of sorting (Zelazo, 2006). NIH toolbox DCCS is a <5‐min test and valid for individuals 3–85 years old (Cole et al., 2021).

Inhibition (Flanker task)

The Flanker task measures both participant's attention and inhibitory control. The test requires the participant to focus on a given stimulus while inhibiting attention to stimuli (fish for ages 3–7 or arrows for ages 8–85) flanking it. Sometimes the middle stimulus is pointing in the same direction as the “Flanker” (congruent) and sometimes in the opposite direction (incongruent). Scoring depends on accuracy and reaction time, and the test takes approximately 3 min to administer. This test is suitable for ages 3–85.

Working memory (List Sorting)

This task assesses the working memory using List Sorting. Participants are required to sort pictures (from smallest to largest in terms of size) of different food items and/or animals. Images provide a written name as well as a sound clip. The task has two conditions: (1) list requires sorting only one series of items (either food or animals), whereas the (2) list condition requires sorting both, one after other (food items followed by animals).

Neuroimaging variables

The neuroimaging variables of brain gray matter employed in our study are the FreeSurfer (FS) (https://www.sciencedirect.com/topics/medicine-and-dentistry/freesurfer) output. FS is a suite of tools for analyzing neuroimaging data that provides various algorithms to quantify the human brain's functional, connectional, and structural properties (http://surfer.nmr.mgh.harvard.edu) (Fischl et al., 1999; Fischl & Dale, 2000). We projected the Destrieux cortical atlas (Destrieux et al., 2010) from FS for brain parcellation. Destrieux atlas is gyri and sulci‐based and has 148 regions.

Total gray matter volume

Total gray matter volume is the sum of the volume at each vertex for the whole brain, and the measurement unit is mm3. Volume is the product of cortical surface area and gray matter thickness at each vertex.

Anatomical variables for cingulate cortex

Destrieux cortical atlas parcellates the CC into five subregions: two from the anterior cingulate cortex (ACC) and three posterior regions (Destrieux et al., 2010; Mashhoon et al., 2014). For this study, we have combined the five CC subregions into two main areas, anterior and posterior CC.For illustration purposes only, Figure 2 shows a parcellation of the CC. Morphometry statistics generated by FS are the physical brain feature measurements like total surface area, mean cortical thickness, and total volume. The surface area is the sum of all vertex in a brain parcel, and the measuring unit is mm2. Thickness is measured as the mean/average distance to the gray‐CSF boundary and measured in mm. Furthermore, we have used curvature‐corrected thickness for the HCP data set, which is available at https://balsa.wustl.edu/sceneFile/L6vp9. Average surface curvature was regressed out, resulting in curvature‐corrected thickness, as explained in Glasser and Van Essen (2011) and Sigalovsky et al. (2006). Volume is the product of area and thickness at each vertex and measured in mm3. Most studies on brain structure have focused on brain volume. However, volume is the composition of area and thickness. Studies that analyze all three measures simultaneously are scarce. It is also not well established which component of volume (area or thickness or both) is more significant (Cox et al., 2018).

FIGURE 2

Parcellation of cingulate cortex (CC) in five subregions from Destrieux cortical atlas (Destrieux et al., 2010). Parcellated anterior CC is divided into anterior cingulate gyrus and sulcus (ACC) and middle anterior cingulate gyrus and sulcus (aMCC). Parcellated posterior CC is divided into middle posterior cingulate gyrus (pMCC), dorsal posterior cingulate gyrus (dPCC), and ventral posterior cingulate gyrus (vPCC).

Anterior CC

Anterior cingulate gyrus and sulcus (ACC)

Middle anterior cingulate gyrus and sulcus (aMCC)

Posterior CC

Middle posterior cingulate gyrus (pMCC)

Dorsal posterior cingulate gyrus (dPCC)

Ventral posterior cingulate gyrus (vPCC)

Traditional statistical analysis

Handling incomplete or missing data

As the first step, we check for any missing values and omit cases with incomplete data.

Outlier detection and elimination

In the second step, we screen for the multivariate outliers. For multivariate outlier detection, we utilized all the measured variables. A multivariate outlier is any subject for which the combination of values for different variables is different from most subjects. It is essential to detect and handle such cases to avoid bias in the results. We applied Mahalanobis distance, Leverage, and Cook's distance to detect the outliers. Mahalanobis distance assigns a score to each subject, which is the distance between the subject and the centroid of the distribution. The distance cut‐off value was selected using a p‐value <0.001 (Tabachnick et al., 2007). We implement Leverage to measure how much one subject influences the slope of regression line β (Welsch, 1982).

Furthermore, we administered Cook's distance which accounts for individual Leverage and residual (Cook & Weisberg, 1982). We carried out this analysis using “Mahalanobis,” “hatvalues,” and “cook.distance” functions in R (R Core Team, 2020). Any subject beyond the cut‐off limit in at least two measures is considered an outlier and excluded from the analysis.

Checking data for normality, homoscedasticity, and linearity

Checking for normality, linearity, and homoscedasticity is essential for any inference‐based analysis because results can degrade if data is skewed or curvilinear. Therefore, we examined multivariate normality using the standardized regression residuals for Flanker's score on all anatomical variables. A histogram and a quantile–quantile (QQ) plot are employed to check for the multivariate normality of data regarding skewness and kurtosis (Tabachnick et al., 2007). Homoscedasticity defines that residuals have the same variance for all the values of independent variables. Homoscedasticity and normality ensure the linear trends in data. We have used a scatter plot for standardized residuals versus fitted values to check homoscedasticity (Tabachnick et al., 2007).

Regression analysis between anatomical variables and EF‐based tasks

To simplify and establish a causal diagram (directed acyclic graph [DAG]), we implemented a linear regression to observe simple associations between the EFs and anatomical variables. It is evident from the literature that the association of intelligence with cortical thickness and area varies across different age groups and racial backgrounds (Frangou et al., 2022; Lett et al., 2020; Menary et al., 2013). Thus, we have included age and race as covariates for the initial estimate. We choose the significant variables for further analysis at a p‐value <0.05.

Examining the multicollinearity of anatomical markers

We tested anatomical variables for multicollinearity to identify the relationships between the three morphological metrics (volume, area, and thickness). We implemented bivariate correlations between the three anatomical markers. In the case of high correlation, we only selected the more simple markers to simplify the causal DAG.

Causal inference under the counterfactual framework

Causal directed acyclic graph

After pre‐processing, cleaning, and more traditional statistical analysis, we established a simplified causal DAG to identify the causal relationship between the EFs and cingulate anatomical markers.

Marginal structural models

The MSM is a subclass of statistical models. MSMs implement the marginal distribution of counterfactual random variables or potential outcome models for causal inference (Breskin et al., 2018). Ideally, for a dichotomous treatment, the true causal effect (CE) is the difference between the two outcomes for the same individual (1) when the treatment was applied and (2) when the treatment was not applied . Treatment (area, volume, or thickness) is continuous in this study. However, we will use a dichotomized version of the treatment variable to establish an initial model and theoretical basis. If is the mean of the continuous treatment for HCP young adults' data set, then for subject is Furthermore , there is a CE if the outcome for the individual is different under treatment and no treatment However, only one of these outcomes is observable, while the other is a counterfactual or potential outcome. For example, the subject can only have the ACC surface area greater than or less than M; however, they cannot simultaneously have both. As the individual CE is not identifiable, the average causal effect (ACE) or the average treatment effect (ATE) is calculated in a population of individuals. Furthermore, for a sample (which is a subset of the population), ATE is estimated using the conditional averages is the observed average of the sample, is the outcome for treated individuals, and is the outcome for untreated individuals. is the unbiased estimator for true ATE under identifiability conditions (exchangeability/ignorability, positivity, consistency); more details can be found in James and Hernán (2020). However, these conditions can be violated by confounding factors. Therefore, confounding factors should be appropriately handled, and the approaches to deal with measured confoundings are discussed in the next section.

Inverse probability of treatment weighting: A marginal method for confounding control

Measured confounding in the observational studies is commonly dealt with using conditional approaches (stratification) and marginal approaches (G‐methods). Stratification‐based methods estimate the association between treatment and outcome in the sample under study (conditional). In contrast, G‐methods estimate the CE on the entire population (marginal) (Ji et al., 2020).

IPTW comes under G‐methods (James & Hernán, 2020). We implemented IPTW because we are interested in the CE on the population. IPTW creates a pseudo‐population with no confounding. where T is 1 for treated or 0 for untreated individuals, and is the probability of being treated given the confounding levels l (Pezzi et al., 2016). We balance the confounders across treatment levels by assigning weights to each sample.

Problem formulation for marginal structural model

MSM estimation via IPTW is a two‐step process. Firstly, IPT weights are computed. We use these weighted samples to fit an outcome model in the second step. Moreover, the unmeasured confounding is handled using a sensitivity analysis, which estimates the amount of unmeasured confounding required to nullify the estimated ATE (VanderWeele & Ding, 2017).

We implemented the marginal models, using cingulate anatomy as a continuous treatment variable and EFs as an outcome , age, and race as covariates . The measured confounding is age, race, gender, handedness, father ID, mother ID, and total gray matter volume. Father and mother ID account for any unmeasured common genetic and environmental factors within families. L also includes age and race to ensure exchangeability within levels of covariates (James & Hernán, 2020). We applied stabilized or augmented IPTW to avoid biased results due to large weights for the measured confounding variables L. The outcome model for ATE on Y is We evaluated the outcome model for the MSM using the GeePack R package (Halekoh et al., 2006). The R code can be found at https://github.com/CCC-members/HCP-Anat-R.

Sensitivity analysis

Sensitivity analysis assesses the robustness of CEs against unmeasured confounding. We implement “E‐value” as a sensitivity parameter administered for the stability of causal estimates in observational studies. E‐value defines the association between the unmeasured confounder and treatment/outcome in the risk ratio scale, which is needed to completely explain the estimated ATE. Thus, the larger the E‐value, the more stable the estimated ATE (VanderWeele & Ding, 2017). However, E‐value should not be used alone and should be administered cautiously (Ioannidis et al., 2019; VanderWeele et al., 2019). We will only consider the E‐value for this study if the p‐value is significant at alpha 0.05. We applied the E‐value function for continuous treatment and outcome (Mathur et al., 2018) based on estimates from MSM. E‐value takes a parameter delta to dichotomize the exposure. Our choice for delta is the mean value for each anatomical variable. For dichotomization of the continuous outcome, the E‐value function utilizes the effect size (Chinn, 2000; VanderWeele, 2017).

RESULTS

Data pre‐processing and cleaning

Our data consist of 899 samples, out of which 231 subjects had missing values. We excluded those subjects from the analysis, and 668 subjects remained in the study. We implemented outlier detection based on Mahalanobis, Leverage, and Cook's distance; and excluded 17 more subjects as outliers from the analysis based on the distance cut‐off value (if the subject was an outlier in two out of three distances). The outlier elimination resulted in n = 668 valid samples with which we carried out the rest of the analysis. We checked for multivariate normality, linearity, and data homogeneity using standardized results from regression . Figure 3 summarizes the results. The QQ plot (Figure 3a) for standardized residuals from the regression shows a linear trend. The histogram for standardized residuals (Figure 3b) shows a normal distribution, so no corrective measures were needed.

FIGURE 3

Results for data's multivariate normality, linearity, and homogeneity. The plots show the result of regression “Flanker~anatomical variables.” (a) A quantile–quantile plot for the standardized residual versus a theoretical normal distribution shows a linear trend. (b) Histogram for standardized regression residuals to check for multivariate normality. (c) A scatter plot between standardized residuals and z‐scored predicted values. The resultant spread is homogeneous.

Furthermore, to check for homogeneity/homoscedasticity, we have plotted z‐scored fitted values versus standardized residuals in a scatter plot (Figure 3c). The results are homogenous as the spread across zero for both axes [−2, 2] is the same. The solid red line in the plot marks zero for both axes.

Linear trend

Before estimating the ATE using MSM, we used linear regression with covariates (age, race) to check for simple associations between the EFs and measured anatomical variables. Results are in scatter plots (Figures 4, 5, 6) for Flanker, Card Sort, and List Sort, respectively. Each figure shows six subplots for anatomical measurements and test scores; the first column is for surface area, the second for curvature‐corrected mean thickness, and the third for volume. Each column has two subplots, one for each subregion in the CC; anterior and posterior. The x‐axis shows the z‐scored area/thickness/volume, and the y‐axis shows test scores (Flanker in Figure 4, Card Sort in Figure 5, List Sort in Figure 6). The solid red line shows the linear trend between x and y, and the plots also exhibit a p‐value for the regression .

FIGURE 4

Six scatter plots for measured anatomical variables versus Flankers. The first row from the top is for anterior cingulate cortex (CC; area, followed by curvature‐corrected thickness and volume). The second row is for posterior CC. The y‐axis is test scores for the Flanker and the x‐axis for standardized anatomical measurements. The solid red line shows the linear trend between x and y with respective p‐values for “Flanker~Anatomy + Age + Race” regressions. Area and volume are significant at alpha 0.05.

FIGURE 5

Six subplots each for one measured anatomical variable. The first row is for anterior cingulate, and the second row is for posterior cingulate. Each row includes three subplots for anatomical measures. The y‐axis is test scores for Card Sort and the x‐axis for anatomical measurements. The solid red line shows the linear trend for y~x.

FIGURE 6

Linear trend between List Sort and anatomical variables. Columns represent the area, thickness, and volume, respectively, whereas rows are for subregions anterior cingulate cortex (CC) and posterior CC.

Results showed that the variables for area/volume have a significant association with Flanker and Card Sort with p‐values ≤0.05. In contrast, variables for curvature‐corrected‐mean thickness only showed an association for anterior cingulate with List Sort at a significance level of 0.05. We do not adjust these results for multiple comparisons.

Multicollinearity in morphological metrics

We have also tested anatomical variables for multicollinearity. To check for the relationship between the three morphological metrics (volume, area, and thickness), we have used bivariate correlations and found that the surface area and volume are highly correlated. Figure 7 shows a correlation plot between anatomical variables; the dark blue indicates highly correlated variables with the highest correlation value, 0.96. Volume is a composite variable based on area and thickness. Based on high correlation values of volume with surface area, we have decided to exclude volume and use only the surface area and cortical thickness for causal analysis.

FIGURE 7

Correlation plot for anatomical measurements. The blue color shows a positive correlation with the darkest value of 0.96 between area volume variables for any specific region.

Inverse probability treatment weighting

The stabilized IPTW resulted in a nominal range with no extreme weights. The weight range was [0, 26] for surface area and [0, 4] for thickness. Figure 8 shows the density plot for weights. The x‐axis is the weights in the log scale, and the y‐axis is the density for each weighting value, whereas different colors represent the CC‐specific weights.

FIGURE 8

Inverse probability of treatment weighting (IPTW) density plot, the x‐axis is log weights for the pseudo‐population adjusted for age, gender, race, handedness, father and mother ID, and total gray matter volume.

Marginal structural modeling

After analyzing the linear trends and collinearity between the anatomical variables, we modified our initial MSM model based on the findings; volume was excluded from the causal analysis. Figure 9 illustrates the updated model for MSM.

FIGURE 9

Directed acyclic graph (DAG) for inverse probability of treatment weighting (IPWT) and marginal structural model (MSM). DAG shows the average treatment effect (ATE) direction from treatment variable T as cingulate surface area/to outcome Y (cognitive test scores); age, race, gender and handedness, and total gray matter volume are confounding factors.

We tested four marginal models for each cognitive test score, two for the area and two for the thickness. We tested 12 marginal models and summarized the results in Table 2. Four models were significant at a significance level of 0.05. None of the List Sort tests were significant. ATE for Flanker and posterior surface area was highly significant with a p‐value <0.001 and a causal estimate of 0.008 (95% CI [0.003–0.012]). PCC area and Card Sort also resulted in significant CE with ATE 0.009 (95% CI [0.003–0.015]) and p‐value = 0.005. Anterior surface area also showed a significant CE for both Flanker and Card Sort with a p‐value of <0.001 and 0.001, respectively. ATE of ACC area on Flanker was 0.004 (95% CI [0.002–0.006]), and on Card Sort was 0.005 (95% CI [0.002–0.008]).

TABLE 2

Results for marginal structural model (MSM) between executive functions test scores (Card Sort, Flanker, List Sort) and CC surface area (anterior and posterior)

Y	X (CC)	p‐value	ATE	95% CI‐lower	95% CI‐upper	E‐value
Card Sort	ACC area	0.001**	0.005	0.002	0.008	21.455
	PCC area	0.005**	0.009	0.003	0.015	20.901
	ACC thickness	0.280	−3.874	−10.912	3.164	4.87
	PCC thickness	0.503	−3.003	−11.804	5.796	3.421
Flanker	ACC area	<0.001***	0.004	0.002	0.006	12.915
	PCC area	<0.001***	0.008	0.003	0.012	16.639
	ACC thickness	0.289	3.966	−3.372	11.304	5.216
	PCC thickness	0.483	2.904	−5.227	11.036	3.427
List Sort	ACC area	0.958	0.000	−0.006	0.006	1.347
	PCC area	0.805	0.001	−0.009	0.011	2.107
	ACC thickness	0.340	−4.301	−13.139	4.537	4.999
	PCC thickness	0.286	5.827	−4.884	16.539	6.385

Note: First columns show p‐values, followed by average treatment effect (ATE), 95% upper and lower bounds, and E‐value for sensitivity analysis. ***p = 0; **p = 0.001; *p = 0.01

Abbreviations: ACC, anterior cingulate cortex; ATE, average treatment effect; CC, cingulate cortex; CI, confidence interval; PCC, posterior cingulate cortex.

Sensitivity analysis

To check for the stability of MSM results, we have used E‐values as sensitivity analysis. The significant models from MSM at alpha 0.05 resulted in the largest E‐values. The PCC area resulted in 16.639 for the Flanker and 20.901 for Card Sort (Table 2). E‐values for ACC models were 21.455 for Card Sort and 12.915 for Flanker (Table 2).

DISCUSSION

We have analyzed healthy young adults' data to estimate the CEs of cingulate anatomy (total surface area, curvature‐corrected mean thickness, and total volume) on the brain's EFs for cognitive flexibility, control inhibition, and working memory. We found that the CC surface area and volume were associated with EFs (Card Sort and Flankers task) with a p‐value <0.01 when age and race were included as covariates. On the other hand, curvature‐corrected mean thickness for anterior cingulate showed association with List Sort with a p‐value of 0.0348. These results are supported by previous studies, which state that thickness and cortical surface area are phenotypically different and evolve independently (Chen et al., 2015; Geschwind & Rakic, 2013; Winkler et al., 2010) and have different associations with cognition and brain disease (Schnack et al., 2015; Vuoksimaa et al., 2016). In a study of human intelligence and its relationship with the cortical thickness (structural) and neural activations (functional), (Choi et al., 2008) found that ACC is functionally associated with intelligence. Moreover, a recent study by Lett et al. (2020) found that general intelligence scores were associated with ACC surface area and thickness (Lett et al., 2020).

Furthermore, as the volume is a composite measure based on area and thickness, we found volume is highly associated with the surface area with a positive correlation value of 0.96. That suggests that the cortical surface area is CC's most relevant morphological metric for the EFs. The importance of the cortical surface area has also been stressed by Lett et al. (2020) as a mediator of the relationship between polygenic scores for intelligence and general intelligence.

Based on the findings mentioned above, we analyzed only the surface area and average cortical thickness for the causal inference. We have implemented IPTW and MSMs to estimate the CEs under the potential outcome/counterfactual paradigm. The results showed that the surface area of the posterior cingulate has a positive ACE on test scores for Flanker (control inhibition) and Card Sort (cognitive flexibility). A unit increase (+1 mm2) in posterior cingulate surface area will cause a 0.008% (Flanker), 0.009% (Card Sort) increase in the test scores from the NIH normative mean score of 100, after accounting for age, gender, handedness, total gray matter volume for the brain, race, and family (indirectly accounting for common genetic and environmental factors).

E‐values test for the robustness of the results under any unmeasured confounding factors, and the results showed that posterior CC has robust ATE. E‐value used dichotomization of the independent/treatment variable, and results are in risk ratio scale. We used the average surface area for dichotomization. The results showed an E‐value of 16.639 for Flanker and 20.901 for Card Sort. This suggests that for any unmeasured confounding factor to change the ATE of PCC, it would need to have a 16‐ to 20‐fold stronger association with both exposure and outcome. Some recent associational studies between posterior cingulate and general intelligence scores support our findings (Basten et al., 2015; Nomi et al., 2018). An extensive recent study from UK BIOBANK data found some of the strongest associations between a latent variable of intelligence “g” (integrated by Matrix Reasoning, Trail Making, Symbol‐Digit, and Verbal Numerical Reasoning) and the volumes of the insula, posterior cingulate/precuneus r < 0.20 (Deary et al., 2021). Zhao et al. (2019) found in another large‐scale study that posterior cingulate volume and generic intelligence have shared genetic variants.

Moreover, the anterior surface area also showed a significant CE for Flanker and Card Sort, with a p‐value of ≤0.001. A unit increase (+1 mm2) in the anterior cingulate surface area will cause a 0.004%, 0.005% increase in Flankers and Card Sort. E‐values for ACC models were 21.455 for Card Sort and 12.915 for Flanker. Lett and colleagues found similar results where general intelligence scores were associated with ACC surface area (Lett et al., 2020).

In addition, even though we found an association between anterior cingulate thickness and List Sort (working memory), no CE was found, suggesting that their relationship is not causal. The association might be due to some common underlying cause (Altman & Krzywinski, 2015). One of the explanations could be the age range of the samples, as the HCP data set is restricted to young adults. Several studies where thickness measures were linked to intelligence (Karama et al., 2011; Luders et al., 2011; Schnack et al., 2015) and other cognitive performances are focused on children and adolescents populations, where neurodevelopmental factors are present in the critical age range for brain maturation and plasticity.

Furthermore, no studies compare causal relationships among all three cortical measures in the same sample, and also, most of them rely upon global measures of intelligence. Our findings provide a new perspective for examining EFs and their relationship with cingulate anatomy.

LIMITATIONS

Even though the sensitivity analysis protects against measured confounding, a more elaborated causal DAG can be tested to understand the chain of causality better. For example, the individual differences in the cognitive performance of tasks implying EFs such as working memory, flexibility, and inhibitory control, are not solely related to the brain gray matter measures studied here. We have implemented a simple sensitivity analysis; a more robust alternative can be administered. Moreover, a more heterogeneous sample in terms of age and race for healthy subjects can allow us to understand better the contribution of covariates in the relationship between anatomy and EFs. Here, we are focused only on the CC, but we will address all the circuits involved in EFs, especially the pre‐frontal cortex, in future studies. Moreover, we restricted our analysis to mean corrected thickness. A complete vertex‐wide study using whole‐brain structural information, connectivity, and genetic information can also be done in the future.

CONCLUSION

Many associational studies exist that describe the role of the CC in EFs. However, to the best of our knowledge, this is the first study using counterfactual‐based causal analysis to map the relationship between morphology (surface area, curvature‐corrected mean thickness, and volume) of the CC and EF‐based tasks (Flanker, List Sort, Card Sort). We identified that volume and surface area are highly associated with the test scores. Cortical volume is the product of surface area and cortical thickness and is positively correlated with surface area (Pearson correlation >0.90). We also found no causal relationship even though mean ACC thickness was associated with EF‐based tasks. These findings suggest that the CC surface area is the primary morphological metric that has any causal relationship with EF‐based tasks. The causal analysis between CC area and tasks resulted in a positive CE of anterior and posterior cingulate on inhibition (Flanker task) and cognitive flexibility (Card Sort). A unit increase (+1 mm2) in the posterior cingulate surface area will cause a 0.008% and 0.009% increase from the NIH normative mean in Flankers (p‐value <0.001) and Card Sort (p‐value 0.005), respectively. Moreover, a unit increase (+1 mm2) in the anterior cingulate surface area will cause a 0.004% and 0.005% increase in Flankers and Card Sort.