Assessing brain structural associations with working-memory related brain patterns in schizophrenia and healthy controls using linked independent component analysis.

Schizophrenia (SZ) is a psychotic disorder with significant cognitive dysfunction. Abnormal brain activation during cognitive processing has been reported, both in task-positive and task-negative networks. Further, structural cortical and subcortical brain abnormalities have been documented, but little is known about how task-related brain activation is associated with brain anatomy in SZ compared to healthy controls (HC). Utilizing linked independent component analysis (LICA), a data-driven multimodal analysis approach, we investigated structure-function associations in a large sample of SZ (n = 96) and HC (n = 142). We tested for associations between task-positive (fronto-parietal) and task-negative (default-mode) brain networks derived from fMRI activation during an n-back working memory task, and brain structural measures of surface area, cortical thickness, and gray matter volume, and to what extent these associations differed in SZ compared to HC. A significant association (p < .05, corrected for multiple comparisons) was found between a component reflecting the task-positive fronto-parietal network and another component reflecting cortical thickness in fronto-temporal brain regions in SZ, indicating increased activation with increased thickness. Other structure-function associations across, between and within groups were generally moderate and significant at a nominal p-level only, with more numerous and stronger associations in SZ compared to HC. These results indicate a complex pattern of moderate associations between brain activation during cognitive processing and brain morphometry, and extend previous findings of fronto-temporal brain abnormalities in SZ by suggesting a coupling between cortical thickness of these brain regions and working memory-related brain activation.

Introduction

Schizophrenia (SZ) is a debilitating illness characterized by delusions, hallucinations and disorganized thought. Impairments of cognitive functions, such as working memory, are also considered core features of the disorder (Green, 2006; Kahn and Keefe, 2013; Park and Gooding, 2014) and have been linked to genetic liability (Agnew-Blais and Seidman, 2013; Reichenberg and Harvey, 2007). However, little is known about how cognitive dysfunction is related to underlying brain anatomy and brain function, and improving our understanding of these associations may help uncover the neuronal substrates of the disease (Kahn and Keefe, 2013; Schultz et al., 2012a).

Neuroimaging studies have identified brain regions and networks involved in cognitive processing (Alnaes et al., 2015; Cabeza and Nyberg, 2000; Cole et al., 2014), and both hypo- and hyperactivation (Brandt et al., 2014; Glahn et al., 2005; Kraguljac et al., 2013; Ragland et al., 2007), as well as brain network dysconnectivity (Brandt et al., 2015; Fornito et al., 2012b; Kaufmann et al., 2015; Pettersson-Yeo et al., 2011) have been reported in SZ. In particular, fMRI studies have demonstrated abnormal neuronal recruitment during working memory processing (Anticevic et al., 2013; Callicott et al., 2003; Glahn et al., 2005; Henseler et al., 2009; Karlsgodt et al., 2007; Kim et al., 2010; Potkin et al., 2009; Schneider et al., 2007; Thormodsen et al., 2011) reflecting dysfunction of a fronto-parietal network including the lateral prefrontal cortex and posterior parietal cortex, as well as anterior cingulate and other regions related to execution of challenging cognitive tasks (Owen et al., 2005). This network shows increased activation during cognitive processing and is thought to reflect focus on task-relevant information, like updating information in working memory. Similar networks have been referred to as central executive (Bressler and Menon, 2010; Sridharan et al., 2008), executive control (Seeley et al., 2007), and task-positive (Fox et al., 2005) networks. The function of such a task-positive fronto-parietal network and its interactions with a task-negative default-mode network (DMN) are considered to play important roles in cognition (Cocchi et al., 2013; Fornito et al., 2012a). The DMN shows task-related deactivation and has been hypothesized to reflect the suppression of task-irrelevant internal activity to optimize goal-directed cognition (Anticevic et al., 2012). Cognitive impairments in SZ may not only be associated with the failure to recruit the fronto-parietal network, but also a relative lack of DMN deactivation (Whitfield-Gabrieli and Ford, 2012). In line with this, a failure of DMN deactivation may reflect a key feature of SZ (Broyd et al., 2009; Landin-Romero et al., 2015; Pomarol-Clotet et al., 2008), and the reciprocity between these two networks may be affected in patients (Anticevic et al., 2013; Nygard et al., 2012; Whitfield-Gabrieli et al., 2009).

Structural brain abnormalities have been consistently reported in patients with SZ, including reductions in cortical thickness, surface area, gray matter volumes, and gyrification (Ellison-Wright and Bullmore, 2010; Gupta et al., 2015; Nesvag et al., 2014; Rimol et al., 2010; Rimol et al., 2012), global and subcortical volumes (van Erp et al., 2015), as well as white matter microstructure as measured by diffusion tensor imaging (DTI) (Ellison-Wright and Bullmore, 2009). Gray matter reductions have been identified in widespread brain regions, but are particularly pronounced in fronto-temporal regions (Glahn et al., 2008; Gupta et al., 2015; Nesvag et al., 2008; Rimol et al., 2010; Schultz et al., 2010), including the insula (Glahn et al., 2008; Shepherd et al., 2012). It has been suggested that these structural brain alterations are associated with abnormal integration of information between frontal and temporal cortical areas (Friston and Frith, 1995; Schultz et al., 2012b). Similarly, functional neuroimaging studies have reported disrupted fronto-temporal connectivity in SZ related to cognitive processing including working memory (Cocchi et al., 2014; Crossley et al., 2009; Meyer-Lindenberg et al., 2001; Wolf et al., 2007), indicating that cognitive impairment is a result of underlying brain dysconnectivity.

