A task-general connectivity model reveals variation in convergence of cortical inputs to functional regions of the cerebellum

  1. Maedbh King  Is a corresponding author
  2. Ladan Shahshahani
  3. Richard B Ivry
  4. Jörn Diedrichsen
  1. Department of Psychology, University of California, Berkeley, United States
  2. Western Institute for Neuroscience, Western University, Canada
  3. Helen Wills Neuroscience Institute, University of California, Berkeley, United States
  4. Department of Statistical and Actuarial Sciences, Western University, Canada
  5. Department of Computer Science, Western University, London, Canada

Abstract

While resting-state fMRI studies have provided a broad picture of the connectivity between human neocortex and cerebellum, the degree of convergence of cortical inputs onto cerebellar circuits remains unknown. Does each cerebellar region receive input from a single cortical area or convergent inputs from multiple cortical areas? Here, we use task-based fMRI data to build a range of cortico-cerebellar connectivity models, each allowing for a different degree of convergence. We compared these models by their ability to predict cerebellar activity patterns for novel Task Sets. Models that allow some degree of convergence provided the best predictions, arguing for convergence of multiple cortical inputs onto single cerebellar voxels. Importantly, the degree of convergence varied across the cerebellum with the highest convergence observed in areas linked to language, working memory, and social cognition. These findings suggest important differences in the way that functional subdivisions of the cerebellum support motor and cognitive function.

Editor's evaluation

This study presents the valuable finding that the human cerebellum receives highly convergent (rather than sparse) task-related connectivity from the cortex. A compelling case is made that this convergent connectivity supports a wide variety of cognitive functions, though it will be important for future work to fully verify the cortical-to-cerebellar directionality of the results. The work will be of broad interest to cognitive neuroscientists interested in the role of connectivity in cognition.

https://doi.org/10.7554/eLife.81511.sa0

Introduction

The last 30 years has witnessed a paradigm shift with regard to cerebellar function, with broad recognition that this subcortical structure is engaged in many aspects of human cognition. Since the first report of cerebellar activation in a semantic retrieval task (Petersen et al., 1989), thousands of neuroimaging papers have reported cerebellar recruitment during a broad range of tasks that cannot be attributed to the motor demands of these tasks. Functional interpretations include hypotheses concerning how the cerebellum may facilitate attentional shifts (Allen et al., 1997), stimulus-response mapping (Bischoff-Grethe et al., 2002), higher order rule processing (Balsters et al., 2013), verbal working memory (Marvel and Desmond, 2010), language (Fiez, 2016), and social cognition (Van Overwalle et al., 2015). This body of work has produced functional maps of the cerebellum that depict the association of particular cognitive processes with different subregions of the cerebellum (King et al., 2019).

Given the relative uniform cytoarchitecture of the cerebellar cortex, it is assumed that differences in function mainly arise from variation in the input to the cerebellum. Trans-synaptic tracing methods employed in non-human primates studies have revealed extensive reciprocal connections between many frontal and parietal areas and the cerebellum (Dum and Strick, 2003; Kelly and Strick, 2003). These studies have highlighted the closed-loop nature of these connections, with each (neo-)cortical region projecting to a specific cerebellar region, and receiving input from the same area (Strick et al., 2009). In humans, resting state functional connectivity analyses have revealed a set of cerebellar networks, each one associated with a specific cortical network (Buckner et al., 2011; Ji et al., 2019; Marek et al., 2018).

An important unanswered question is whether each cerebellar region receives input from a restricted cortical region or whether it receives convergent input from multiple cortical regions. Providing an answer to this question has important implications for our understanding of cerebellar function. An architecture marked by a one-to-one relationship between cortical and cerebellar regions would suggest that the function of each cerebellar region is to fine-tune the dynamics in its cortical input. In contrast, a convergent architecture would suggest that subregions within the cerebellum integrate information across disparate cortical regions and may coordinate their interactions. Indeed, recent work in the rodent brain has suggested convergence of mossy fibers from diverse sources onto the same cerebellar region (Henschke and Pakan, 2020; Pisano et al., 2021), or even onto the same granule cells (Huang et al., 2013). Furthermore, the pattern of cortico-cerebellar convergence may vary across the cerebellar cortex, similar to how cortical areas show considerable variation in the degree to which they serve as points of convergence from other cortical regions (Bertolero et al., 2015; Yeo et al., 2014; Yeo et al., 2015).

In the current study, we introduce a novel approach to study cortico-cerebellar connectivity. Using the data from a large battery of tasks, we derived models that could be used to predict the activity in each cerebellar voxel based on the activity pattern in the neocortex. While this approach allowed us to evaluate the degree of convergence of cortical inputs to the cerebellar cortex, we recognize that the model could also be evaluated in the opposite direction, namely, to predict activity in the neocortex based on cerebellar activity patterns. However, we believe that using the model to make directional predictions from neocortex to cerebellum is most appropriate for fMRI data. The BOLD signal in the cerebellar cortex overwhelmingly reflects cortical input (via the pons) with no measurable contribution from the Purkinje cells, the output neurons of the cerebellar cortex (Mathiesen et al., 2000; Thomsen et al., 2004; Alahmadi et al., 2016; Alahmadi et al., 2015; Mapelli et al., 2017; Gagliano et al., 2022). In contrast, the neocortical BOLD signal reflects many sources including inputs from other cortical regions, local activity (input and output), as well as ascending input. While the latter will include cerebellar input via the thalamus, this source is likely to make a relatively small contribution to the overall BOLD response in the neocortex. As such, the relationship between neocortical and cerebellar BOLD signals will be most informative in evaluating cortico-cerebellar connectivity.

We trained multiple models of neocortical-cerebellar connectivity on a fMRI data set obtained while human participants completed a large task battery that was designed to engage cognitive processes across a broad range of functional domains (e.g. visual cognition, memory, attention, cognitive control). The models varied in terms of the degree of convergence of cortical inputs onto each cerebellar area. To evaluate the models in a cross-validated fashion, we examined how well each model predicted cerebellar data obtained from different tasks and/or different participants, using only the corresponding neocortical data. These analyses reveal a novel picture of cortico-cerebellar connectivity, one in which the degree of convergence was higher in cerebellar regions associated with more complex cognitive functions.

Results

Overview

We compared three models of cortico-cerebellar connectivity by using a region-to-region predictive modeling approach (Cole et al., 2016; Mell et al., 2021). For each model, the task-evoked activity in each cerebellar voxel was predicted as a linear combination of task-evoked activity patterns across the entire neocortex (Figure 1). As a model of sparse connectivity, we used a Winner-Take-All (WTA) model which imposes the strong constraint that only a single cortical region is used to predict the activity in each cerebellar voxel. The other two models allowed for some degree of convergence. Models estimated with Lasso Regression (L1 regularization) find sparse solutions, minimizing the number of cortical inputs by setting to zero those weights that make a negligible contribution to the predicted activity. Models estimated with Ridge Regression (L2 regularization) allow for a wide distribution of inputs, keeping each weight as small as possible. For model input, we used a set of cortical parcellations that varied in terms of their level of granularity. We trained the models on cortical and cerebellar fMRI data obtained from 24 participants who performed a task battery with 29 task conditions (Task Set A, Figure 1A). These were acquired over two sessions on separate days, with the tasks identical across sessions. These data were used to estimate the connectivity weights W^separately for each model and participant (see Methods).

Connectivity model training and evaluation.

(A) Models were trained on Task Set A (session 1 and 2) of the multi-domain task battery (MDTB; King et al., 2019). Model hyperparameters were tuned using fourfold cross validation. Three types of models were used (WTA, Lasso, Ridge), each with 7 cortical parcellations of different granularity. (B) Models were tested on an independent Task Set B (sessions 3 and 4), which included both novel and common tasks. Models had to predict the cerebellar activity solely from cortical activity patterns. To avoid the influence of shared noise correlations across cortex and cerebellum, the models were trained and tested by using cortical and cerebellar activity patterns from different sessions within each task set (see Methods). (C) As a test of generalization, the models were used to predict cerebellar activity from cortical data obtained from a separate experiment (King et al., unpublished data).

To compare the models, we tested their prediction performance on two independent datasets. First, we used data from the same 24 participants when tested in two different sessions (Task Set B) that included 18 novel conditions and 14 conditions repeated from Task Set A (Figure 1B). To predict the cerebellar activity patterns, we used the observed cortical activity patterns from Task Set B and the estimated connectivity weights (W^) from Task Set A for each participant. Note that because the model predictions relied on a single set of connectivity weights across all tasks, the input is based only on the cortical activity patterns without reference to any features of the tasks themselves. Second, we also tested how the models would generalize when tested with data from a new group of participants tested on a set of novel tasks (Figure 1C).

The use of separate training and evaluation datasets allowed us to determine the best task-general model of cortico-cerebellar connectivity without overfitting the data. To validate that this approach enabled us to distinguish between different forms of cortico-cerebellar connectivity, we conducted a range of model recovery simulations (see Methods for details). Using the measured cortical activity for each participant, we generated artificial sets of cerebellar data imposing either a one-to-one or a many-to-one mapping. Following the procedure used with the real data, we trained the three models on Task Set A and tested them on Task Set B. The simulations (Figure 2—figure supplement 1) showed that the WTA performed best if each cerebellar voxel was connected to only one cortical parcel, whereas ridge regression performed better if there was substantial convergence, with Lasso performing at an intermediate level. Thus, despite the presence of some degree of collinearity between the cortical parcels, these simulations demonstrate that our modeling approach is able to distinguish between different forms of connectivity.

Cortico-cerebellar connectivity is best captured by models with convergence

We first compared the different models, asking how well they predicted activity patterns obtained when the same participants were tested on Task Set B (Figure 2A). Models allowing for some degree of convergence outperformed the WTA model, and this advantage was observed across all levels of cortical granularity. Indeed, the prediction performance for the Ridge, Lasso, and WTA models was relatively independent of granularity. Averaged across all levels of granularity, the Ridge models outperformed the WTA models (F1,23 = 47.122, p<0.001). Post-hoc tests revealed that this advantage for the Ridge model was consistent across all levels of granularity, starting with a parcellation of 80 regions (t23=5.172, p<0.001). We also found a significant difference in predictive accuracy between the Ridge and Lasso model (F1,23 = 15.055, p<0.001). Post-hoc tests showed that there was no significant difference in predictive accuracy at the lowest level of granularity (t23=1.279, p=0.213), whereas the difference was significant for the finer parcellations (all t23 >10.937, p<0.001).

Figure 2 with 4 supplements see all
Performance of cortico-cerebellar connectivity models.