In order to understand more of the relationship between structural and functional brain abnormalities in SZ, several studies have examined associations between structural MRI and neuropsychological performance (Ehrlich et al., 2012; Gutierrez-Galve et al., 2010; Hartberg et al., 2010), or between brain structure and brain activation during cognitive tasks using different methods and modalities (see Schultz et al., 2012a for review). A common finding has been significant structure–function associations in patients, but not in healthy controls (Fusar-Poli et al., 2011b; Pujol et al., 2013; Schultz et al., 2012b), indicating differential patterns of associations. Previous studies have however often performed separate analyses of brain structure and function without correlating them directly (Pomarol-Clotet et al., 2010; Skudlarski et al., 2010), while only a few have combined structural and functional measures in the same analysis allowing for interpretation of joint features across modalities (Calhoun et al., 2006; Correa et al., 2008; Michael et al., 2011).

Further, a common approach has been to use structural and functional regions of interest (Harms et al., 2013; Pujol et al., 2013; Schultz et al., 2012b), thus restricting the anatomical interpretations to these regions. The sample size has often been small, comprising around 15 patients and 15 controls (Calhoun et al., 2006; Fusar-Poli et al., 2011a; Pujol et al., 2013; Rasser et al., 2005), which have limited the generalizability. Thus, there is a need for studies combining structural and functional modalities in the same analysis in large samples of patients and controls. Lastly, despite the relevance of task-positive and task-negative networks in cognition and SZ, only one regions-of-interest based study comprising a small number of patients has investigated the relationship between these functional brain systems and brain anatomy in SZ (Pujol et al., 2013).

Summarized, functional and structural imaging has documented a plethora of brain abnormalities in SZ. However, more knowledge is needed to map task-related brain activation onto underlying brain structure in SZ and HC, and to assess if patients show differential structure–function relationships compared to HC. These are important questions pertaining to the fundamental aim of delineating associations between structural and functional properties of the brain, and for increasing the understanding of the mechanisms of severe mental illness, such as SZ.

Thus, the main aim of the current study was to determine structure–function relationships in a large sample of SZ (n = 96) and HC (n = 142) by combining patterns of fMRI activation during a working-memory paradigm (n-back) and key brain morphometric properties including vertex- and voxel-based measures of surface area, cortical thickness, and gray matter volume. We used linked independent component analysis (LICA), a data-driven approach in which several imaging modalities may be combined (Groves et al., 2011; Groves et al., 2012), and tested for associations between LICA components reflecting brain structural features and functional brain networks related to task-positive and task-negative activation, respectively, and to what degree these associations were different in SZ compared to HC. Based on the few previous studies examining associations between task-related fMRI activation and gray matter structure in SZ (Calhoun et al., 2006; Fusar-Poli et al., 2011a; Harms et al., 2013; Michael et al., 2011; Pujol et al., 2013; Rasser et al., 2005; Schultz et al., 2012b), we hypothesized that strong task-related activation and deactivation would be associated with brain patterns reflecting increased structural integrity, including cortical thickness and gray matter volume in overlapping brain regions. Further, we expected that any group differences in the structure–function relationships would reflect disruptions of fronto-parietal, default-mode, and fronto-temporal brain regions in SZ. Lastly, likely related to larger between-subject variance, we expected stronger structure–function associations in SZ compared to HC.

Material and methods

Sample

Table 1 summarizes participant demographics and clinical variables. 238 participants, overlapping with Brandt et al. (2014, 2015), were recruited as part of the Thematically Organized Psychosis (TOP) study, comprising 96 DSM-IV-diagnosed patients with schizophrenia spectrum disorders (70 schizophrenia, 15 schizoaffective disorder, 11 schizophreniform disorder), referred to as SZ, and 142 healthy controls (HC).

Table 1

Demographic and clinical characteristics.

	SZ	HC	Test	p
Demographics
Sex (male); n (%)	61 (63.5)	75 (52.8)	χ2 = 2.7	.101
Age (years); mean (s.d.)	32.0 (8.2)	35.2 (9.0)	t = 2.8	.007
Handedness (right); n (%)	83 (86.5)	133 (93.7)	χ2 = 3.5	.060
Education (years); mean (s.d.)a	12.9 (2.5)	14.3 (2.3)	t = 4.3	<.001
IQ score; mean (s d.)b	104.5 (14.9)	114.9 (9.8)	t = 6.4	<.001
Duration of illness (years); mean (s.d.)c	6.5 (6.8)	–	–	–
Comorbid disorders; n (%)
Substance used	24 (25.0)	–	–	–
Somatic illnesse	15 (16.7)	–	–	–
Lifetime episodes; n (%)f
Psychosis	96 (100)	–	–	–
Depression	51 (53.1)	–	–	–
Mania	8 (8.3)	–	–	–
Current symptoms; n (%)g
Psychotic symptoms	29 (30.9)	–	–	–
Elevated mood symptoms	11 (11.7)	–	–	–
Depressive symptoms	25 (26.6)	–	–	–
Medicationh
Antipsychotics
n (%)	65 (73.9)	–	–	–
DDD; mean (s.d.)	1.1 (1.2)
Antiepileptics
n (%)	7 (7.8)	–	–	–
DDD; mean (s.d.)	.06 (.30)
Antidepressants
n (%)	20 (22.7)	–	–	–
DDD; mean (s.d.)	.31 (.67)
Anxiolytics
n (%)	7 (8.1)	–	–	–
DDD; mean (s.d.)	.07 (.29)
Substance usei
Alcohol use (AUDIT score); mean (s.d.)	5.8 (6.4)	5.3 (3.1)	t = 0.9	.361
Illicit drug use (DUDIT score); mean (s.d.)	2.9 (6.4)	.33 (1.5)	t = 4.5	<.001
Smoking; n (%)j	36 (48.6)	n.a.	–	–

SZ, schizophrenia; HC, healthy controls; DDD, defined daily dose; AUDIT/DUDIT, Alcohol/Drug Use Disorders Identification Test.

The total number of years of completed education as reported by the participants.

Wechsler Abbreviated Scale of Intelligence. Missing in SZ group: n = 4.

Number of years between age at onset and age at fMRI scanning. Age at onset was defined as age at first contact with the mental health service due to a primary symptom (n = 97) or age at first experience of symptoms (n = 2).

Lifetime abuse/dependency diagnosis of alcohol/cannabis/other drugs: 17/17/13%.

Lifetime somatic illness, included cardiovascular (2%), respirational (9%), endocrinological (1%), neurological (1%), or cancer (0%). Missing: n = 6.

Lifetime psychotic/depressive/manic episode, based on the SCID-interview (n = 99/90/98), age at first contact with the mental health service due to an episode (n = 0/6/0), or age at first experience of SCID-verified symptoms of an episode (n = 0/1/1).

Missing: n = 2.

Defined daily dose. Missing: antipsychotics n = 8, anti-epileptics n = 6, antidepressants n = 8, and anxiolytics n = 10.

Missing in SZ/HC groups: AUDIT n = 4/1, and DUDIT n = 4/2.

Daily smoking (yes/no) in the previous year. Missing: n = 22.

Patients were recruited consecutively from psychiatric units at four major hospitals in Oslo. Healthy controls were randomly selected from the same catchment area as the patient group using statistical records. 28% replied to the invitation and consented to participate, and of these 37.4% participated in the MR scanning. The study is approved by the Regional Committee for Medical Research Ethics and the Norwegian Data Inspectorate. Written informed consent was obtained from all participants according to the Declaration of Helsinki.

All participants had to be aged 18–65 years and speak a Scandinavian language. Exclusion criteria were presence of a developmental disorder or serious brain damage, and having metal implants, cardiac pacemaker or other MRI contraindications. Healthy controls were screened with a questionnaire about severe mental illness and the Primary Care Evaluation of Mental Disorders (Spitzer et al., 1994), and were excluded if they or a first-degree relative had a lifetime history of SZ, bipolar disorder or major depression.

Patients were characterized through a clinical interview conducted by trained physicians or clinical psychologists, covering diagnostics, symptoms, cognition, drug use and medication status. Diagnosis was established using the Structured Clinical Interview for DSM-IV Axis I Disorders (SCID) (First et al., 1995). Diagnostic reliability was satisfactory, with overall agreement for DSM-IV diagnosis categories of 82% and overall kappa = 0.77 (95% CI: 0.60–0.94).

Current symptoms were assessed on the day of scanning using a brief interview covering psychosis, depression, and elevated mood in the previous week (Aminoff et al., 2011). Symptoms were rated as present, possibly present or absent. Current psychotic symptoms were present in less than one third of the patients, indicating low symptom levels at the time of scanning. Regular medication was also recorded on the scanning day, as well as the use of alcohol and illicit drugs.

Experimental paradigm

The experimental paradigm was an n-back task with consecutive presentations of pairs of numbers between 1 and 9 (Haatveit et al., 2010; Hugdahl et al., 2004). In a 0-back condition, participants were instructed to press a response button when the two numbers were identical. In a 2-back condition, the numbers in each stimulus pair were identical and participants were instructed to press a response button when they were the same as the ones presented two trials earlier. The two conditions were separate runs, and all participants performed 0-back before 2-back. Stimuli were presented using a blocked design with four on-blocks and four off-blocks. On-blocks of 52 s duration comprised 18 stimuli presented in a pseudo-randomized sequence, including 3–4 targets (in total 12 targets in 0-back and 13 targets in 2-back). Stimulus duration was 300 ms and the inter-stimulus interval was 2500 ms. On-blocks were followed by off-blocks consisting of a fixation cross of 26 s duration. The paradigm is identical to the one used in Brandt et al. (2015).

The paradigm was programmed in E-Prime (Psychology Software Tools, Inc., Sharpsburg, USA) and stimuli were presented through goggles (NordicNeuroLab Inc., Bergen, Norway). Responses were collected using the ResponseGrips system (NordicNeuroLab Inc., Bergen, Norway). Participants responded with their right or left index finger or thumb. All participants completed a training procedure to ensure they had understood the task instructions.

MRI acquisition

MRI data were acquired on a 1.5 T Siemens Magnetom Sonata (Siemens Medical Solutions, Erlangen, Germany) using a standard head coil at Oslo University Hospital.

Structural data were acquired using a 3D T1-weighted magnetization prepared rapid acquisition gradient echo (MPRAGE) sequence with the following parameters: TR/TE/TI/FOV/FA/matrix = 2730 ms/3.93 ms/1000 ms/240 mm/7°/192 × 256; voxel size 1.33 × 0.94 × 1 mm, 160 sagittal slices. The sequence was repeated twice.

Functional imaging with 164 whole brain volumes per run was obtained with a T2*-weighted echo-planar imaging (EPI) pulse sequence with the following parameters: TR/TE/FOV/FA/matrix = 2040 ms/50 ms/224 × 224 mm/90°/64 × 64. Each volume consisted of 24 axial slices with a voxel size of 3 mm in the axial plane, and a slice thickness of 4 mm with 1 mm gap between slices. The first seven volumes and the last volume were discarded, leaving 156 volumes for analysis.

Quality control

An initial screening of the data was performed to exclude datasets with obvious movement or scanner artifacts. From the remaining available sample of 115 patients and 163 controls, 16 SZ and 20 HC were excluded due to task performance below chance level (SZ: n = 3), excessive relative head motion (displacement from 1 volume to the next) defined as more than 3 SDs above the mean (SZ: n = 6, HC: n = 3), and poor data quality (structural images; SZ: n = 3, HC: n = 1; signal loss in functional images; SZ: n = 7, HC: n = 16), yielding 96 SZ and 142 HC in the included sample.

Structural MRI processing