(A) Predictive accuracy (Pearson correlation) of the Ridge, Lasso, and WTA regression models for the test data of Task Set B (B) Predictive accuracy normalized to the noise ceiling based on reliability of both cerebellar and cortical data (see Methods). Source data files are available for panels A and B (‘model_evaluation.csv’) (C) Voxelwise map of inter-session reliability of the test data. (D) Voxelwise map of predictive accuracy of the Ridge model (1848 parcels), normalized to the noise ceiling. (E) Observed and predicted activity for a novel task involving spatial working memory (Spatial Map) and a social cognition task (Theory of Mind). The spatial map task was not included in Task Set A. Error bars indicate 95% confidence intervals across participants (n=24).

Figure 2—source data 1

Source data file contains model evaluation predictive accuracies for each of the three methods (Ridge, Lasso, Winner-Take-All), and for each parcellation (80-1848).

https://cdn.elifesciences.org/articles/81511/elife-81511-fig2-data1-v1.csv

The mean predictive accuracy of the Lasso and Ridge models was 0.257 (Figure 2A). There are two issues of note here. First, the models predicted activity of individual voxels, without any smoothing across cerebellar voxels. Thus, the upper bound for the predictive accuracy is limited by the reliability of the cerebellar test data. The correlation of the measured cerebellar activity patterns across the two sessions of Task Set B was, on average, r=0.51 (SD = 0.102). Reliability was fairly consistent across the cerebellum (Figure 2C), with some decreases in lobules I-IV. This dropoff likely reflects the fact that our battery only included hand and eye movements, and did not activate the lower body representation that is prominent in this area. Second, predictive accuracy is also limited by the reliability of the cortical data that was used to make the prediction. This is an especially limiting factor for WTA models that use fine granularity, since the predictions will be based on data obtained from a small cortical region.

Given these issues, we calculated a noise ceiling for each model, taking into account the reliability of the cerebellar data, the reliability of the cortical data, and the effect of granularity (see Methods, Figure 2—figure supplement 2). As an unbiased comparison of the model predictions, Figure 2B re-plots the predictive accuracy of each model normalized by its noise ceiling. As with the original analysis, the Ridge model significantly outperformed the WTA model (F1,23 = 16.49, p<0.01), and in post-hoc tests, the advantage was especially pronounced for finer parcellations (1848 regions; t23=5.073, p<0.001). Overall, the noise ceiling calculation showed that the Ridge model was able to predict approximately 45% of the systematic variance of the cerebellar activity patterns across tasks (average R2=0.672). While the predictive accuracy (Figure 2D) was best in anterior motor regions, it was reasonably high across the entire cerebellar surface.

The predicted and observed activity patterns for two exemplary tasks (Figure 2E) demonstrate the quality of these predictions at the group level. In both of these examples, (spatial working memory and social cognition tasks), the connectivity model predicted the pattern of task activity with a high degree of fidelity. This is especially compelling for the spatial working memory task as the training set did not include a task with similar characteristics.

To ensure that our results were not biased by the use of an essentially arbitrary parcellation of the neocortex, we repeated the analysis using a range of published functional parcellations. For example, we used a 7-network cortical parcellation based on resting state fMRI data (Yeo et al., 2011) to train and test the three models. Here, too, the WTA model was inferior to the Ridge models (t23=2.956, p<0.01), with no performance difference between Ridge and Lasso models (t23=–1.235, p=0.229). The same pattern held when we used other common functional parcellations of the neocortex (Figure 2—figure supplement 3).

In summary, these results demonstrate that models which entail some degree of convergence from the neocortex to the cerebellum outperform a model in which cerebellar activity is based on input from a single cortical region. This conclusion holds across a broad range of cortical parcellations (Zhi et al., 2022).

Convergence of neocortical inputs varies across the cerebellum

To gain insight into where these cortical inputs came from, we visualized the cortical weights for each of the 10 functional regions of the cerebellum (Figure 3—animation 1). As expected given the crossed connectivity between M1 and lobules IV/V (Kelly and Strick, 2003; Krienen and Buckner, 2009; O’Reilly et al., 2010), the input to the hand regions of the cerebellum (regions 1 and 2) was centered around the contralateral primary sensorimotor cortex with some additional input from premotor and parietal cortex. Regions 3 and 4 are the other two cerebellar regions associated with motor function. Activity in region 3, the oculomotor vermis, is predicted by a bilateral set of cortical regions including the frontal eye fields (FEF), regions in the intraparietal sulcus, and extrastriate visual regions. Region 4, an area strongly activated during action observation, is predicted by a bilateral network of regions including premotor cortex, supplementary motor cortex (SMA) and parietal cortex. Cerebellar regions 5–10, the regions associated with more cognitive processes are predicted by a distributed set of cortical regions, with stronger input coming from the contralateral cortical hemisphere. For example, cerebellar region 5, a region restricted to the left cerebellum, receives much stronger input from the right cerebral hemisphere. In summary, these results suggest that most cerebellar regions are best predicted by multiple cortical regions, pointing to some degree of cortical-cerebellar convergence.

Figure 3 with 3 supplements see all
Cortical connectivity weight maps for the Ridge regression model with 1848 cortical parcels for each of 10 functional cerebellar regions.

Each region is denoted by the most important functional term (King et al., 2019). Results are averaged across participants. Regression weights are in arbitrary units. See Figure 3—animation 1 for a gif of the connectivity weight maps, and see Figure 3—figure supplement 1 for the corresponding analysis using Lasso regression. Figure 3 has been adapted from Figure 5 from King et al., 2019.

Visual inspection of the connectivity patterns (Figure 3) also highlights interregional variation of convergence. For example, inputs to hand motor regions MDTB (regions 1 and 2) arise from a relatively small area of the neocortex, while inputs to regions in lobule VII (regions 3–10) come from a larger area. To quantify this observation, we tallied the number of non-zero regression weights for each cerebellar voxel as a measure of its input surface area. For this calculation, we used the Lasso model (80 cortical parcels) since it uses the minimal set of cortical areas necessary for predicting each cerebellar voxel. As the Lasso model forces the estimates of the other weights to zero, it allows for a quantification of input area without applying an arbitrary threshold. This calculation revealed a substantial variation in the degree of convergence across the cerebellar cortex (Figure 4A). For example, predicting the activity pattern of voxels within Crus I required inclusion of up to 10% of the cortical surface whereas predicting the activity of voxels in the anterior lobe required less than 5% of the neocortex.

Figure 4 with 2 supplements see all
Cortico-cerebellar convergence measures using the Lasso model.

(A) Map of the cerebellum showing percentage of cortical parcels with non-zero weights for the Lasso model (n=80 parcels). (B) Percentage of parcels with non-zero weights for functional subregions of the cerebellum. (C) Spherical dispersion of the connectivity weights on the cortical surface for each cerebellar voxel. (D) Average cortical dispersion for each functional subregion of the cerebellum. Error bars indicate standard error of the mean across participants. See Figure 4—figure supplement 1 for the same results using Ridge regression. Error bars indicate 95% confidence intervals across participants (n=24).

To statistically analyze these data, we opted to bin the cerebellar data using a functional 10-region parcellation of the cerebellum (King et al., 2019). This analysis confirmed that the size of the estimated cortical input area differed across functional regions (F9,207 = 7.244, p<0.001, Figure 4B). The lowest level of convergence was observed for regions 1 and 2, the anterior hand regions of the cerebellum. The highest levels of convergence were found in regions of the cerebellum associated with more cognitive functions (e.g. regions 7 and 8, areas engaged during language comprehension). Noteworthy, this pattern holds for cortical parcellations of different levels of granularity (Figure 4—figure supplement 2), as well as for the thresholded coefficients from the Ridge regression models (Figure 4—figure supplement 1A, B).

To assess whether these differences were driven by the collinearity of the cortical data, rather than by cortical-cerebellar connectivity, we ran a simulation (see Methods, model recovery simulations) in which we replaced the cerebellar data with the activity profile of the most similar cortical parcel. In this simulation, we did not find any differences in the area of estimated input across cerebellar regions (F9,207=1.762, p=0.076). Thus, the observed variation in convergence cannot be explained by the collinearity between different cortical input regions.

As an independent method to quantify convergence, we determined the spatial spread of the inputs across the cortical surface. For example, MDTB region 4, which is activated by action observation, is explained by a set of cortical regions with a relatively small surface area. Nonetheless these regions are spread widely across the cortex (e.g. fusiform gyrus, parietal and premotor regions). For a measure of dispersion, we calculated the spherical variance of the non-zero connectivity weights on an inflated model for each cerebral hemisphere (see Methods). This analysis revealed a similar pattern as seen in the surface area measure (Figure 4C, D F9,207=18.322, p<0.001). For example, the cortical inputs to the hand motor regions of the cerebellum were more concentrated, whereas the cortical inputs to lobule VII were more dispersed. Again, this pattern also holds for the Ridge models (Figure 4—figure supplement 1C, D).

In summary, we observed variation in the degree of cortico-cerebellar convergence across the cerebellar cortex. In particular, cerebellar motor areas such as the hand region in lobule V and VIII received relatively focal cortical input whereas cerebellar areas associated with cognitive processes (e.g. language, working memory) in lobule VII exhibited a higher degree of convergence.

Cortico-cerebellar connectivity model predicts new cerebellar data

A strong test of a model is how well it predicts novel data obtained in contexts distinct from that data used to develop the model. We conducted a generalization test using data from a separate experiment involving novel tasks and new participants. The data came from a study in which participants were trained over multiple sessions on five tasks, selected to span a set of distinct cognitive domains.

All three model types (WTA, Lasso, Ridge) were evaluated, using the connectivity weights estimated from Task Set A of the main study. Because the new study had different participants, we averaged the weights for each model across the individuals in the training set. These group-averaged weights were then used to predict cerebellar activity patterns for the new data.

As shown in Figure 5, the connectivity models were generally successful in predicting cerebellar activity patterns for the new dataset. The overall pattern is quite similar to that seen in the initial model tests (in which we had used data from different tasks) but involved training and test data from the same participants. Predictive accuracy was stable across levels of cortical granularity and the Ridge model provided the best overall predictive accuracy (r=0.657) and the WTA (r=0.352) the worst predictive accuracy. These results provide evidence that our cortico-cerebellar connectivity models capture important aspects of connectivity that are stable both across tasks and participants. Moreover, in accord with the earlier analyses based on predictions at the individual level, this generalization analysis again suggests that approximately 43% of the variation of cerebellar activity across tasks can be predicted by cortical activity alone.

Generalization to new dataset.

Models of cortico-cerebellar connectivity are tested in a new experiment. Each model is tested across different levels of cortical granularity. Predictive accuracy is the Pearson correlation between observed and predicted activity patterns, normalized to the noise ceiling. Error bars indicate 95% confidence intervals across participants (n=20).

Discussion

To date, models of connectivity between the human neocortex and cerebellum have been based on fMRI resting-state data (Buckner et al., 2011; Ji et al., 2019; Marek et al., 2018). This work demonstrates that each region of the cerebellum receives input from a distinct set of cortical regions. For example, anterior and posterior cerebellar motor regions show correlated activity with contral-lateral sensorimotor cortex and non-motor or ‘cognitive’ cerebellar regions show correlated activity with specific parietal and frontal association areas.