T1-weighted scans were processed using FreeSurfer (http://surfer.nmr.mgh.harvard.edu) including 3D surface reconstruction and full brain segmentation (Fischl et al., 2002). The segmented volume was used in order to produce high quality brain masks for co-registration.

FreeSurfer was also used to estimate vertex-wise cortical thickness and surface area measures (Dale et al., 1999; Fischl et al., 1999). Cortical thickness was obtained by reconstructing the gray/white matter boundary and the pial surface (Dale et al., 1999) and then calculating the distance between the surfaces at each vertex. Surface area was estimated by registering each subject's reconstructed surfaces to a common template, and the relative amount of expansion or compression at each vertex was used as a proxy for regional arealization. Surface maps were resampled and mapped to a common coordinate system (fsaverage5, 10,242 vertices) using a non-rigid high-dimensional spherical averaging method to align cortical folding patterns (Fischl et al., 1999).

FSL-VBM (Douaud et al., 2007), a voxel-based morphometry (VBM) analysis (Ashburner and Friston, 2000; Good et al., 2001), was used to derive a gray matter volume (GMV) map. The registered map was divided by the Jacobian map of the deformation field to account for local contraction and expansion.

fMRI processing

fMRI data were processed using FEAT, part of FMRIB's Software Library (FSL; http://www.fmrib.ox.ac.uk/fsl) (Jenkinson et al., 2012). Preprocessing included motion correction (Jenkinson et al., 2002), non-brain removal (Smith, 2002), spatial smoothing using a Gaussian kernel of full width of half maximum (FWHM) = 6 mm, and high-pass temporal filtering with a 90 s window. Registration from fMRI to structural space was carried out using FLIRT (Jenkinson et al., 2002), and fMRI data were warped to MNI space via the high-resolution structural volume using FNIRT (http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/FNIRT).

In a lower-level general linear model (GLM) analysis for each run (0-back and 2-back), the onset and duration of the on-blocks were modeled with the off-blocks as implicit baseline. The design matrix was filtered and convolved with a hemodynamic response function (HRF) before the model fit. Motion (translation and rotation) parameters as estimated by the motion correction during preprocessing were included in the model to remove residual motion effects that were left after correction, and a temporal derivative was added to adjust for regional differences in the timing of the HRF. This first-level analysis resulted in individual contrast parameter estimate (COPE) maps reflecting the on vs. off contrast within 0-back and 2-back. Then, an intermediate-level fixed-effects analysis was performed to calculate individual COPEs reflecting the 2-back vs 0-back contrast.

Downsampling and smoothing of structural and functional measures

Cortical thickness maps were smoothed using a Gaussian kernel with an FWHM of 15 mm, yielding an effective smoothness of 20.8 mm as estimated based on Hagler et al. (2006). To match the smoothness across measures, which may be beneficial in order to match the dimensionality across features, smoothing of the surface area, GMV, and COPE maps was performed such that the resulting effective smoothness was comparable to the thickness maps. The surface area maps were smoothed with an FWHM of 11 mm. GMV and COPE maps were resampled to have an isotropic voxel size of 4 mm using a trilinear interpolator, and subsequently smoothed using FreeSurfer's tool mri_fwhm with an FWHM of 10.5 and 9.8 mm, respectively.

Whereas the main analysis was performed using matched smoothness across modalities, the necessity and utility of assuming a common smoothing level across modalities in a LICA context is unresolved. Therefore, LICA and the main between-subjects analyses were rerun with no additional smoothing of the COPE maps and less smoothing of GMV maps (using the same FWHM value of 9.4 mm as described by Groves et al. (2012)). The purpose was to assess the robustness of the results across smoothing levels, and specifically test whether significant findings from the main analysis were replicated in the additional analysis.

Linked independent component analysis

Multimodal analysis was performed using LICA (Groves et al., 2011, 2012) including three structural and three functional modalities: cortical surface area, thickness, and GMV, as well as full-brain and unthresholded individual COPE maps for 0-back (task vs. rest), 2-back (task vs. rest), and 2-back > 0-back.

Dimensionality selection was automatically performed by LICA (Groves et al., 2011), yielding 66 components (83 components for the additional analysis with less smoothing). In LICA, each component has a shared individual subject weight (loading) across modalities, indicating the degree to which each subject contributes to that component. In addition, all components have corresponding spatial maps for each modality. The spatial maps of all components were inspected, and the relative weight of each modality in each component was assessed.

No components showed substantial fusing between structural and functional modalities, i.e. more than 10% involvement of at least one structural and one functional modality. 27 components were defined as structural (<10% total involvement of functional modalities) and 39 as functional (<10% total involvement of structural modalities).

The structural components, including the modalities and brain regions involved, are summarized in Table A.1. Based on the spatial maps, two functional components were selected for correlation with these structural components: A task-positive fronto-parietal network (FPN; ic2) overlapping with a central executive network (Bressler and Menon, 2010; Sridharan et al., 2008), and a task-negative default-mode network (DMN; ic12) (Cole et al., 2010; Smith et al., 2009). Both components showed strong contributions from 2-back and 2-back > 0-back COPE maps, and were chosen based on their spatial similarity with canonical task-positive and task-negative patterns of brain activation (Cocchi et al., 2013; Raichle et al., 2001; Spreng et al., 2010), as well as the variance explained by these two components within the 2-back modality across subjects. The spatial maps of the FPN and DMN LICA components overlapped highly with brain regions showing positive and negative correlations with the task, respectively (Fig. 1), which was confirmed by the correlations between the spatial maps of FPN/DMN components and 2-back COPE-maps, as estimated using fslcc in FSL (FPN vs. 2-back > rest: r = .76; DMN vs. rest > 2-back: r = .66). FPN explained 17.1% (ranked first, i.e. being the component explaining the largest amount of variance) and DMN explained 6.1% (ranked fourth) of the variance, i.e. together these components explained 23.2% of the total variance across the brain in the 2-back modality. Other functional components related to 2-back reflected brain networks such as left (6.3%) and right (4.9%) fronto-parietal, cingulo-opercular (4.3%), parietal (dorsal attention: 3.9%), and occipital (primary visual: 4.2%; secondary visual: 3.6%) networks.

Fig. 1

The task-positive fronto-parietal (FPN; red) and task-negative default-mode (DMN; blue) functional components from LICA (z > 3, uncorrected), and GLM-maps of activation (2-back > rest) and deactivation (rest > 2-back; z > 2, uncorrected), overlapping with functional components from LICA.

Statistical analysis

Task performance

N-back performance was assessed using d-prime (d′) (see Haatveit et al., 2010) and response time (RT) on correct responses. Repeated measure ANOVAs were performed with load (0-back, 2-back) as within-subject factor and diagnosis (SZ, HC) as between-subject factor to test for main effects of load and diagnosis and their interactions on d′ and RT.

Functional components

In all analyses, the functional and structural components identified using LICA were investigated by using the individual subject weights as input. Linear regressions were performed to investigate group differences within FPN and DMN while covarying for sex and age. Associations between the two functional components were assessed using partial correlations across groups (covarying for sex, age, and diagnosis) and within groups (covarying for sex and age). We tested for associations with task performance (d′ and RT during 2-back) and relative subject motion during 2-back on the two functional components while covarying for sex, age, and diagnosis.

Structure–function associations

The functional components showed no clear multi-modal fusion with structural modalities, and structure–function relationships were assessed based on the correlations between the subject weights of the two functional components and the structural components.

In order to investigate the associations between structure and function across groups, linear regressions were performed with each structural component as independent variable and each functional component as dependent variable, using two different models: (1) while covarying for sex, age, and diagnosis, and (2) while covarying for sex, age, diagnosis, and the interaction between structural component and diagnosis. For associations that did not show a significant interaction effect, the results are reported and interpreted from the across groups analysis without interaction term. For cases with a significant interaction effect, results are reported from the across groups analysis with interaction term, but the structure–function association within each group is used for interpretation of the results.

Between groups (case–control) differences in structure–function associations were inferred from the interaction term (structural component × diagnosis) from the same statistical model. Associations within groups are also reported (covarying for sex and age).

Comparison with conventional analysis

To assess the methodological generalizability of the findings, we tested the associations between the two functional components and several conventional structural measures: (1) surface area (vertex-based analysis using FreeSurfer), (2) cortical thickness (vertex-based analysis using FreeSurfer), (3) GMV (voxel-based analysis using FSL-VBM), and (4) selected global and regional brain volumes, including intracranial volume, total cortical volume, total subcortical gray matter volume, as well as bilateral volumes of the lateral ventricles, hippocampus, amygdala, thalamus, caudate, putamen, pallidum, corpus callosum, and cerebellar cortex. Here, using linear regressions, we tested the effect of each functional component on surface area, thickness and GMV, while we tested the effect of each brain structural volume on the functional components, both across groups (covarying for sex, age, and diagnosis) and within groups (covarying for sex and age). For the associations with surface area and structural volumes, intracranial volume was also included as covariate.

Associations with possible confounders

For structure–function associations surviving correction for multiple comparisons within patients, post-hoc tests were performed to examine the effects of possible confounders, including current symptoms (psychosis, depression, elevated mood), lifetime symptoms (depression, mania), medication level (antipsychotics, antiepileptics, antidepressants, anxiolytics), duration of illness, and substance use (alcohol use, illicit drug use, drug disorder, smoking). Each covariate was included sequentially in the within-group model in addition to sex and age.

For associations surviving correction for multiple comparisons across and between groups, as well as within SZ, the analyses were also performed after excluding patients with schizoaffective (n = 15) and schizophreniform (n = 11) disorders. Further, effects of task performance (d′ and RT during 2-back) and relative subject motion during 2-back were tested by including these variables in the model.

Statistical threshold

For the main analysis of structure–function associations across, between, and within groups, correction for multiple comparisons was performed using permutation testing (Nichols and Holmes, 2002). The subject weights of each of the two functional components (FPN, DMN) were permuted together with the covariates (sex, age, diagnosis) 10,000 times with respect to the subject weights of the structural components. In each iteration, the maximum positive and minimum negative t statistics across the structural components were computed and stored, resulting in a two-tail null distribution of the t-statistics. The corrected p-value was obtained by comparing the t-statistics obtained from the original analyses to the null distribution.

The main focus was on associations that were significant at a corrected p-level, but effects with a nominal p < .05 (uncorrected) are also reported for full transparency and in order to facilitate comparisons with previous and future studies. For other analyses, the statistical threshold was set to p < .05 (uncorrected).

All statistical analyses were carried out using IBM SPSS Statistics version 22 and R (http://www.r-project.org).

Results

Task performance

Table 2 summarizes task performance. In line with previous reports from an overlapping sample (Brandt et al., 2014; Brandt et al., 2015), significant group differences in d′ and RT were found in both load conditions, indicating reduced target discrimination and increased response times in SZ. In addition to main effects of load and diagnosis on d′ (load: F = 204.8; diagnosis: F = 37.4; p < .001) and RT (load: F = 98.2; diagnosis: F = 12.6; p < .001), there was also a significant load × diagnosis interaction on both measures (d′: F = 38.0, p < .001; RT: F = 8.7, p = .004), indicating larger group differences during 2-back compared to 0-back.

Table 2

Task performance.

	SZ	HC	t	p
0-back
Accuracy — % hits (s.d.)	99.4 (1.4)	99.7 (0.9)	2.1	<.05
d-prime — mean (s d.)a	3.98 (.30)	4.04 (.18)	2.0	<.05
RT hits, ms — mean (s d.)b	551.1 (118.7)	521.0 (81.5)	2.2	<.05
RT total, ms — mean (s d.)b	550.6 (120.5)	521.0 (81.9)	2.2	<.05
2-back
Accuracy — % hits (s.d.)	93.9 (6.3)	97.7 (3.0)	6.1	<.001
d-prime — mean (s d.)c	3.06 (.90)	3.68 (.59)	6.4	<.001
RT hits, ms — mean (s d.)b	683.5 (213.6)	592.7 (150.8)	3.7	<.001
RT total, ms — mean (s d.)b	697.3 (210.5)	602.6 (150.6)	3.9	<.001

SZ, schizophrenia; HC, healthy controls; RT, response time.

Max score: 4.13.

Missing in SZ/HC groups: 14/5.

Max score: 4.16.

Functional components

Fig. 1 shows spatial maps of the FPN and DMN components, as well as the 2-back vs. rest COPE maps. There was no significant (t = 1.3, p = .196) difference in FPN subject weights between SZ (mean = −.090, SD = .99) and HC (mean = .057, SD = 1.00) covarying for age and sex. Similarly, there was no significant (t = −.08, p = .443) difference in DMN between SZ (mean = .035, SD = .99) and HC (mean = −.016, SD = 1.01). FPN and DMN subject weights were negatively correlated (across groups: r = −.46, p < .001; within SZ: r = −.40, p < .001, within HC: r = −.51; p < .001), indicating increasing task-positive activations with increasing deactivation of the DMN. d′ in 2-back was associated with both functional components (FPN: t = 5.5, p < .001; DMN: t = −3.4, p = .001, covarying for sex, age, diagnosis), reflecting increasing FPN activation and DMN deactivation with increasing performance. RT in 2-back was associated with FPN (t = −3.0, p = .004, covarying for sex, age, diagnosis), but not DMN (t = .80, p = .424), indicating increasing FPN activation with slower responses. Subject motion showed a significant association with both functional components (FPN: t = −2.4, p = .017; DMN: t = 2.2, p = .027), indicating reduced FPN activation and reduced DMN deactivation with increasing motion. These results (d′, RT, motion) were highly similar also when not including diagnosis in the model.

Structure–function associations

Across groups

Structure–function associations across groups are shown in Table A.2. Permutation testing revealed a significant (t = 3.0, uncorrected p = .003, corrected p < .05) association across groups between the FPN component and a component reflecting fronto-temporal thickness (ic60), indicating stronger FPN activation with increasing thickness.

Eight other structure–function associations were nominally significant (p < .05, uncorrected) across groups: FPN was associated with components reflecting global surface area (ic4; t = 2.2, p = .031), global thickness (ic5; t = 2.2, p = .029), and precentral surface area (ic62; t = 2.7, p = .009), while DMN was associated with cerebellar GMV (ic22; t = −2.3, p = .023), frontal thickness (ic45; t = 2.2, p = .031), temporal GMV (ic50; t = 2.5, p = .013), precuneal GMV (ic58; t = −2.1, p = .041), and fronto-temporal thickness (ic60; t = −2.0, p = .047).

Between groups

Table A.2 also shows the case–control differences in structure–function associations, as indicated by the interaction term. None of these remained significant after permutation testing, but the association between DMN and a structural component reflecting frontal GMV (ic59) was at the border of significance (t = −2.8, uncorrected p = .005; corrected p = .066), indicating reduced DMN deactivation with increased frontal GMV in controls (t = 2.4, uncorrected p = .020), and a tendency towards an opposite pattern in SZ (t = −1.8, uncorrected p = .078).

Four other structure–function associations showed nominally significant (p < .05, uncorrected) differences between SZ and HC: FPN and fronto-temporal thickness (ic60; t = 2.6, p = .009), DMN and cerebellar GMV (ic22; t = −2.2, p = .026), DMN and fronto-cingulate surface area (ic35; t = −2.1, p = .040), and DMN and precuneal GMV (ic58; t = −2.0, p = .046). These group differences were mainly driven by significant associations within patients, indicating increased FPN activation with increased fronto-temporal thickness (ic60), and increased DMN deactivation with increased cerebellar (ic22) and precuneal (ic58) GMV. DMN and fronto-cingulate surface area (ic35) were not significantly associated in any of the groups.

Within groups

Table A.3 shows structure–function associations within each group. In SZ, the association between FPN and fronto-temporal thickness (ic60) was significant at a corrected p-level threshold (t = 3.0, uncorrected p = .003, corrected p = .038).

Nominal associations (p < .05, uncorrected) were found between FPN and global surface area (ic4; t = 2.3, p = .027), global thickness (ic5; t = 2.1, p = .041), and occipital surface area (ic66; t = −2.0, p = .048). DMN was nominally associated with global surface area (t = −2.4, p = .017), cerebellar GMV (ic22; t = −2.3, p = .024), frontal thickness (ic45; B = .28, t = 2.5, p = .015), temporal GMV (ic50; t = 2.1, p = .043), cerebellar-putamen GMV (ic54; t = −2.0, p = .049), precuneal GMV (ic58; t = −2.1, p = .035), and fronto-temporal thickness (ic60: t = −2.1, p = .035).

Within HC, no associations survived permutation testing. Nominal associations (p < .05, uncorrected) were found between FPN and temporal GMV (ic50; t = −2.0, p = .045) and precentral surface area (ic62; t = 2.2, p = .028), while DMN was associated with temporo-parietal surface area (ic39; t = −2.0, p = .048) and frontal GMV (ic59; t = 2.4, p = .020).

Comparison with conventional analysis

Across groups, the analysis of conventional structural measures showed associations (p < .05, uncorrected) between FPN and surface area in small clusters mainly in the left precentral, postcentral and inferior temporal gyri; cortical thickness mainly in the precentral, postcentral, and superior temporal gyri (Fig. A.1); and GMV in small clusters in the postcentral gyrus, precuneus, and lateral occipital cortex. DMN was associated (p < .05, uncorrected) with surface area in small clusters mainly in the medial frontal pole, precuneus, and superior parietal cortex; cortical thickness in the inferior frontal, precentral, postcentral, and superior temporal gyri; and GMV in small clusters in the frontal pole, occipital cortex, and middle frontal gyrus. As shown for cortical thickness in Fig. A.1, these associations were generally more pronounced in patients than in controls, with a main pattern of positive associations with FPN and negative associations with DMN.

The analyses of brain volumes (Table A.4) showed that across groups FPN was associated with total cortical volume (t = 2.0, p = .046, uncorrected), while there were no associations with any of the other brain structural volumes, and no associations between DMN and brain volumes. Within SZ, FPN was associated with total cortical volume (t = 2.2, p = .033, uncorrected) and total subcortical gray matter volume (t = 2.2, p = .034, uncorrected), while DMN was associated with cerebellar volume (t = −3.0, p = .003, uncorrected). In HC, no associations were found between functional components and brain volumes.

Associations with possible confounders

The association between FPN and fronto-temporal thickness (ic60), which was significant within SZ also at a corrected p-level, remained significant at a nominal level when controlling for the effect of all clinical variables (FPN–ic60: p < .01, uncorrected), and when including d′ (FPN–ic60: p = .013, uncorrected), RT (FPN–ic60: p = .010, uncorrected), and relative subject motion (FPN–ic60: p = .002, uncorrected) in the model. When excluding patients with schizoaffective (n = 15) and schizophreniform (n = 11) disorders from the analysis the FPN-ic60 association within patients (n = 70) was not significant, but showed a tendency towards a positive association (B = .15, t = 1.4, p = .178).

Consistency across smoothing levels

In line with the results from the main analysis using matched smoothness, the additional analysis with less smoothing on the 3D GMV and COPE maps showed a significant association between the FPN resulting from this LICA run (ic2) and a component reflecting fronto-temporal thickness (ic64) in SZ, which also survived correction for multiple comparisons (t = 3.2, uncorrected p = .002, corrected p = .038). See Table A.5 for all within groups associations from this analysis.

As shown in Fig. 2, the functional FPN components were highly similar across analyses (spatial map correlation: r = .93; subject weight correlation: r = .95), including the frontal pole, middle frontal gyrus, inferior frontal gyrus, insula, precentral gyrus, paracingulate gyrus, supplementary motor area, lateral occipital cortex and superior parietal lobule in both hemispheres. The somewhat more extended pattern in the original analysis was due to the higher smoothing level. The structural fronto-temporal thickness components were also similar across analyses (spatial map correlation: r = .59; subject weight correlation: r = .59), including the superior temporal gyrus, insula, and postcentral gyrus in both hemispheres. The inferior frontal and precentral gyri showed a stronger involvement in the original analysis, and the insula in the additional analysis, but all these regions were present in both structural components at a z-score threshold <3, indicating that the components captured the same spatial variability.

Fig. 2

The task-positive fronto-parietal network (FPN) and fronto-temporal thickness (FTT) components showing a significant association in schizophrenia, from (A) the main analysis using matched smoothness across modalities, and (B) the additional analysis using less smoothing of COPE and GMV maps (3 < z < 10).

Discussion

Using a multivariate approach for multimodal fusing of brain imaging data, we have documented a robust association between task-positive fronto-parietal brain activation and cortical thickness of fronto-temporal regions in SZ. This finding was also confirmed in an additional analysis using less smoothing of fMRI and GMV maps. Whereas the remaining structure–function associations were in general stronger in patients compared to controls, they were moderate both across, between and within groups, with nominally significant associations. The characteristics and implications of these novel findings will be detailed below.

Fronto-parietal brain activation and fronto-temporal thickness

The strongest and most consistent structure–function relationship was found in SZ, between the functional fronto-parietal network (FPN) and a structural component reflecting cortical thickness of frontal and temporal brain regions, including the insula. This association was specific to SZ, indicating a disrupted integration of brain structure and function in the implicated brain regions.

The current finding is in line with previous reports of fronto-temporal and insular gray matter reductions in SZ (Glahn et al., 2008; Nesvag et al., 2008; Rimol et al., 2010; Schultz et al., 2012b; Shepherd et al., 2012), which may reflect alterations at the neuronal and synaptic level, with consequences for cognitive networks and processing. Evidence from DTI studies suggests disruptions of fronto-temporal white matter bundles in SZ, including the uncinate fasciculus (Kubicki et al., 2002; Samartzis et al., 2014). Similarly, functional imaging studies have reported connectivity alterations between frontal and temporal areas (Cocchi et al., 2014; Crossley et al., 2009; Meyer-Lindenberg et al., 2001; Meyer-Lindenberg et al., 2005; Wolf et al., 2007), as well as regional abnormalities in frontal (Callicott et al., 2003; Glahn et al., 2005; Thormodsen et al., 2011), temporal (Hugdahl et al., 2009), and insular (Palaniyappan and Liddle, 2012; Uddin, 2015) regions during cognitive processing. However, the current finding extends previous reports of fronto-temporal abnormalities in SZ by suggesting a coupling between cortical thickness in these brain regions and working memory-related brain activation. The regions involved in the structural and functional components were partly overlapping, e.g. in the insula which shows consistent abnormalities in SZ and has been suggested a role in hallucinations (Palaniyappan and Liddle, 2012). However, the main pattern emerging was a lack of anatomical overlap, indicating that activation of the FPN in patients was related to brain structural variability in other regions. Yet the causal relationship is unclear (see Harms et al., 2013 for discussion).

A few other studies have also reported structure–function alterations in SZ in similar regions using different methodological approaches. Correlating gray matter and fMRI voxels across the whole brain, Michael et al. (2011) reported a differential and aberrant association in SZ between working memory-related brain activation and fronto-temporal anatomy, though this correlation was negative, indicating decreased activation with increased gray matter volume (Michael et al., 2011). In line with the current results, Pujol et al. (2013) reported that cortical thickness in the left inferior frontal gyrus and insula explained 57% of the variability in task-positive network activation and task-negative deactivation in SZ. However, this study comprised a small number of patients (n = 14) and restricted the analyses to regions showing cortical thinning in the patient group. The current positive association between FPN and fronto-temporal thickness indicates increased activation with increased thickness. In line with this, a positive correlation between anterior cingulate working memory-related activation during a Sternberg task and fronto-temporal thickness in SZ has been reported (Schultz et al., 2012b).

The observed association between FPN and fronto-temporal thickness in SZ may be related to several factors which should be taken into consideration when interpreting the results, including behavioral effects, clinical variables and patient subgroups. The structure–function association remained significant when controlling for task performance and in-scanner subject motion. Still, behavioral differences within patients may have influenced the findings indirectly, since differences in task performance may reflect illness severity and other aspects of the disorder. Similarly for clinical variables, the association was also present when controlling for the effect of illness duration, current and lifetime symptoms, medication level, and substance use. However, since this is a naturalistic study, there is an inherent association between clinical severity, symptoms, and medication status, which is difficult to disentangle. Also, since all patients were medicated, it is not possible to isolate effects of disease from effects of medication.

In line with this, a recent study (Lesh et al., 2015) showed aspects of this complexity by reporting that cortical thinning, including in fronto-temporal regions, was identified in patients receiving antipsychotic treatment, but not in unmedicated patients, as compared to healthy controls. The medicated group also presented with higher prefrontal activation and better behavioral performance than the unmedicated group in response to a continuous performance task (Lesh et al., 2015). Due to this complex relationship between medication status and brain alterations in SZ it is unclear to which degree the current findings reflect brain abnormalities related to vulnerability and underlying pathophysiology, or secondary disease- and treatment-related effects.

When excluding patients with schizoaffective and schizophreniform disorders, the structure–function association in patients largely disappeared, implying that the effect was partly driven by these subgroups. However, since these groups consisted of 15 and 11 patients only, no clear conclusions should be drawn from this, and further studies are needed to assess the reliability of the present findings.

Other structure–function associations

Across and between groups, no associations survived correction for multiple comparisons using permutation testing. The results suggest a complex relationship between brain structure and function, which does not support a direct anatomical overlap. Instead, as also shown in the main finding, functional networks such as FPN and DMN may rely on brain structural variability in distant regions (Harms et al., 2013), which is in line with a system-level perspective of brain function. Further, whereas all effect sizes were moderate, patients in general showed more numerous and stronger associations than HC, in line with recent reports (Fusar-Poli et al., 2011a; Pujol et al., 2013; Schultz et al., 2012b). This may be related to factors such as stronger variance in patients and clinical characteristics that are specific to the patient group.

Methodological aspects and possible confounders

In addition to clinical variables, several factors should be mentioned in relation to the current results. Compared to HC, SZ performed significantly worse on the 2-back task, implying reduced target discrimination and slower responses, while there was no group difference in FPN and DMN activations. Increased target discrimination was associated with increased FPN activation and DMN deactivation across groups, while slower responses were associated with increased FPN activation. There was also an effect of subject motion on the functional components, indicating reduced FPN activation and DMN deactivation with increasing subject motion. These behavioral effects on brain activation may also have influenced the observed structure–function associations and the lack thereof.

Further, it is not clear to what degree varying spatial smoothing levels of imaging modalities influence LICA results. Therefore, as a supplement to the main analysis involving comparable smoothness across modalities, we performed an additional analysis using less smoothing on GMV and fMRI maps, in order to assess consistency across analyses. The results showed that the main finding of an association between FPN and fronto-temporal thickness in patients was reliable across approaches and independent of smoothing level, suggesting that robust associations are not sensitive to such variation. However, since issues related to spatial smoothing are unresolved, future studies should assess the impact of spatial smoothing on the results, and also consider other approaches such as surface-based alignment, which has been shown advantageous in SZ (Anticevic et al., 2008).

An important question related to the generalizability of the findings is whether LICA is sensitive to the same variability as other approaches. First, this question is relevant to the functional LICA components used here and to what degree these are comparable to networks resulting from, e.g. independent component analysis of task- and resting-state fMRI data, and from conventional fMRI analysis. Encouragingly, even in absence of any time series information in the COPE maps submitted to LICA, the FPN and DMN components were isolated as separate components that were independent of components reflecting other brain networks (e.g. visual and lateralized fronto-parietal networks). FPN and DMN also showed a high spatial overlap with the group COPE maps reflecting task vs. rest. Second, although univariate and multivariate approaches may in principle target different sources of variability, conventional analyses of structural imaging data yielded an overlapping pattern of nominally significant structure–function associations and slightly stronger effects in SZ than HC. Similar patterns and directions (positive/negative) were found for several structural components, including the FPN association with fronto-temporal thickness in SZ. Together, these results show that the current approach detects patterns of structure–function relationships which are consistent with conventional structural approaches. The advantages of LICA include that the variation in other components is taken into account when a particular component is investigated, which allows for stronger interpretations regarding specificity. Also, it provides a substantial down-sampling of the multidimensional imaging data, which increases statistical power and facilitates further associations with other variables, including clinical and genetic variability.

Conclusion

Using a multimodal fusion approach to investigate structure–function associations in a large sample of SZ and HC, we found a significant and robust relationship between a task-positive fronto-parietal network and fronto-temporal thickness in SZ. This finding extends previous reports of fronto-temporal abnormalities in SZ by suggesting a coupling between cortical thickness in these brain regions and working memory-related brain activation.