Despite these important insights, previous work has not been designed to examine the patterns of convergence between the neocortex and cerebellum. Resting state connectivity maps are generally produced by assigning each cerebellar voxel to a single cortical network, a de facto winner-take-all model. In the present study, we quantified and compared models of cortico-cerebellar connectivity. The models differed in the degree of convergence of cortical inputs to each cerebellar region, ranging from an architecture constrained to a strict one-to-one mapping to architectures that allowed for distributed inputs. We evaluated these models in terms of how well they could predict cerebellar activity patterns on a novel Task Set. For nearly the entire cerebellum, models that allowed for some convergence predicted cerebellar activity better than the WTA model.

Convergence differs across cerebellar circuits

Importantly, the amount of convergence differed across the cerebellar cortex. Specifically, regions in anterior (lobules I-V) and inferior cerebellum (lobules VIII-X) were well predicted by a relatively small and concentrated area of the neocortex. In contrast, regions in lobule VI and especially lobule VII required input from a larger and spatially more dispersed set of cortical regions that were primarily located in association areas of prefrontal, parietal, and cingulate cortex. This finding underscores the heterogeneity of cortico-cerebellar connectivity with some cerebellar areas functioning in nearly a 1:1 relationship with a single cortical region, whereas other areas integrate input from a more diverse set of cortical regions.

This variation bears some resemblance to a motor/cognitive gradient identified from resting state data (Guell et al., 2018). However, there are a few notable exceptions. First, based on our evaluation metrics, Region 3 (oculomotor vermis) was best explained by a large and relatively dispersed set of cortical regions, including intraparietal sulcus, the frontal eye fields (FEF), and extrastriate visual areas (Figure 4). Second, sub-regions of lobules IX, functionally associated with non-motor tasks (Buckner et al., 2011), are best explained by a relatively restricted set of cortical regions (Figure 4A and C). Thus, rather than following a functional distinction, it appears that regions with strong convergence are anatomically restricted to lobules VI and VII, two areas that have disproportionately increased in size during primate evolution (Balsters et al., 2010).

Variation in cortico-cerebellar connectivity has important implications for theories of cerebellar function (De Zeeuw et al., 2021; Diedrichsen et al., 2019; Henschke and Pakan, 2020). A cerebellar region that forms a closed-loop circuit with a single cortical region can only fine-tune the computational processes that are ongoing in its connected cortical region. In contrast, a cerebellar region that receives convergent input would be able to coordinate the interactions between dispersed cortical regions. Indeed, using a rodent model, Pisano et al., 2021 have identified cerebellar areas that project back to multiple cortical areas, providing a substrate by which the cerebellum can influence neuronal dynamics not only in a focal cortical area but also interactions within networks of cortical areas (Pisano et al., 2021).

Task-based vs. resting-state connectivity analyses

While our modeling approach could be applied to resting-state data, we opted to use the data from the multi-domain task battery for a variety of reasons. First, the breadth of tasks included in the battery assured that there would be reliable variation in activity across most of the neocortex and cerebellum, a necessary precondition for building a complete connectivity model. Second, the dataset allowed us to avoid biases that arise from correlated noise across abutting cortical and cerebellar regions (Buckner et al., 2011). We used the cortical activity patterns from one session to predict cerebellar activity patterns in a different session, relying on the fact that measurement noise is independent across sessions. Third, the MDTB data set allowed us to test the model across a broad set of tasks and mental states. While resting-state correlations are predictive of task-based activation patterns (King et al., 2019; Tavor et al., 2016; Zhi et al., 2022), resting in an fMRI scanner is arguably a relatively restricted situation. Our generalization test shows that the connectivity model can robustly predict activity patterns for new tasks and participants.

Directionality of the model

We opted in the present study to model the correlations between neocortical and cerebellar activity in a directional manner, making inferences from these data about connectivity from the neocortex to the cerebellum. We recognize that the two structures communicate with each other in a closed loop. Our decision to focus on cortico-cerebellar connectivity is based on studies of cerebellar blood flow in rodents. Increases in cerebellar blood flow, a major contributor to the BOLD signal, are the result of increases in mossy fiber activity, the primary input to the cerebellar cortex (Alahmadi et al., 2015; Alahmadi et al., 2016; Gagliano et al., 2022; Mapelli et al., 2017; Mathiesen et al., 2000; Thomsen et al., 2004). In contrast, even a dramatic increase in the (complex or simple spike) firing rate of Purkinje cells does not produce a measurable change in blood flow within the cerebellar cortex (Thomsen et al., 2004; Thomsen et al., 2009). Thus, the cerebellar BOLD signal provides a relatively clear image of the cortical inputs to the cerebellum; importantly, this signal does not provide information about the output of the cerebellar cortex (or cerebellum in general). In contrast, the BOLD signal in the neocortex reflects a combination of signals including local inputs/outputs across cortical layers, inputs from other cortical regions, and ascending inputs such as from the thalamus which will include input from the cerebellum. Given this asymmetry, a priori, we should expect that correlations in the fMRI signal between the neocortex and cerebellum will more strongly reflect the flow of information from the cortex to the cerebellum than the reverse.

Nonetheless, it is possible that some of the convergent inputs identified in our modeling work may be due to divergent projections from a single cerebellar region to multiple cortical regions. Indeed, recent viral tracing work has shown that some cerebellar areas project to the reticular nuclei in the thalamus, which in turn projects to a wide array of cortical regions (Kelly and Strick, 2003; Pisano et al., 2021). A complete analysis of the entire cortico-cerebellar circuit will require the methods such as viral tracing techniques (Kelly and Strick, 2003; Pisano et al., 2021) or functional activation measures that target key nodes in the ascending pathway (deep cerebellar nuclei, thalamus).

Methodological limitations

Inter-region correlations of fMRI data can of course only provide indirect evidence of true functional connectivity. As such, it is important to consider methodological limitations that may influence the validity of our conclusions. From a statistical perspective, it was not clear, a priori, that we would have the power to distinguish between models of connectivity given that there can be substantial collinearity between different cortical regions. The model-recovery simulations (Figure 2—figure supplement 1A, B) suggest that the present dataset was suitable to make such inferences, namely, activity patterns in different cortical regions were sufficiently de-correlated, likely reflecting the use of a broad task battery. Thus, we were able to recover the correct model used to simulate the data, regardless of whether we assumed one-to-one connectivity or different degrees of convergence.

However, the simulations also indicated that the approach (and data) was not sufficient to determine the absolute degree of convergence with high confidence. For example, the size of the cortical input area for each cerebellar region differed substantially between the Ridge and Lasso regression models (Figure 3 vs. Figure 3—figure supplement 1). Nonetheless, the two models result in similar predictive accuracy when using real data. Therefore, the true extent of the cortical input to a cerebellar region likely lies somewhere between the extremes provided by the Ridge and Lasso models. Importantly, this ambiguity does not impact the core observation that the degree of convergence varies systematically across the cerebellar cortex. Whether using measures based on cortical surface area or dispersion, the general picture of variation in connectivity holds for both the Ridge and Lasso models, as well as when using cortical parcellations of different granularity. Thus, we are confident that the variation in convergence reflects a stable, method-independent characteristic of the cortico-cerebellar system.

Future directions

Our approach was designed to identify the best task-general model of cortico-cerebellar connectivity. The connectivity model was able to take cortical activity patterns to make fairly accurate predictions of cerebellar activity, even when the data were obtained in a completely separate experiment from that used to build the model. This finding suggests that a substantial proportion of the cortico-cerebellar communication can be described by fixed information flow (Cole et al., 2016); that is, each area of the cerebellum simply receives, independent of task, a copy of the neural activity occurring in the corresponding cortical areas.

However, the task-general model did not provide a perfect prediction of cerebellar activity. While these predictions may improve with more training data, we believe that there will likely be some systematic failures of the model. Such task-dependent deviations from the model prediction may offer important insights into cerebellar function, indicating that the cerebellum is more active than predicted by the cortical activity for some tasks and less for others. This pattern would suggest a ‘gating’ of cortical input to the cerebellum, perhaps within the pontine nuclei, which is the main relay station of cortical inputs to the cerebellum. Rather than just transmitting cortical input to the cerebellar cortex, these nuclei may serve as an adaptive gate (Schwarz and Thier, 1999), amplifying or attenuating information as a function of the relative importance of cerebellar processing for the current task.

The models presented here allow us to detect such deviations by comparing the observed cerebellar activity for any task against the activity that is predicted by a task-invariant connectivity model. As such, the model provides a potent new tool to test hypotheses concerning cerebellar function.

Methods

Multi-domain task battery

To build models of cortico-cerebellar connectivity, we used the publicly available multi-domain task battery dataset (MDTB; King et al., 2019). The MDTB includes fMRI data from two independent Task Sets (A and B). Each set consists of 17 tasks, eight of which were common to both sets. The 26 unique tasks were designed to sample activity during a broad range of task domains including motor (e.g. sequence production), working memory (e.g. 2-back task), language (e.g. reading), social (e.g. theory of mind), cognitive control (no-go, Stroop), and emotion (e.g. facial expression; see Table S1 in King et al., 2019).

Twenty-four participants (16 females, 8 males, mean age = 23.8) were scanned during four sessions, with Task Set A employed in sessions 1 and 2 and Task Set B in sessions 3 and 4. Within a session, there were eight 10-min runs, with each run composed of 17 blocks of 35 s each. Each block involved a unique task and included an initial 5 s instruction period. Most tasks contained multiple conditions (i.e. different levels of difficulty), resulting in a total of 47 unique task conditions.

Image acquisition and preprocessing

All fMRI and MRI data were collected on a 3T Siemens Prisma located at the Center for Functional and Metabolic Mapping at Western University, Canada. The protocol used the following parameters: 1 s repetition time; field-of-view measuring 20.8 cm; P-to-A phase encoding direction; 48 slices; 3 mm thickness; in-plane resolution 2.5×2.5 mm2. In order to localize and normalize the functional data, a high-resolution anatomical scan (T1-weighted MPRAGE, 1 mm isotropic resolution) of the whole brain was acquired during the first session.

Functional data were realigned for head motion artifacts within each session and different head positions across sessions using a six-parameter rigid body transformation. The mean functional image was co-registered to the anatomical image and this rigid-body transformation was applied to all functional images. No smoothing or group normalization was applied. SPM12 ( spm/doc/spm12_manual.pdf) and custom-written scripts in MATLAB were used to conduct data pre-processing.

General linear model (GLM)

To generate estimates of the activity related to each task condition, a general linear model (GLM) was fitted to the time series data of each voxel (ti). This was done separately for each imaging run and Task Set (A and B). The Design matrix for the GLM (Z1) of the GLM consisted of the regressors for the different conditions, plus separate regressors for each of the 5-s task-specific instructions period. For each Task Set, the beta weight for each task condition was averaged across the 8 runs within a session, resulting in the average activity estimate for that session (b¯j) for each voxel. The regressor for Task Set A resulted in 46 activity estimates (17 Instructions + 29 condition regressors) per session and Task Set B resulted in 49 activity estimates (17 instructions + 32 conditions regressors) per session.

Neocortex surface reconstruction

For each of the 24 participants, the anatomical surfaces of the cortical hemispheres were reconstructed using the standard recon-all pipeline of the FreeSurfer package (Fischl, 2012; v. 5.0). The pipeline included brain extraction, generation of white and pial surfaces, inflation, and spherical alignment to the symmetric fsLR-32k template (Van Essen et al., 2012). Individual surfaces were re-sampled to this standard grid, resulting in surfaces with 32,492 vertices per hemisphere.

Spatial normalization of cerebellar data

The cerebellum was isolated and normalized to the high-resolution Spatially Unbiased Infratentorial Template of the cerebellum using the SUIT toolbox (Diedrichsen, 2006). This non-linear transformation was applied to both the anatomical and functional data. Task condition activity estimates (i.e. the beta weights) were resampled to a resolution of 3 mm isotropic and resliced into SUIT space. The cerebellar isolation mask was edited to remove voxels in the superior cerebellum that abutted voxels in the primary visual cortex. Functional images were masked with the cerebellar isolation mask resulting in activation signals that originate only from the cerebellar cortex. The cerebellar data were visualized using a surface-based representation of the cerebellar gray matter included in the SUIT toolbox (Diedrichsen and Zotow, 2015). This flat map is not intended to represent a true unfolding of the cerebellar cortex, but rather provides a convenient way to visualize volume-averaged cerebellar imaging data in a 2d representation.

Connectivity models

All of the task-based connectivity models were, in essence, multiple-regression models in which the data for each of cerebellar voxels (yi) was modeled as a linear combination of Q cortical parcels (X, see next section).

(1) yi=Xwi+ei

By combining the data (Y), the connectivity weights (W), and the measurement error (E) for each cerebellar voxels as columns into a single matrix, the entire model could be written as:

(2) Y=XW+E

The connectivity approach that we used relied only on task-evoked activity: model training and testing were done on the fitted time series rather than the full or residual time series (as done in many other connectivity approaches). We multiplied Zj from the first-level GLM with the 46 (for training) or 49 (for testing) activity estimates (b¯j) from the first-level GLM with the first-level design matrix (Z) to obtain the fitted time series the first-level design matrix for each session (Z¯j) from the first-level GLM with tfittedi=Zb¯+E for each voxel (see General Linear Model). This was done for both the cortical and the cerebellar voxels.

In calculating the variances and covariances between cortical and cerebellar data (which are computed in the multiple regression models), this procedure reweighted the activity estimates to account for the fact that estimates for the instructions were based on 5 s of fMRI data, while estimates for the conditions were based on 10–30 s of data. For example, the covariance between the cortical and cerebellar time series is:

(3) YTX=b¯XTZTZb¯Y

Therefore, if we had used the task-evoked activity estimates (b¯) for cerebellar and cortical data, but reweighted the influence of each regressor by the diagonal of, ZTZ we would have obtained similar results. Because of the off-diagonal terms, this method also mitigates problems that may arise in event-related designs due to correlations between regressors (Mumford et al., 2012), one example is the estimation covariance between the instruction period and the task-related regressor that follows immediately afterwards.

We normalized the fitted time series by dividing them by the standard deviation of the residuals from the first-level GLM. This procedure emphasized voxels with large signal-to-noise ratios over voxels with low signal-to-noise ratios. The normalized time series were then averaged within all Q cortical parcels, resulting in a TxQ matrix (X) of cortical features. For the cerebellum, we modeled each of the p=6937 cerebellar voxels separately (SUIT atlas space with 3 mm resolution), resulting in a TxP data matrix (Y). The cortical data were further z-standardized (as in a partial correlation analysis) before entering them into regression models (see Model Estimation).

When estimating connectivity models, correlated fMRI noise across the cortex and cerebellum (e.g. head movement, physiological artifacts) can lead to the erroneous detection of connectivity between structures. This is especially problematic for the superior cerebellum and the directly abutting regions of the occipital and inferior temporal lobes. To negate the influence of any noise process that was uncorrelated with the tasks, we used a ‘crossed’ approach to train the models. Because there were two sessions (same tasks, different order) in each task set (A and B), we were able to predict the cerebellar fitted time series for the first session by the cortical time series from the second session, and vice-versa (see Figure 1).

(4) Z1b1XZ1b2Y

We did this separately for each task set, and given that the sequence of tasks was randomized across sessions, we could conclude from this approach that this prediction was based on task-related signal changes rather than fluctuations attributable to noise processes given that noise is not shared across scanning sessions.

Cortical parcels

We created a set of models that used different levels of granularity to subdivide the cortex. Specifically, each hemisphere was divided into a set of regular hexagonal parcels, using icosahedrons with 42, 162, 362, 642, or 1002 parcels per hemisphere [see (Zhi et al., 2022) for details]. Parcels that were completely contained within the medial wall were excluded, resulting in 80, 304, 670, 1190, or 1848 parcels combined for the left and right hemisphere. Activity within a parcel was averaged across voxels for each condition, defining the Q regressors for each model.

The regular icosahedron parcellations are arbitrary in the sense that they do not align to functional boundaries in the human neocortex more than expected by chance (Zhi et al., 2022). We therefore also employed functionally-defined parcellations, repeating the main analyses using 12 different parcellations derived from resting state fMRI data (see Figure 2—figure supplement 3).

Model estimation

We used three regression methods to estimate the connectivity matrix W^ at the individual participant level. In the core analyses, each of these methods was combined with the five arbitrary cortical parcellations, resulting in 15 connectivity models. Each method was selected to favor a specific pattern of cortico-cerebellar connectivity. For the winner-take-all models, we assumed that the time series of each cerebellar voxel is explained by one and only one cortical region. To build this model, we simply chose the cortical region with the highest correlation with the cerebellar voxel in question. The connectivity weight corresponding to the rest of the cortical regions were set to 0.

The other two methods allowed for some degree of convergence in the cortical input to each cerebellar voxel. Lasso regression [least absolute shrinkage and selection operator (Tibshirani, 1996)], seeks to explain activity in each cerebellar voxel by a restricted set of cortical features. Specifically, Lasso minimizes the squared-error, plus an additional penalty calculated as the sum of the absolute values (L1-norm) of the regression coefficients:

(5) W^=argminwYXW22+λW1

In contrast, Ridge regression penalizes large connectivity weights, attempting to find a broad set of cortical regions to explain each cerebellar voxel. Ridge regression minimizes the following loss function:

(6) W^=argminwYXW22+λW22

Both Lasso and Ridge models include a hyperparameter λ that must be tuned in model estimation. The value of these regularization parameters was determined using gridsearch with ourfold cross-validation on the training set. We divided the conditions of the training set into four non-overlapping subsets, reconstructed the time series using tasks from three of the four subsets and evaluated the model using the left-out subset. Figure 2—figure supplement 4 depicts the average predictive accuracy for the training set and the optimal value of λ for the Lasso and Ridge models.

Model testing

After model training with Task Set A, Task Set B was used as the test set. We applied the same procedure to construct the X and Y matrices. We then used the estimated weights for each model to generate predicted voxel-wise cerebellar activity patterns from the observed cortical time series. These predictions were generated at each level of granularity for each participant. Model performance was measured by correlating predicted and observed cerebellar time series for each voxel. As with the procedure for model training, the cortical and cerebellar time series were crossed: Cortical time series from session 3 were used to predict cerebellar time series from session 4 and vice versa. Model predictions were calculated separately for each cerebellar voxel and visualized on voxelwise maps of the cerebellar cortex.

In sum, we evaluated 15 models in which the cortical parcels were based on an arbitrary icosahedron, based on the combination of three regression methods and five levels of granularity. For the 12 functionally-defined parcellations (Figure 2—figure supplement 3), we limited the evaluation to the WTA and Ridge models.

Noise ceiling

For each model evaluated in the current study, the noise ceiling quantifies the expected performance of each model, under the assumption that the estimated weights reflect the true connectivity between the neocortex and cerebellum. When predicting the cerebellar activity patterns, our models are limited by two factors: Measurement noise associated with the cerebellar data and measurement noise associated with the cortical data.

To estimate these sources of noise, we used the fitted cerebellar time-series for the two sessions, assuming that each session’s data are composed of the true time series Y* plus noise ε: Yi=Y+ϵ. If XiW signifies the variance of the true time series, and σϵ2 the variance of the noise, then the correlation between the two sessions, i.e. the split-half reliability, is:

(7) r(Y1,Y2)=σY2σY2+σϵ2.

Similarly, the reliability of the prediction, even if we knew the true connectivity weights W, is limited by the noise on the cortical data X. The prediction XiW has noise variance, σXW2 which can be estimated by calculating the split-half reliability across the two sessions:

(8) r(X1W,X2W)=σY2σY2+σXW2

Thus, for the true model the expected correlation would be:

(9) rceil=r(XiW,Yi)=σY2σY2+σXW2σY2+σϵ2=r(Y1,Y2)r(X1W,X2W)

Because we do not know the true weights, we use the estimated weights from each model to estimate the noise ceiling. Thus, the noise ceiling is model dependent, specifying what the correlation would be if the current model was the true model.

Cortico-cerebellar convergence

Two measures were used to assess the amount of convergence of cortical inputs to each cerebellar voxel. For the first measure, we calculated the percentage of the cortical surface contributing to the prediction in each cerebellar voxel. We used the estimated weights from the Lasso regression model with the best unnormalized performance (80 cortical parcels). We opted to focus on the Lasso model since it sets weights that make a negligible contribution to the predictions to zero, thus identifying the sparsest cortical input that explains the cerebellar data. To determine the likely input area, we determined the number of non-zero coefficients for each cerebellar voxel, and expressed this as a percentage of the cortical surface. We also repeated these analyses using Ridge regression (Figure 3—figure supplement 1). Because Ridge does not result in zero coefficients, we applied a threshold, counting the number of connectivity weights with a value one standard above the mean of the weight matrix.

For the second measure, we calculated the dispersion (spherical variance) of the input weights. Using the spherical representation of the cerebral hemispheres from FreeSurfer, we defined unit vectors, vi , pointing from the center of the sphere to the center of each of the cortical parcels. Each vector was then weighted by the cortico-cerebellar connectivity weight for that region wi. All negative connectivity weights were set to zero. The spherical variance for hemisphere h across all Qh parcels of this hemisphere was defined as:

(10) varh=1|1QhQhwivi|

To obtain a composite measure, we averaged the variances for the two hemispheres, with each weighted by the size of the summed (non-negative) connectivity weights. For summary and statistical testing, we then averaged the size of the input area, as well as the dispersion, across cerebellar voxels within each of the 10 functionally defined MDTB regions.

Model recovery simulations

Even though we are using data from a broad task battery, there is substantial collinearity between the cortical parcels. This raises the questions of whether our approach can reliably distinguish between models of sparse and broad connectivity, and whether we can detect differences in convergence across regions in an unbiased fashion.

To validate our model-selection approach we conducted a series of model recovery simulations. Using the observed cortical activity patterns from each individual participant (with 80–1848 parcels), we simulated artificial cerebellar data for 2000 voxels. The artificial connectivity weights were chosen under two scenarios. To simulate one-to-one connectivity, each cerebellar voxel was connected to only one randomly chosen cortical parcel. To simulate broad convergence, the connectivity weights were drawn from a normal distribution with mean = 0, and a variance of 0.2. We then generated artificial data following Equation 2, both for the training data (using the cortical data from Task Set A) and for the testing data (using the cortical data from Task Set B). The measurement noise (E) was drawn from a normal distribution with mean = 0. The variance (0.25) was adjusted to approximate the signal-to-noise level observed in the empirical data. Model estimation (including hyperparameter-tuning) and model testing were then conducted following the procedure used for the empirical data (see above).

We also tested whether differences in cortical input area (Figure 4b) could be due to differences in collinearity in the associated cortical parcels. To address this question, we conducted a simulation with 80 cortical parcels. For each MDTB region and participant, we simulated cerebellar data using one-to-one connectivity with the cortical parcel that had the highest average lasso-regression coefficient for that region and participant. Thus, we replaced the average activity profile in each cerebellar region by the profile of the most similar cortical parcel. We submitted these data to a Lasso regression with log-lambda set to –3 and calculated the percentage of the cortical surface with non-zero connectivity weights in the same manner as done with the real data.

Generalization to new participants

Although the models were trained and tested on distinct datasets with more than half the tasks unique to each set, the data for training and testing were obtained from the same set of participants. As a stronger test of generalization, we evaluated the models with data from a new experiment involving naive participants (n=20, 11 females, 9 males, mean age = 21.3). These participants were trained to perform a 5-task battery (cognitive control, social prediction, action prediction, visual search, semantic prediction). For each task, there was an easy and hard condition (e.g., 2-back vs 0-back for cognitive control, high vs low predictability), along with the task’s associated instruction. Thus, there were a total of 16 conditions when including rest (fixation). Each participant completed five sessions, with fMRI data obtained in the first, third, and fifth sessions. Within a scanning session, there were six 11-min runs, and each run included three repetitions of each task (35 s) with a 10 s rest period prior to each 5 s instruction period.

The participants were scanned on a Siemens MAGNETOM TrioTim syngo MR B17 located in the Henry Wheeler Jr. Brain Imaging Center at the University of California, Berkeley. The protocol used the following parameters: 1 s repetition time; field-of-view measuring 20.8 cm; A-to-P phase encoding direction; 48 slices; 3 mm thickness; in-plane resolution 2.5×2.5 mm2. To localize and normalize the functional data, a high-resolution, whole-brain anatomical scan (T1-weighted MPRAGE, 1 mm isotropic resolution) was acquired during the first session. fMRIPrep (Esteban et al., 2019; Esteban et al., 2020; https://fmriprep.org/en/stable/) was used to preprocess the anatomical and functional data, following the same analysis procedures as conducted for the main experiment.

To generate estimates of activity (beta weights), a General Linear Model (GLM) was fitted for each voxel to the time series data. Each of the 16 conditions was modeled as a separate regressor, with this repeated separately for each imaging run. Nipype (https://nipype.readthedocs.io/en/latest/) and custom-written scripts in Python were used to estimate the beta weights from a first-level GLM.

To evaluate model generalization to this new dataset, we used the full set of cortico-cerebellar connectivity models estimated with the original data set, limiting the analysis to the regular cortical parcellations. Since the participants were different from the original data set, we averaged the model weights for the original 24 participants (i.e. 3 model types x 5 levels of granularity). Each model was then used to predict cerebellar activity from the cortical activity obtained in the new study, with these predictions compared to actual cerebellar activity at the individual level. Given that the experiment was a training study, we expected that both the true activation patterns, as well as the signal-to-noise ratio would change across sessions. We therefore determined the noise ceiling for each session separately by calculating the reliability across odd and even runs within each session. The normalized predictive accuracy was determined separately for each of the three imaging sessions, with the results across sessions averaged.

Source data files

Model training (‘model_training_ridge.csv’ and ‘model_training_lasso.csv’) and evaluation (‘model_evaluation.csv’) summary results are included as additional source data in the manuscript.

Data availability

The multi-domain task battery (MDTB) published in King et al., 2019 was analyzed in the current study. The raw behavioral and whole-brain imaging data generated by this task battery are available on the openneuro data sharing repository (https://openneuro.org/datasets/ds002105/share). The cortical parcellations used in the current study are available to download from https://github.com/DiedrichsenLab/fs_LR_32/tree/master (copy archived at Diedrichsen et al., 2022) and the cerebellar parcellations and contrast maps are available to download from https://www.diedrichsenlab.org/imaging/mdtb.htm. The cerebellar maps are also available for interactive visualization in an atlas viewer (https://www.diedrichsenlab.org/imaging/AtlasViewer/). The experimental code is publicly available at https://github.com/maedbhk/cerebellum_connectivity (copy archived at King, 2023). Model training and evaluation summary results are included as additional source data in the manuscript.

The following previously published data sets were used

References

Decision letter

  1. Michael W Cole
    Reviewing Editor; Rutgers, The State University of New Jersey, United States
  2. Chris I Baker
    Senior Editor; National Institute of Mental Health, United States

Our editorial process produces two outputs: (i) public reviews designed to be posted alongside the preprint for the benefit of readers; (ii) feedback on the manuscript for the authors, including requests for revisions, shown below. We also include an acceptance summary that explains what the editors found interesting or important about the work.

Decision letter after peer review:

Thank you for submitting your article "A task-general connectivity model reveals variation in convergence of cortical inputs to functional regions of the cerebellum" for consideration by eLife. Your article has been reviewed by 3 peer reviewers, one of whom is a member of our Board of Reviewing Editors, and the evaluation has been overseen by Chris Baker as the Senior Editor. The reviewers have opted to remain anonymous.

The reviewers have discussed their reviews with one another, and the Reviewing Editor has drafted this to help you prepare a revised submission.

Essential revisions:

1) Claims of directionality are overstated and should be revised accordingly.

2) The authors should address the differences in model complexity, ideally through an additional analysis (e.g., involving AIC or BIC).

3) More details on the connectivity method used, especially how it handles the task-evoked variance inflation problem.

4) The authors should discuss the possibility that motor tasks are inherently simpler than non-motor tasks, which may have biased their results regarding motor vs. non-motor regions.

Reviewer #2 (Recommendations for the authors):

The authors should compare model performance while accounting for model complexity. One solution is for the Lasso/Ridge models to have the "WTA" weights to be fixed while randomly permuting all other weights. If the Ridge model is truly better, you would expect this random model to perform similarly to the WTA model but worse than the original Ridge/Lasso models. This would suggest that more inputs do not significantly explain more variances in the cerebellum just by chance.

The authors should expand the discussion on the closed-loop anatomical relationship between the cerebellum and the cortex and how its multi-synaptic relationship can allow converging inputs from the cortex to the cerebellum.

For figures 3 and 4, it would be nice to expand and show physically where and how widespread those cortical inputs are. That can help readers contextualize the physical level of "convergence" of the cortex onto cerebellar regions.

Instead of interpreting the results as cortico-cerebellar inputs, discuss how they can be bidirectional interactions.

Reviewer #3 (Recommendations for the authors):

– The authors note that a "convergent architecture would suggest that subregions in the cerebellum integrate information across disparate cortical regions and may coordinate their interactions." It would be interesting to evaluate the degree to which the "disparate" regions that may provide convergent inputs to the cerebellum, as identified in this study, are functionally/anatomically related.

– The authors discuss findings from Pisano et al. 2021 as anatomical evidence for cerebellar outputs from specific regions reaching a variety of cortical targets. However, it is important to note that in most species, and especially in the mouse, it is quite unlikely that any one anatomical tracer injection into the cerebellar cortex can target a functionally distinct region of the cerebellum. See, for example, Fujita et al. 2010 J Comp Neurol where small injections into marmoset cerebellar cortex are difficult to localize to specific longitudinal bands in the cerebellar cortex.

– A map of the cerebellum showing the MDTB regions would be useful to include as a guideline in Figure 3

[Editors' note: further revisions were suggested prior to acceptance, as described below.]

Thank you for resubmitting your work entitled "A task-general connectivity model reveals variation in convergence of cortical inputs to functional regions of the cerebellum" for further consideration by eLife. Your revised article has been evaluated by Chris Baker (Senior Editor) and a Reviewing Editor.

The manuscript has been improved but there are some remaining issues that need to be addressed, as outlined below:

Reviewer #1 (Recommendations for the authors):

The authors did an excellent job improving their manuscript based on previous reviews, except for one major point described below.

R1 point 3 was responded to, but not meaningfully addressed. The authors seemed to misunderstand the issue that regressing out mean task-evoked activity is meant to deal with, resulting in the authors effectively ignoring the issue with regard to their own study. Indicative of misunderstanding the issue, the authors state: "Returning to the papers cited by the reviewer, these are designed to look at connectivity not related to the task-evoked activity." This is a common misunderstanding in the field: Rather than looking at connectivity not related to task-evoked activity, regressing out mean task-evoked activity is meant to deal with a causal confound (common inputs from stimuli), while still allowing estimation of connectivity related to the task-evoked activity. This interest in the impact of the task-related activity on connectivity is clear from the simple logic of most studies utilizing mean task-evoked activity regression, wherein connectivity estimated in the context of condition-specific task-evoked activity is compared. (If all task-evoked activity were removed then there would be no point in comparing task conditions.) What is left over for estimating functional connectivity estimation is trial-by-trial (and moment-to-moment) variation in task-evoked activity (along with spontaneous activity). So, the authors have not dealt with the possibility that common stimulus inputs – such as centrally presented visual stimuli simultaneously activating left and right hemispheres, generating false inter-hemispheric functional connections – are driving some of their connectivity estimates. Ideally, this confound would be dealt with, but at the very least the authors should acknowledge this limitation of their methods.

Beyond this issue (and perhaps more importantly), it is difficult to reconcile the present study's approach with more standard task-related functional connectivity approaches. How can connectivity values estimated via the present study's approach be interpreted? How do these interpretations differ from those of more common task-related functional connectivity approaches? Is there enough advantage to using the present study's unconventional task-related functional connectivity method to make it worth the extra effort for readers to make sense of it and learn to relate it to more standard approaches? The presence of all of these unanswered questions makes it difficult to properly evaluate the manuscript and its results.

Further, in my attempt to answer the above questions myself, I have read and re-read the description of the connectivity approach, and yet I have not been able to make sense of it. I doubt this is my own lack of knowledge or intelligence, but even if that were the case it is likely that many (perhaps most) readers would have the same issue. Therefore, I recommend rewriting the "Connectivity Models" section (within "Methods") to include more details. Some specific details that would be helpful to include:

– The meaning of "fitted time series" is ambiguous. Does this mean the original time series multiplied by the fit β estimate for each condition, resulting in the original time series modulated by the fit to the task regressor for each condition? Or does this mean the regressor time series multiplied by the fit β estimate for each condition? If it is the first option I can imagine interpreting the resulting functional connectivity estimates in a fairly similar way to standard task-related functional connectivity. However, if it is the second option, then I think the interpretation would be quite different from standard task-related functional connectivity. In such a case I think much more details on how to interpret the resulting connectivity values would be warranted. Of course, it could be a third option I haven't considered, in which case it would be helpful to have that option detailed.

– It is unclear what this detail means (pg. 24): "This procedure reweighted the activity estimates in an optimal manner, accounting for the fact that estimates for the instructions were based on 5 s of fMRI data, while estimates for the conditions were based on 10-30 s of data." In what sense is this accounting for the duration of task conditions? This makes me wonder whether the connectivity estimates were based on using cortical β estimates (rather than time series) to predict cerebellar β estimates (rather than time series), which could account for condition duration quite well. But the description is ambiguous with respect to these kind of details, and this possibility doesn't seem to fit well with earlier descriptions of the method.

– It is unclear what this detail means (pg. 24): "This method also mitigates problems that may arise in event-related designs due to correlations between regressors (Mumford et al., 2012)." How exactly? Also, what exactly are the problems that may arise in event-related designs due to correlations between regressors? It sounds like this is related to collinearity, but it is unclear which aspects of the method help with collinearity.

– The details regarding the "crossed" approach need to be fleshed out. It is stated (pg. 25): "we used a "crossed" approach to train the models: The cerebellar fitted time series for the first session was predicted by the cortical time series from the second session, and vice-versa". This procedure is also presented in an ambiguous manner in Figure 1. Right now it is stated in the manuscript that the session 2 cortical time series are used to predict the session 1 cerebellar time series, and yet it is unclear how that would be possible when there are different task conditions across the sessions. Is the connectivity model built using both cortical and cerebellar time series during session 2 (betas from cortex predicting cerebellum), then both cortical and cerebellar time series are used from session 1 to test the session 2 connectivity model (session 2 betas multiplied by session 1 cortical time series to predict session 1 cerebellar time series)? Put more simply: is it the case that the session 2 connectivity model is tested using session 1 data, and that the session 1 connectivity model is tested using session 2 data? If so, the description should change to better describe this. In any case, the description should be clarified with more details.

https://doi.org/10.7554/eLife.81511.sa1

Author response

Essential revisions:

1) Claims of directionality are overstated and should be revised accordingly.

We agree that the correlation structure of BOLD signals in the neocortex and cerebellum is shaped by the closed-loop (bi-directional) interactions between the two structures. As such, some of the observed convergence could be caused by divergence of cerebellar output. We have added a new section to the discussion on the directionality of the model (Page 18).

That said, there are strong reasons to believe that our results are mainly determined by how the neocortex sends signals to the cerebellum, and not vice versa. An increasing body of physiological studies (and this includes newer papers, see response to reviewer #1, comment #1 for details) show that cerebellar blood flow is determined by signal transmission from mossy fibers to granule cells and parallel fibers, followed by Nitric oxide signaling from molecular layer interneurons. Importantly, it is clear that Purkinje cells, the only output cell of the cerebellar cortex, are not reflected in the BOLD signal from the cerebellar cortex. (We also note that increases in the firing rate of inhibitory Purkinje cells means less activation of the neocortex). Thus, while we acknowledge that cerebellar-cortical connectivity likely plays a role in the correlations we observed, we cannot use fMRI observations from the cerebellar cortex and neocortex to draw conclusions about cerebellar-cortical connectivity. To do so we would need to measure activity in the deep cerebellar nuclei (and likely thalamus).

The situation is different when considering the other direction (cortico-cerebellar connections). Here we have the advantage that the cerebellar BOLD signal is mostly determined by the mossy fiber input which, at least for the human cerebellum, comes overwhelmingly from cortical sources. On the neocortical side, the story is admittedly less clear: The cortical BOLD signal is likely determined by a mixture of incoming signals from the thalamus (which mixes inputs from the basal ganglia and cerebellum), subcortex, other cortical areas, and local cortical inputs (e.g., across layers). While the cortical BOLD signal (in contrast to the cerebellum) also reflects the firing rate of output cells, not all output cells will send collaterals to the pontine nuclei. These caveats are now clearly expressed in the Discussion section2.

On balance, there is an asymmetry: Cerebellar BOLD signal is dominated by neocortical input without contribution from the output (Purkinje) cells. Neocortical BOLD signal reflects a mixture of many inputs (with the cerebellar input making a small contribution) and cortical output firing. This asymmetry means that the observed correlation structure between cortical and cerebellar BOLD activity (the determinant of the estimated connectivity weights) will be determined more directly by cortico-cerebellar connections than by cerebellar-cortical connections. Given this, we have left the title and abstract largely the same, but have tempered the strength of the claim by discussing the influence of connectivity in the opposite direction.

2) The authors should address the differences in model complexity, ideally through an additional analysis (e.g., involving AIC or BIC).

We apologize for not having been clearer on the topic of model complexity in the original submission. Our approach (Figure 1) is designed to compare models of different complexity in an unbiased fashion. Parameters for all models were estimated on the training set (Task set A, Figure 1A). This includes the large number of connectivity weights (W_hat) and the single regularization coefficient for Ridge and Lasso regression. The models were then tested on the data from two independent imaging sessions (Task set B, Figure 1B). Using independent training and test data sets is the gold-standard in statistical model comparison (Stone, 1974), and fully corrects for model complexity i.e., any model that overfits the data will, on average, perform worse on an independent test dataset (i.e., Hastie et al. 2009, Chapter 7).

Training and testing was performed within-subjects to account for the considerable inter-individual variation in cortical and cerebellar functional organization. However, to test how well each model generalizes to new subjects, we also averaged the estimated weights across subjects and re-tested the model on a dataset with new subjects performing a new set of tasks (Figure 1C).

AIC is a well-established approximation to the expected cross-validation error in the setting of Gaussian linear regression (Stone 1977). While extensions to WTA (i.e. best-subset selection), Ridge regression and Lasso regression exist, the comparison of models across these approaches using AIC has, to our knowledge, not been well-validated. Thus, we think that using separate training and test sets is the safest and most transparent way to compare models with different structures.

Another important reason to compare models based on their ability to predict independent datasets, rather than using AIC, is that we can assess how well the models generalize to new task conditions. AIC (and related methods) approximates the expected cross-validation error in predicting Y (the cerebellar activity patterns) given a fixed X (the cortical activity patterns). However, in our case the number of tasks is smaller than the number of cortical parcels. Therefore, there will be strong linear dependencies between the different cortical parcels, such that a model may associate a cerebellar region with the wrong combination of cortical parcels (see pt 3). To protect against this, we really want to test on a new random X (e.g., on new tasks). New task sets will lead to somewhat different co-dependencies between cortical regions, and a model that learned the wrong combination of cortical parcels should perform poorly on these new (random) X’s, but not on the old (fixed) X’s. We therefore think that the out-of-sample test on a new taskset is a much stronger criterion for model comparison.

We now point this logic out more clearly in the Introduction and Results. We hope this is a more appropriate and convincing way of dealing with the problem of model comparison than AIC.

References:

Hastie, T., Tibshirani, R., and Friedman, J. (2009). The Elements of Statistical Learning. Second Edition.

Stone, M. (1974). Cross-validatory choice and assessment of statistical predictions, Journal of the Royal Statistical Society Series B 36: 111–147.

Stone, M. (1977). An asymptotic equivalence of choice of model by cross validation and Akaike’s criterion, Journal of the Royal Statistical Society Series B. 39: 44–7.

3) More details on the connectivity method used, especially how it handles the task-evoked variance inflation problem.

We have expanded the description of our connectivity approach in the revised manuscript. It is important to stress that our approach only relied on task-evoked activity. Rather than removing the task-related activity from the time series (as for example, Cole et al., 2019) or using the entire time series, we used only the fitted (and hence task-evoked) activity. All other fluctuations were removed using the crossed approach of predicting cerebellar activity from one half of the data based on cortical activity from the other half of the data for each task set (see response to reviewer #1, comment #4). Models estimated based on the fitted (task-evoked) activity were able to predict the cerebellar activity patterns in Task Set B better than models based on the entire time series or models based on the residuals after removing task-evoked activity.

A task-based approach comes with two important concerns: First, for a small number of tasks, collinearity between different cortical areas may prevent us from determining the correct connectivity model. Second, the particular choice of tasks may bias the estimate of connectivity, and the models may not generalize well to other tasks.

In this paper, we addressed these problems by deliberately choosing a very broad task-set, one designed to activate different cortical areas as independently as possible. Additionally, we evaluated connectivity models in how well they predicted activity patterns from tasks that were not included in task set A (see also pt 2) – this ensures that we pick a model that generalizes well.

To determine how successful we were in picking tasks that allowed us to estimate a comprehensive connectivity model, we calculated the variance-inflation-factor (VIF) for each cortical parcel for the Ridge regression model (Marquardt, 1970). The VIF tells us how variable our estimate of the connectivity weight to a particular cortical parcel is due to collinearity with other cortical parcels. All cortical activity profile (columns of X) were z-standardized before entering into Ridge regression, such that:

VIF=(XTX+Iα)1XTX(XTX+Iα)1

Author response image 1 shows the average VIF for all 1848 cortical parcels (log-α = 8). The collinearity was most severe in anterior and inferior temporal regions, orbitofrontal cortex and some regions of the insula. This indicates that our task battery did not systematically activate and/or dissociate activity in these regions. Importantly, however, this limitation is not observed in the rest of the brain. The VIF in visual, motor, and associative regions was fairly uniform. Thus, the differences in convergence that we found between hand-motor regions and more cognitive regions of the cerebellum does not appear to be caused by different degrees of collinearity in these regions.

Author response image 1

To further address this issue, we conducted an additional simulation analysis using the Lasso model and 80 cortical parcels (as in Figure 4b). For each MDTB region and participant, we simulated artificial cerebellar data using the activity profile from each cortical parcel that was most similar to that parcel and added Gaussian noise (see methods, model recovery simulation). We then estimated a connectivity model (Lasso, log λ = -3, 80 parcels) and calculated the percentage of cortical area with non-zero connectivity weights.As can be seen from Author response image 2, the cerebellar regions showed comparable input areas. Thus, the increased convergence observed for cerebellar region 3-10 was not an artifact of the collinearity of the cortical regressors.

Author response image 2

References:Cole, Michael W., Takuya Ito, Douglas Schultz, Ravi Mill, Richard Chen, and Carrisa Cocuzza. 2019. "Task Activations Produce Spurious but Systematic Inflation of Task Functional Connectivity Estimates." NeuroImage 189 (April): 1-18.

Marquardt, D. W. (1970). Generalized Inverses, Ridge Regression, Biased Linear Estimation, and Nonlinear Estimation. Technometrics, 12(3), 591. https://doi.org/10.2307/1267205

4) The authors should discuss the possibility that motor tasks are inherently simpler than non-motor tasks, which may have biased their results regarding motor vs. non-motor regions.

We do not believe that the difference in convergence between hand-motor and cognitive regions of the cerebellum are biased by complexity differences between motor and non-motor tasks. As we discuss in pt. 3, almost all regions of the cortex (including primary, premotor, and supplementary motor regions) were systematically activated and dissociated by our broadband set of tasks.

The above is an empirical argument. Theoretically we recognize there are many ways of defining complexity, a non-trivial exercise. We included two basic motor tasks in our battery, simple finger tapping and finger sequence tasks. Both activated primary and somatosensory motor cortices (as did most of our cognitive tasks since many required a left- or right-hand response). The main difference between the two motor tasks was that the sequence task elicited higher activation across these regions, as well as elicited activation in premotor and parietal areas. We do think it will be interesting in future research to include a broader set of motor tasks (e.g., include foot and tongue movements). Nonetheless, we do think we have a reasonable manipulation of complexity in the motor domain and yet, still find a difference in convergence between the motor and cognitive regions. Moreover, based on our statistical analysis of collinearity across the cortex (pt. 3), we show that differences in convergence between motor and non-motor regions cannot be attributed to different degrees of collinearity in these regions.

Reviewer #2 (Recommendations for the authors):

The authors should compare model performance while accounting for model complexity. One solution is for the Lasso/Ridge models to have the "WTA" weights to be fixed while randomly permuting all other weights. If the Ridge model is truly better, you would expect this random model to perform similarly to the WTA model but worse than the original Ridge/Lasso models. This would suggest that more inputs do not significantly explain more variances in the cerebellum just by chance.

Our approach strongly controls for the differences in complexity between the WTA, Lasso and Ridge models. We elaborate on this in our response to General Comment #2.

The authors should expand the discussion on the closed-loop anatomical relationship between the cerebellum and the cortex and how its multi-synaptic relationship can allow converging inputs from the cortex to the cerebellum.

This is an important point. We have expanded on this in the revised Discussion section and provide an extended response in General Comment #1.

For figures 3 and 4, it would be nice to expand and show physically where and how widespread those cortical inputs are. That can help readers contextualize the physical level of "convergence" of the cortex onto cerebellar regions.

As noted above, we have reversed the order of Figures 3 and 4, using the examples to now set the stage for the quantitative metrics concerning “physical spread” of cortical inputs. We hope that this reordering provides a clearer picture of connectivity between cortex and cerebellar regions, and allows for a better contextualization of the dispersion metrics used to quantify convergence.

Instead of interpreting the results as cortico-cerebellar inputs, discuss how they can be bidirectional interactions.

We have added an extended discussion of this important point (page 18, see response to General Comment #1).

Reviewer #3 (Recommendations for the authors):

– The authors note that a "convergent architecture would suggest that subregions in the cerebellum integrate information across disparate cortical regions and may coordinate their interactions." It would be interesting to evaluate the degree to which the "disparate" regions that may provide convergent inputs to the cerebellum, as identified in this study, are functionally/anatomically related.

This is an interesting issue. Many cerebellar regions appear to receive convergent inputs from disparate cortical regions that are functionally and/or anatomically connected. For example:

  • MDTB Region 3 (oculomotor vermis) was best explained by a large and relatively dispersed set of cortical regions, including intraparietal sulcus, the frontal eye fields (FEF), and extrastriate visual areas. There are established anatomical connections between FEF and extrastriate visual regions in the macaque (Schall et al. 1995). Moreover, neurostimulation studies have demonstrated that activity in the human FEF has a direct effect on sensitivity of extrastriate visual regions, suggesting a causal relationship between these regions (Silvanto et al. 2006).

  • MDTB regions 7 and 8, which include left and right lobule IX and Crus I/11 respectively, are explained by a relatively restricted set of cortical regions, many of which overlap with subregions of the default mode network (angular gyrus; middle temporal cortex, middle frontal gyrus and inferior frontal gyrus; Smallwood et al. 2021). Previous work has suggested a similar relationship between lobule IX and subregions of the DMN (Guell et al. 2018). Further, based on our previous work, we find that this combination of cerebellar regions is often activated in DMN-like tasks such as rest and movie watching (King et al. 2019).

References:

Schall, J. D., Morel, A., King, D. J., and Bullier, J. (1995). Topography of visual cortex connections with frontal eye field in macaque: convergence and segregation of processing streams. Journal of Neuroscience, 15(6), 4464-4487.

Silvanto, J., Lavie, N., and Walsh, V. (2006). Stimulation of the human frontal eye fields modulates sensitivity of extrastriate visual cortex. Journal of neurophysiology, 96(2), 941-945.

Smallwood, J., Bernhardt, B. C., Leech, R., Bzdok, D., Jefferies, E., and Margulies, D. S. (2021). The default mode network in cognition: a topographical perspective. Nature reviews neuroscience, 22(8), 503-513.

Guell, X., Schmahmann, J. D., Gabrieli, J. D., and Ghosh, S. S. (2018). Functional gradients of the cerebellum. ELife, 7, e36652.

King, M., Hernandez-Castillo, C. R., Poldrack, R. A., Ivry, R. B., and Diedrichsen, J. (2019). Functional boundaries in the human cerebellum revealed by a multi-domain task battery. Nature neuroscience, 22(8), 1371-1378.

– The authors discuss findings from Pisano et al. 2021 as anatomical evidence for cerebellar outputs from specific regions reaching a variety of cortical targets. However, it is important to note that in most species, and especially in the mouse, it is quite unlikely that any one anatomical tracer injection into the cerebellar cortex can target a functionally distinct region of the cerebellum. See, for example, Fujita et al. 2010 J Comp Neurol where small injections into marmoset cerebellar cortex are difficult to localize to specific longitudinal bands in the cerebellar cortex.

We appreciate the reviewer’s point that tracing methods have their own limits in terms of drawing conclusions about connectivity, and thank the reviewer for pointing us to the Fujita paper as providing an example of the challenges in targeting an injection site. In one sense, these limitations make an argument for the utility of our approach -- namely, to use functional data as a way to look at the convergence question (especially since this method makes no claims about whether the cortical input is direct) (via the pons) or indirect (via projections to other cortical areas and then to the pons). We have opted to not make any changes to the manuscript with respect to this point since it seemed secondary to the inferences we make in this paper.

– A map of the cerebellum showing the MDTB regions would be useful to include as a guideline in Figure 3

As discussed in our response to Reviewer 2, we have reversed the order of Figures 3 and 4 in the revised manuscript. We did this to provide examples of differences in the physical spread of cortical regions before turning to the quantitative analysis shown in Figure 4 (formerly Figure 3). Given that the MDTB regions are included in Figure 3, the reader is now provided with this map as a guideline for Figure 4.

[Editors' note: further revisions were suggested prior to acceptance, as described below.]

The manuscript has been improved but there are some remaining issues that need to be addressed, as outlined below:

Reviewer #1 (Recommendations for the authors):

The authors did an excellent job improving their manuscript based on previous reviews, except for one major point described below.

R1 point 3 was responded to, but not meaningfully addressed. The authors seemed to misunderstand the issue that regressing out mean task-evoked activity is meant to deal with, resulting in the authors effectively ignoring the issue with regard to their own study. Indicative of misunderstanding the issue, the authors state: "Returning to the papers cited by the reviewer, these are designed to look at connectivity not related to the task-evoked activity." This is a common misunderstanding in the field: Rather than looking at connectivity not related to task-evoked activity, regressing out mean task-evoked activity is meant to deal with a causal confound (common inputs from stimuli), while still allowing estimation of connectivity related to the task-evoked activity. This interest in the impact of the task-related activity on connectivity is clear from the simple logic of most studies utilizing mean task-evoked activity regression, wherein connectivity estimated in the context of condition-specific task-evoked activity is compared. (If all task-evoked activity were removed then there would be no point in comparing task conditions.) What is left over for estimating functional connectivity estimation is trial-by-trial (and moment-to-moment) variation in task-evoked activity (along with spontaneous activity). So, the authors have not dealt with the possibility that common stimulus inputs – such as centrally presented visual stimuli simultaneously activating left and right hemispheres, generating false inter-hemispheric functional connections – are driving some of their connectivity estimates. Ideally, this confound would be dealt with, but at the very least the authors should acknowledge this limitation of their methods.

We apologize for the confusion caused by our response. When we said “task-evoked” activity in our last response, we meant the averaged (across repeated trials) task-evoked activity, not the trial-by-trial fluctuations in the task evoked activity that you are talking about here. The connectivity model training and evaluation in our paper is relying exclusively on this averaged task-evoked activity (see below).

As you state, this of course raises the immediate concern that a high connectivity weight between region A and B could be caused by a task that activates both regions, without the regions being connected. For example, the verbal N-back task is always associated with responses of the left hand, so if we were estimating connectivity based on this task alone, we could find high “connectivity” between the right hand-region and the working memory region. However, our task set was designed to be very broad, optimally dissociating the activity profiles across cortical regions (see Figure R1 in our last response). For example, the task set also included an object N-back task that was associated with right-hand responses. Thus, we addressed this concern by deliberately choosing a very broad task-set. Using data from this task set should strongly reduce the influence of spurious associations. Additionally, the evaluation on an independent set of tasks ensures that connectivity models that do reflect spurious between-region correlations that are induced by the specific choice of tasks, are penalized and will perform worse.

Thus, our approach follows a similar logic as the papers discussed in the last response.

The standard approach is to rely on the fact that the neural fluctuations captured in the residuals (or in the trial-by-trial variations of the regression weights) are rich and varied enough, such that two areas that are not connected would not always correlate. Here, we are relying on the fact that our task set is rich and varied enough to ensure that two areas that are not connected would not always be co-activated. To be clear, both approaches do not guarantee that the estimated connectivity between region A and region B is not caused by a consistent indirect influence from a third region.

Beyond this issue (and perhaps more importantly), it is difficult to reconcile the present study's approach with more standard task-related functional connectivity approaches. How can connectivity values estimated via the present study's approach be interpreted? How do these interpretations differ from those of more common task-related functional connectivity approaches? Is there enough advantage to using the present study's unconventional task-related functional connectivity method to make it worth the extra effort for readers to make sense of it and learn to relate it to more standard approaches? The presence of all of these unanswered questions makes it difficult to properly evaluate the manuscript and its results.

Our approach indeed deviates in another aspect from other connectivity approaches. Rather than using correlation between time-series to derive a network structure, it runs a multiple regression analysis to explain the task-evoked activity in the cerebellum from the task-evoked activity from the neocortex. This approach is justified by the particular neuroanatomical and physiological properties of the system under study: (1) BOLD signal in the cerebellar cortex reflects almost exclusively the incoming input to the cerebellum, which in the human is dominated by the input from neocortex, and (2) there is very little recurrent activity and there are few association fibers between different regions of the cerebellum.

We therefore do not believe that our approach could not be used to model cortico-cortico connectivity. However, our chosen approach is perfectly suited to model the cortico-cerebellar system. We have now pointed out this difference, and the justification for our approach, more clearly in the manuscript.

Further, in my attempt to answer the above questions myself, I have read and re-read the description of the connectivity approach, and yet I have not been able to make sense of it. I doubt this is my own lack of knowledge or intelligence, but even if that were the case it is likely that many (perhaps most) readers would have the same issue. Therefore, I recommend rewriting the "Connectivity Models" section (within "Methods") to include more details. Some specific details that would be helpful to include:

– The meaning of "fitted time series" is ambiguous. Does this mean the original time series multiplied by the fit β estimate for each condition, resulting in the original time series modulated by the fit to the task regressor for each condition? Or does this mean the regressor time series multiplied by the fit β estimate for each condition? If it is the first option I can imagine interpreting the resulting functional connectivity estimates in a fairly similar way to standard task-related functional connectivity. However, if it is the second option, then I think the interpretation would be quite different from standard task-related functional connectivity. In such a case I think much more details on how to interpret the resulting connectivity values would be warranted. Of course, it could be a third option I haven't considered, in which case it would be helpful to have that option detailed.

– It is unclear what this detail means (pg. 24): "This procedure reweighted the activity estimates in an optimal manner, accounting for the fact that estimates for the instructions were based on 5 s of fMRI data, while estimates for the conditions were based on 10-30 s of data." In what sense is this accounting for the duration of task conditions? This makes me wonder whether the connectivity estimates were based on using cortical β estimates (rather than time series) to predict cerebellar β estimates (rather than time series), which could account for condition duration quite well. But the description is ambiguous with respect to these kind of details, and this possibility doesn't seem to fit well with earlier descriptions of the method.

Your second interpretation is correct. To be perfectly clear, the GLM decomposes the time series of each voxel (t) for each run (r) and session (s) into a fitted time series (Z bsr) and a residual time-series (rsr). The fitted activity can be further decomposed into the part of the time series that is explained by the activity for each condition (averaged across runs), and a part that is explained by the run-to-run variations around the mean activity for each condition:ts,r=Zs,rb¯s+Z(bs,rb¯s)+rs,r

Where Zs,r is the (n_timepoints x n_regressors) design matrix for session s and run r. bs,r is the vector of activity estimates for a specific run, and is the average vector of β weight sacross runs for session s. For the connectivity analysis, we used only the first term – the part of the response that was explained by the average task response. We have now clearly pointed out in methods section (see connectivity models) how this approach differs diametrically from other connectivity approaches that remove the fitted time series (Z b) and use the residual time series (r) only to estimate connectivity.

– It is unclear what this detail means (pg. 24): "This method also mitigates problems that may arise in event-related designs due to correlations between regressors (Mumford et al., 2012)." How exactly? Also, what exactly are the problems that may arise in event-related designs due to correlations between regressors? It sounds like this is related to collinearity, but it is unclear which aspects of the method help with collinearity.

We apologize for not being clearer here. As you point out, our approach of using the fitted time-series (based on averaged β weights) is conceptually very similar to using the vector of averaged activity estimates to estimate the connectivity model (indeed this is what we did first). The difference is that the former approach reweights the activity estimates based on the impact the activity estimate has on the predicted time series.

This can be clearly seen when we calculate the variances and covariances between cortical and cerebellar data, which are the only quantities that matter for the connectivity model. For example, when we compute the inner product between the data for a cortical (x) and a cerebellar voxel, we have:

yTx=b¯xTZTZb¯y

If we had used only the average activity estimates (b), we would have:yTx=b¯xTb¯y

Comparing these two equations shows that the procedure simply reweights the activity estimates by ZTZ. The diagonal of this matrix contains the squared length of the support of the regressor (i.e. it accounts for the fact that estimates for the instructions were based on 5 s of fMRI data, while estimates for the conditions were based on 10-30 s of data).

The off-diagonal of ZTZ contains the covariance between the regressors. For example, due to the overlap of the instruction period and the following condition, these two regressors are positively correlated. Consequently, the activity estimates for the instruction and following condition are negatively correlated, a problem that arises especially in event-related designs (Mumford et al., 2012). This estimation uncertainty, however, is fully cancelled out when multiplying these estimates with the design matrix again.

Overall, this is a very minor and technical point, and we fear that it would confuse the reader rather than clarify our approach. We have therefore chosen to remove this point from the methods.

– The details regarding the "crossed" approach need to be fleshed out. It is stated (pg. 25): "we used a "crossed" approach to train the models: The cerebellar fitted time series for the first session was predicted by the cortical time series from the second session, and vice-versa". This procedure is also presented in an ambiguous manner in Figure 1. Right now it is stated in the manuscript that the session 2 cortical time series are used to predict the session 1 cerebellar time series, and yet it is unclear how that would be possible when there are different task conditions across the sessions. Is the connectivity model built using both cortical and cerebellar time series during session 2 (betas from cortex predicting cerebellum), then both cortical and cerebellar time series are used from session 1 to test the session 2 connectivity model (session 2 betas multiplied by session 1 cortical time series to predict session 1 cerebellar time series)? Put more simply: is it the case that the session 2 connectivity model is tested using session 1 data, and that the session 1 connectivity model is tested using session 2 data? If so, the description should change to better describe this. In any case, the description should be clarified with more details.

We have now added more detail explaining the “crossed approach” (see Connectivity Models). We have also made it clear throughout the document that there are two sessions (same tasks, separate days) for each task set: the same tasks are performed across sessions, with the order of the tasks randomized across sessions (this is true for both task sets). However, using the averaged activity estimates for each session, and the task-related regressors from the other session, we can predict the time series we should have observed if the task activated the brain in the same way across sessions. Therefore, it is possible to predict the cerebellar time series from one session using the cortical time series based on the activity estimates from another session.

https://doi.org/10.7554/eLife.81511.sa2

Article and author information

Author details

  1. Maedbh King

    Department of Psychology, University of California, Berkeley, Berkeley, United States
    Contribution
    Conceptualization, Data curation, Software, Formal analysis, Validation, Visualization, Methodology, Writing – original draft, Project administration, Writing – review and editing
    Contributed equally with
    Ladan Shahshahani
    For correspondence
    maedbhking@gmail.com
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0001-5374-1011
  2. Ladan Shahshahani

    Western Institute for Neuroscience, Western University, London, Canada
    Contribution
    Software, Formal analysis, Validation, Visualization, Methodology, Writing – original draft, Project administration, Writing – review and editing
    Contributed equally with
    Maedbh King
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0001-6189-9994
  3. Richard B Ivry

    1. Department of Psychology, University of California, Berkeley, Berkeley, United States
    2. Helen Wills Neuroscience Institute, University of California, Berkeley, Berkeley, United States
    Contribution
    Conceptualization, Resources, Supervision, Funding acquisition, Investigation, Writing – original draft, Writing – review and editing
    Competing interests
    is a co-founder with equity in Magnetic Tides, Inc
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-4728-5130
  4. Jörn Diedrichsen

    1. Western Institute for Neuroscience, Western University, London, Canada
    2. Department of Statistical and Actuarial Sciences, Western University, London, Canada
    3. Department of Computer Science, Western University, London, Ontario, Canada
    Contribution
    Conceptualization, Resources, Data curation, Software, Formal analysis, Supervision, Funding acquisition, Investigation, Methodology, Writing – original draft, Project administration, Writing – review and editing
    Competing interests
    Reviewing editor, eLife
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-0264-8532

Funding

Canadian Institutes of Health Research (PJT 159520)

  • Jörn Diedrichsen

National Institutes of Health (NS092079)

  • Richard B Ivry

National Institutes of Health (NS105839)

  • Richard B Ivry

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

This work was supported by the Canadian Institutes of Health Research (PJT 159520 to JD), the Canada First Research Excellence Fund (BrainsCAN to Western University), and the National Institute of Health (NS092079, NS105839 to RBI). Thanks to Da Zhi for providing cerebral cortical parcellations.

Ethics

Human subjects: All participants gave informed consent under an experimental protocol approved by the institutional review board at Western University (Protocol number: 107293). Undergraduate and graduate students were recruited (via posters) from the larger student body at Western University. The sample was biased towards relatively high-functioning, healthy and young individuals. The final sample consisted of 24 healthy, right-handed individuals (16 females, 8 males; mean age=23.8 years old, SD=2.6) with no self-reported history of neurological or psychiatric illness.

Senior Editor

  1. Chris I Baker, National Institute of Mental Health, United States

Reviewing Editor

  1. Michael W Cole, Rutgers, The State University of New Jersey, United States

Version history

  1. Preprint posted: May 8, 2022 (view preprint)
  2. Received: June 30, 2022
  3. Accepted: March 31, 2023
  4. Version of Record published: April 21, 2023 (version 1)

Copyright

© 2023, King, Shahshahani et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 892
    Page views
  • 121
    Downloads
  • 4
    Citations

Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Maedbh King
  2. Ladan Shahshahani
  3. Richard B Ivry
  4. Jörn Diedrichsen
(2023)
A task-general connectivity model reveals variation in convergence of cortical inputs to functional regions of the cerebellum
eLife 12:e81511.
https://doi.org/10.7554/eLife.81511

Share this article

https://doi.org/10.7554/eLife.81511

Further reading

    1. Computational and Systems Biology
    2. Neuroscience
    Tony Zhang, Matthew Rosenberg ... Markus Meister
    Research Article

    An animal entering a new environment typically faces three challenges: explore the space for resources, memorize their locations, and navigate towards those targets as needed. Here we propose a neural algorithm that can solve all these problems and operates reliably in diverse and complex environments. At its core, the mechanism makes use of a behavioral module common to all motile animals, namely the ability to follow an odor to its source. We show how the brain can learn to generate internal “virtual odors” that guide the animal to any location of interest. This endotaxis algorithm can be implemented with a simple 3-layer neural circuit using only biologically realistic structures and learning rules. Several neural components of this scheme are found in brains from insects to humans. Nature may have evolved a general mechanism for search and navigation on the ancient backbone of chemotaxis.

    1. Neuroscience
    Frances Skinner
    Insight

    Automatic leveraging of information in a hippocampal neuron database to generate mathematical models should help foster interactions between experimental and computational neuroscientists.