Impact Factor 4.677 | CiteScore 5.4
More on impact ›


Front. Neurosci., 21 September 2021 |

Simulating Local Deformations in the Human Cortex Due to Blood Flow-Induced Changes in Mechanical Tissue Properties: Impact on Functional Magnetic Resonance Imaging

  • 1Department of Neurophysics, Max Planck Institute for Human Cognitive and Brain Sciences, Leipzig, Germany
  • 2Methods and Development Group Neural Data Science and Statistical Computing, Max Planck Institute for Human Cognitive and Brain Sciences, Leipzig, Germany
  • 3Institute for Medical Informatics and Biometry, Carl Gustav Carus Faculty of Medicine, TU Dresden, Dresden, Germany
  • 4Department of Radiology, Charité – Universitätsmedizin Berlin, Berlin, Germany
  • 5Berlin Center for Advanced Neuroimaging, Charité – Universitätsmedizin Berlin, Berlin, Germany
  • 6Berlin Center for Computational Neuroscience, Berlin, Germany
  • 7Faculty of Physics and Earth Sciences, Felix Bloch Institute for Solid State Physics, Leipzig University, Leipzig, Germany

Investigating human brain tissue is challenging due to the complexity and the manifold interactions between structures across different scales. Increasing evidence suggests that brain function and microstructural features including biomechanical features are related. More importantly, the relationship between tissue mechanics and its influence on brain imaging results remains poorly understood. As an important example, the study of the brain tissue response to blood flow could have important theoretical and experimental consequences for functional magnetic resonance imaging (fMRI) at high spatial resolutions. Computational simulations, using realistic mechanical models can predict and characterize the brain tissue behavior and give us insights into the consequent potential biases or limitations of in vivo, high-resolution fMRI. In this manuscript, we used a two dimensional biomechanical simulation of an exemplary human gyrus to investigate the relationship between mechanical tissue properties and the respective changes induced by focal blood flow changes. The model is based on the changes in the brain’s stiffness and volume due to the vasodilation evoked by neural activity. Modeling an exemplary gyrus from a brain atlas we assessed the influence of different potential mechanisms: (i) a local increase in tissue stiffness (at the level of a single anatomical layer), (ii) an increase in local volume, and (iii) a combination of both effects. Our simulation results showed considerable tissue displacement because of these temporary changes in mechanical properties. We found that the local volume increase causes more deformation and consequently higher displacement of the gyrus. These displacements introduced considerable artifacts in our simulated fMRI measurements. Our results underline the necessity to consider and characterize the tissue displacement which could be responsible for fMRI artifacts.


Understanding the complex anatomical and functional architecture of the human brain has been a major research topic for more than a century. However, little is known about the dynamic interaction between mechanical tissue properties and the shape and microscopic composition of cortical substructures. For example, recent studies suggest that the gyrification of the developing brain can be explained by different growth rates of gray and white matter (Tallinen et al., 2016) that causes mechanical stress resulting in the characteristic pattern of gyri and sulci (Armstrong et al., 1995). This early stage of the nervous system’s development is followed by axonal growth patterns that are guided by mechanical cues, i.e., the tissue’s local stiffness gradients (Koser et al., 2016; Oliveri et al., 2021). Abnormalities in cortical folding, neuronal growth, and viscoelastic tissue properties have been associated with various neurological and cognitive disorders (Bayly et al., 2014) such as hydrocephalus, Alzheimer’s disease, multiple sclerosis, epilepsy, schizophrenia, autism, and brain tumors (Hardan et al., 2004; Lin et al., 2007; Wuerfel et al., 2010; Murphy et al., 2011; Streitberger et al., 2011, 2020; Li et al., 2016; Garcia et al., 2018b; Huesmann et al., 2020) which demonstrates a clear need for understanding the underlying tissue mechanics.

The rapid development of magnetic resonance imaging (MRI) methods has enabled the non-invasive observation of the complex structure-function relationship within the living brain at different timescales, i.e., by detecting morphometric or microstructural changes during development, disease, learning (Maguire et al., 2000; Draganski et al., 2004), or during rapid functional stimulation (Watanabe et al., 2017). Recent structural and functional MRI (fMRI) methods are approaching the spatial resolutions needed for resolving single cortical layers (Fukunaga et al., 2010; Trampel et al., 2019; Weiskopf et al., 2021) with a typical thickness around 400 μm (Wagstyl et al., 2020) that are specialized in processing specific neuronal information (Huber et al., 2017; Fracasso et al., 2018). The raised energy demand during elevated neuronal activity induces vasodilation in the activated brain tissue resulting in an increased inflow of freshly oxygenated blood. The primary fMRI contrast mechanisms rely on this vasodilatory response evoked by neuronal activity (Goense et al., 2016). In the context of high-resolution fMRI, the mechanical implications of vasomotion at the scale of gyral substructures have not yet been investigated (Hetzer et al., 2018a; Polimeni et al., 2018). It should be noted that minimal tissue displacements in the range of some percent of the fMRI voxel size can result in erroneous relative fMRI signal changes of the same magnitude, i.e., up to 10% for a 50-μm displacement within a 500 μm voxel at high contrast interfaces. For example, varying partial volume contributions of cerebrospinal fluid (CSF) can significantly impact the cortical fMRI signal due to its long transverse relaxation time as compared to cortical tissue (Piechnik et al., 2009). These mechanical effects are expected to depend on the baseline biomechanical properties of the tissue and may thus depend on pathologies accompanied by tissue degeneration like protein deposition, dystrophic calcification, gliosis, etc. (Yin et al., 2018).

To our knowledge, the extent to which the fMRI signal is influenced by tissue deformations inherent to the mechanical behavior of cortical brain tissue has never been systematically investigated. Therefore, we here combine simulations of tissue biomechanics and the physics of fMRI to study whether subtle mechanical property changes due to vasodilation can result in artifacts in ultra-high resolution fMRI. Two effects of vasodilation on tissue mechanics will be explored.

First, the viscoelastic properties locally change due to vasodilation (Parker, 2017). Recent studies using MR elastography (MRE) (Muthupillai et al., 1995) detected slight changes in viscoelastic properties of human brain tissue during hypercapnia-induced vasodilation (Hetzer et al., 2018b). This is partially predicted by the well-known Fåhræus–Lindqvist effect (Fåhræus and Lindqvist, 1931), which describes the increase of blood viscosity with vessel diameter. Embedded in folded cortical tissue that is under stress (Xu et al., 2010; Tallinen et al., 2016), a focal increase of tissue stiffness due to vasodilation will result in slight displacements of the mechanically interconnected surrounding tissue.

Second, a focal volume change is expected. Within a constant intracranial volume as assumed in the Monro–Kellie doctrine, any variation in cerebral blood volume (CBV) must be compensated by opposite changes in other compartments, or must affect the intracranial pressure (ICP). While clinically induced vasoconstriction is routinely used to reduce ICP in patients with severe head injury (Laffey and Kavanagh, 2002), the mechanical effects of a focal vasodilatory volume increase of activated tissue and the associated displacement of interconnected cortical tissue is expected to be fully compensated by CSF redistributions without significantly elevating the ICP in accordance with the Monro–Kellie doctrine (Piechnik et al., 2009; Jin and Kim, 2010). With a typically observed relative blood volume increase of 50% in proximity of neuronal activation (Lu et al., 2003) the increase of the resting cortical blood partition fraction of 5% (Hua et al., 2011) is expected to translate into a local tissue volume increase of a couple of percent and concomitant shifts.

To estimate the tissue displacements due to vasodilation, we created a two-dimensional geometric model of a representative gyrus from ultra-high resolution histology data of Amunts et al. (2013) and used a finite-element model to simulate the respective tissue mechanics. Following recent works realistically modeling the nearly incompressible brain tissue (Franceschini et al., 2006; Bayly et al., 2013; Budday et al., 2014) our model includes neo-Hookean material properties for white matter (WM) and cortical gray matter (GM) resembling six neocortical layers. For simplicity, neither CSF compartment changes, CSF redistributions due to focal deformations of cortical tissue nor gravitational forces were explicitly modeled as the brain suspended in CSF exists in neutral buoyancy. A quasi-static simulation with hyperelastic material was chosen as an appropriate model for slow deformations or the equilibrium state.

To predict the sensitivity of typical high-resolution fMRI sequences to the analyzed biomechanical effects, we simulated the MR imaging process at different image resolutions on an idealized 7-Tesla scanner with a high-performance gradient system which will be realistically available within the next years and used to study the intracortical functional neuroanatomy (Brain Initiative, 2021).

Materials and Methods

Anatomical Model

From ultra-high resolution histology data at nearly cellular resolution of 20 μm from the HumanBrainProject (Amunts et al., 2013), we created a two-dimensional geometric model of a coronal slice of the right middle frontal gyrus (at MNI coordinate 24, 30, and 41) of 8.5 mm width and 15.3 mm length (Figures 1A,B).


Figure 1. (A) Original microscopy image of the right middle frontal gyrus at MNI coordinate x = 24 mm, y = 30 mm, z = 41 mm taken from the HumanBrainProject. (B) The geometric model extracted from (A) showing gray matter (GM) (gray) with different anatomic layers (different gray levels), white matter (WM) (cream) and cerebrospinal fluid (CSF) (blue). Simulation of local changes in tissue mechanics: the stiffness and volume of a segment of interest (here superficial layer L1 highlighted in red) were increased to simulate the effect of local increase in blood flow.

We traced the cortical boundaries of the gyrus in the histological image and imported the boundaries for subsequent finite element simulations in COMSOL Multiphysics (V.5.5, COMSOL Multiphysics Göttingen). We modeled the geometry with an adaptive triangular mesh with elements size ranging from 0.38 to 190 μm. To ensure high numerical accuracy, the extremely fine mesh size was used for the cortical gyrus region in our study. The six cortical layers in GM were manually traced by an expert with more than 10 years of experience in histology and neuroanatomy (CJ) and later incorporated into the model. The geometry also includes the surrounding GM, WM, and CSF. Note that we did not explicitly model the biomechanics of CSF in this study. To separately specify parameters locally in different parts of the layered geometry during simulation and analysis, we divided the model further into segments (Figure 1B).

To assess potential simulation resolution effects, the simulations were repeated using different mesh element sizes from normal to extremely fine. The mesh size did not influence our results much: for fine with average mesh element size of 206 μm and extremely fine mesh with average mesh element size of 143 μm, the difference in resulting displacement was less than 0.05% (20 nm) indicating that convergence had been attained.

Biomechanical Simulations

We modeled the brain tissue as an isotropic and hyperelastic material as suggested in Xu et al. (2010) and Tallinen et al. (2016). We assumed a quasi-static regime disregarding any dynamic effects, such as viscosity or hydrodynamic interactions between parenchyma and CSF. Accounting for dynamic effects would dramatically increase the complexity of the physical model and the number of coupling parameters and is therefore left for future studies. We use the strain-energy density function W of a nearly incompressible neo-Hookean material (Holzapfel, 2000) defined as


with μ and κ being the shear and bulk modulus of the material, F the deformation gradient and


Using this, the constitutive stress-strain relation is given by


where nearly incompressible tissue properties are achieved when κ≫μ.

The cortical gyrus is modeled and simulated in COMSOL for characterizing blood flow-induced motion and volume change by means of biomechanical tissue properties. Its material properties were taken from Budday et al. (2015) and Mcintosh and Anderson (2010) for the GM with a shear modulus of 1.4 kPa, approx. 35% stiffer WM (1.9 kPa) and assuming a nearly incompressible material with κ three orders of magnitude larger than μ (Xu et al., 2010). An overview of the specific values used for simulation are given in Table 1.


Table 1. An overview of the material parameters used for mechanical simulation.

In order to gauge the effect sizes and contributions of different tissue characteristics to the fMRI signal, we studied three computational models simulating different biomechanical mechanisms: stiffness change, volume change, and a combination of both. Based on recent research using direct measurements, we know that there is tension in WM and compressive stresses in cortical GM in both smooth mouse brain and folded ferret brain (Xu et al., 2009, 2010; Garcia et al., 2018a). This implies that the brain scanned by in vivo MRI methods is under stress already. We estimated the compressive stress distribution in the gyrus based on the assumption that flattening the folded cortex results in a realistic stress distribution. The underlying simplifying assumption is that the folding pattern results from simple mechanical effects without tissue growth or other influences. In order to obtain the flattened geometry, the input geometry was deformed back to a smooth and unfolded shape. The inner (lower) and outer (upper) gyrus boundaries were subjected to a prescribed displacement in y direction (top-down direction on Figure 1B); it means they were straightened into two parallel lines with a given distance. This was done in small steps in which the gyrus became more and more deformed until these boundaries were flat and even (see Supplementary Video 1). During this procedure, the stress inside the gyrus, especially in curved regions, would increase as well as the overall length, as one would expect. The stresses arising from this deformation were mapped back to the original gyrus shape. The resulting stress distribution was mapped onto the folded input geometry and used for the starting configuration for the simulation studies involving stiffness change.

As a first computational model, solely stiffness changes were taken into account. Then, only focal volume changes were considered. Finally, we combined both volume and stiffness changes in our third model. The following steps and parameters were used for these three simulations:

1. To simulate the effect of stiffness change, we locally increased the shear modulus μ1 of one specific segment of a single layer (highlighted in red in Figure 1B) by up to 10% of its initial value μ0, i.e., μ1 = 1.1 μ0. To trace the dynamics of the resulting shape deformation, we swept the stiffness value from μ0 to μ1 in equidistant steps of 0.5%.

2. The same segment as used in the previous model (Figure 1B) underwent a volume increase by 10%. The volume change occurred along two orthogonal directions in plane (i.e., the extent in the through plane direction was implicitly kept constant). The volume increase was governed by expanding the lateral dimension uniformly in equal steps of 0.5%.

3. In the final model, both volume and stiffness change were combined in an individual simulation. The two parameters were changed simultaneously up to 10%.

To assess the sensitivity of our model on the chosen parameter values and to reasonably determine the tissue mechanical properties, we set up a grid simulation that captures a large range of the parameters κ and μ obtained from literature (Budday et al., 2016, 2020; Ganpule et al., 2018).

For quantifying the resulting tissue displacements, a displacement field with the coordinates u, v was defined. At a given coordinate (X, Y), e.g., after deformation of the geometry, the values u and v quantify the displacement in x- and y-direction, respectively, relative to its original position. The total displacement was then defined asdtot=u2+v2.

Functional Magnetic Resonance Imaging Simulations

To predict the sensitivity of different typical high-resolution fMRI acquisitions to biomechanical tissue displacements due to vasodilatory effects, the gyral geometries and their changes resulting from the mechanical simulations were transferred to MRiLab,1 an open-source software simulating the entire MR imaging process (Liu et al., 2017).

The three most widely used fMRI variants (Goense et al., 2016; Huber et al., 2019) were simulated, i.e., Vascular-Space-Occupancy (VASO), blood-oxygen-level-dependent (BOLD), and proton density contrast as used in arterial spin labeling (ASL), for the next generation 7-Tesla scanner with a high-performance gradient system (maximal gradient amplitude and slew rate: 200 mT/m, 1,000 T/m/s), which is tailored for studies of the mesoscopic cortical functional neuroanatomy, and idealized radio frequency (RF) coils and B0 homogeneity employing a gradient-echo readout with optional inversion recovery preparation.

Proton density values and relaxation times at 7 Tesla for CSF, WM and cortical layers (see also Table 2) were taken from the literature (Rooney et al., 2007; Tofts, 2010; Carey et al., 2018; Caan et al., 2019; McColgan et al., 2021) and mapped to the mechanical gyrus model.


Table 2. Values for proton density (PD), T1, T2* for different tissue types and cortical layers (L1–6) at 7 Tesla used in the simulations.

For gyral geometries resulting from the analyzed focal stiffness and volume changes, magnitude MR images were reconstructed from the simulated k-space data at different isotropic image resolutions (voxel sizes of 0.25, 0.5, and 1.0 mm) of the following fMRI sequence types:

(a) Proton density contrast was created at a short echo time TE = 5 ms to reflect the effort of blood flow measurements to suppress unwanted concomitant BOLD signal, i.e., in ASL (Hetzer et al., 2011; Huber et al., 2019).

(b) BOLD imaging was simulated by calculating the image contrast at TE = T2 of GM where BOLD sensitivity is optimal (Deichmann et al., 2002).

(c) VASO contrast (TE = 5 ms) was simulated by an inversion recovery preparation module for nulling blood with T1b = 2,100 ms (Zhang et al., 2013) at an inversion time of TI = ln(2)T1b = 1,456 ms.

Voxel-wise maps of the relative signal change resulting from tissue displacement were calculated from the corresponding MR images of the modeled gyrus in the resting and “activated” conditions, i.e., by subtracting the rest from the activated condition image.


Mechanical Simulations

Cortical Stress Distribution

From our biomechanical simulation employing the flattening approach, we obtained an estimate of the stress distribution shown in Figure 2. The stress distribution indicates parts with negative stress which correspond to compression in high curvatures and positive stress distribution which corresponds to tension in low curvatures. The arrow shows the stress gradient from highest to lowest stress value.


Figure 2. (A) Stress distribution results of the first simulation step of the geometry and (B) a zoomed region where positive curvature results in a positive stress gradient across layers. The stress was normalized by the stiffness value μ of the GM.

Stiffness Change

Using the model set up in “Cortical Stress Distribution,” we conducted a simulation study to estimate the effect of local stiffness change (e.g., due to blood flow) and the resulting deformation of the cortical tissue. Figure 3A shows the simulation results. One can observe that the largest displacement took place in the gyrus crown, and the displacement decreased with the distance from the gyrus peak. The displacement of the gyrus correlated with the amount of stiffness change. We observed a maximum total displacement of 41.3 μm for the gyrus crown with a 10% stiffness change.


Figure 3. Resulting tissue movement for a 10% increase of (A) tissue stiffness, (B) tissue volume of the activated layer segment, and (C) the combination of both effects. Displacement fields are shown by black arrows and color scale indicates magnitude of displacement.

Volume Change

The displacement pattern arising from a local volume increase of the activated gyral segment was similar to the displacements induced by a local stiffness increase (Figure 3B), with an estimated 69.5 μm total displacement for a volume increase of 10% in the gyral crown. This effect was hence larger than the local stiffness change (Figure 3A).

Combined Stiffness and Volume Change

Figure 3C shows the results for the combined effect of simultaneous increase in stiffness and volume by 10% each, resulting in a larger displacement compared to the individual effects. The total displacement calculated in the gyrus crown was 110.7 μm.

The total displacement of all three models as a function of the parameter change is shown in Figure 4. The plot shows that the total displacement increases linearly with increasing stiffness or volume change with the local volume increase resulting in a larger displacement. Our simulations show that the combination of these two mechanisms are additive (to within an accuracy of <0.1 μm).


Figure 4. (A) The total displacement for the gyrus crown as a function of parameter change for the three simulation studies. (B) The spatial map of the resulting displacement given in per percent volume change. Color map shows the magnitude of local displacement.

Robustness of the Simulation and Sensitivity to Tissue Parameter Choices

We conducted a simulation in which we investigated a broader range of mechanical parameters summarized in Table 3. The results are shown in Figure 5 and Supplementary Figure 2. The stiffness change was kept constant (2 and 10%, respectively) and the total displacement dtot as a function of the bulk modulus κ and the baseline stiffness μ of the gyrus was calculated (Supplementary Figure 1). The values for dtot ranged from 4.69 to 23.90 μm for a change in stiffness of 2% and from 23.00 to 117.50 μm for a stiffness change of 10%.


Table 3. Parameter range for mechanical properties based on the literature (Xu et al., 2010; Budday et al., 2020).


Figure 5. (A) Surface plot of the total displacement for the studied region as a function of baseline stiffness and stiffness change in percent for bulk modulus κ = 5 GPa. The baseline stiffness μ = 1,400 Pa is used for all the previously discussed studies. (B) The effect of scaling the stress field on tissue displacement for the three simulation cases. The dotted line indicates the stress resulting from the flattening procedure employed in this work.

The influence of κ on the displacement of the gyrus is negligible since the material model is nearly incompressible. In Figure 5 the bulk modulus was kept at a fixed value of κ = 5 GPa and dtot of the gyrus as a function of its stiffness value μ and its stiffness change (in %) is shown.

Furthermore, in order to gauge how the baseline stress distribution impacts the results, the simulations were repeated with different scales of obtained stress distribution from the flattening step. The results in Figure 5B show that the displacement values depend on the global scaling of the stress distribution. The dotted line shows the scaling factor of 1 which was used for all other simulations. The internal stress distribution was not relevant for the volume change simulation.

To get a deeper understanding of how the layer depth influences the biomechanical response of the cortex, we repeated the simulations for the same lateral segment across all layer depths (outermost layer i = 1 to innermost layer i = 6). The results for a change in stiffness μ, volume and the combined effect by 10% are shown in Figure 6. Figure 6A shows how dtot changes with increasing index i, i.e., outer to inner layer, when changing the segment from outer to inner layer. One can observe that the larger displacement occurred in layer 3 for all three cases, as expected with the maximum for the combined effect. To rule out the effect of the segment size we normalized the resulting dtot of a segment i by its area Ai as shown in Figure 6B. The total displacement dtot per area Ai clearly decreases with increasing segment index i, with a factor of more than 2 between the innermost (i = 6) and outermost (i = 1) layer. This effect was largest for the combined stiffness and volume change.


Figure 6. Influence of layer position within the cortical segment on (A) absolute displacement for the three different scenarios and (B) displacement normalized by the layer segment volume.

To investigate how the effects depend on the size of the region experiencing blood flow changes, we extended the simulation region from one segment to four. In all three scenarios, the extended activation of segments resulted in a larger displacement. For example, when the activated segment was extended up to layer 4, the total displacement for volume increase simulation increased from 69.5 to 319.7 μm. Also, when the position of the broad simulated segments was moved to other regions, it resulted in larger displacements compared to a single segment.

Changing the position of the activated segment for the first layer in tangential direction, the total displacement in the gyrus crown decreased for all the three cases (Supplementary Figure 3A). By extending the length of the activated region from a single segment to four segments, we observed an increase in the displacement value for all the three simulation cases. For example the displacement value for segment one of layer one for stiffness change simulation was 29.9 μm. However, the displacement value increased as we increased the segment number, reaching 74.4 μm (Supplementary Figure 3B).

Spurious Functional Magnetic Resonance Imaging Signal Changes Due to Displacement

Raw images of the gyrus in the resting condition for the three analyzed fMRI sequence variants at three image resolutions are shown in Figure 7A. The corresponding maps of spurious fMRI signal changes due to a tissue displacement following a 2% volume increase of the “active” layer segment (marked red in Figure 1B) can be seen in Figure 7B.


Figure 7. (A) Simulated MR images for three typical functional magnetic resonance imaging (fMRI) sequence types with proton density contrast (e.g., ASL), BOLD contrast, VASO contrast. The biomechanically modeled gyrus is shown below the corresponding image of a digital brain phantom of MRiLab for illustrating the typical whole-brain MR contrasts seen at 7 Tesla and the relative spatial dimensions (red box). (B) Maps of signal change due to tissue displacement following a 2% increase of both volume and stiffness of the activated layer segment (see also Figure 1B) at three different resolutions, i.e., voxel sizes of 1, 0.5, and 0.25 mm.

Signal changes due to tissue displacement were strongest in fMRI methods based on images with high contrast between neighboring tissues and cortical layers, i.e., with effect sizes decreasing from VASO to BOLD and proton density contrast (e.g., ASL). Further, the artifacts due to displacements increased with the image resolution.

For all fMRI variants, the shift of gyral tissue from left to right (Figure 3) resulted in strongest signal changes along the pial surface, i.e., the boundary of cortical tissue to CSF, with a signal increase on the left side and a decrease of signal on its right side. Similar effects of smaller magnitude can be seen along the boundary between the cortex and WM – with inverted sign in the case of VASO fMRI.

Note the typical signal oscillations of the Gibb’s ringing artifacts (Kellner et al., 2016) unfolding orthogonally to the CSF/GM border within gyral tissue and CSF.


Our study addressed the impact of mechanical effects of blood flow/volume changes on the spatial resolution and artifacts in fMRI by mechanical modeling paired with MRI signal models. We generated a mechanical tissue model to estimate the range of mechanical displacement that is expected in cortical layers of the human brain. This model will be briefly discussed in the following paragraph.

Mechanical Simulations

Since there are not yet accurate methods to estimate the exact values of stress within brain tissue in vivo, we used the geometry itself as the input variable. The shape of the tissue is a result of all forces acting on the brain during development and we regard the tissue as a confined geometry subjected to mechanical stresses during development. We estimate these stresses by numerically inverting the developmental process and bringing the geometry to a smooth, unfolded shape. This stage is typically considered as the starting point for established models of cortical folding (Bayly et al., 2013; Tallinen et al., 2016). Our simulation of the initial stress distribution for the cortical gyrus (Figure 2) shows a good agreement with the results presented by other studies [See Figure 7 in Xu et al. (2010)], which use the forward growth model by Rodriguez et al. (1994). Based on this initial stress distribution and assuming a local change of tissue stiffness due to blood flow, our simulations predict a considerable displacement of the tissue in the gyral crown of 4 μm per percent stiffness change. Independently of the initial distribution of stresses in the tissue, our results show that local volume changes also induce tissue displacements. Here, the effect would even be higher leading to a displacement of 7 μm per percent volume change in our simulation. Our study indicates that the combination of the two different effects, i.e., displacements due to local stiffness and volume change are cumulative. Overall, our simulations predict a linear increase in displacement with increasing local stiffness or volume. This sets a lower bound to the expected effect in real experimental settings, since our presented results are based on conservative estimates of parameters as the values are taken from ex-vivo measurements for tissue stiffness via indentation, which underestimate stiffness of perfused tissue. Furthermore, the volume and stiffness changes from literature were likely underestimated due to partial volume effects of low-resolution data.

Functional Magnetic Resonance Imaging Simulations

Our fMRI simulations show that the predicted tissue displacements following focal vasodilatory activity in a gyrus can create erroneous fMRI signal changes on the same order of magnitude as typically reported for classical fMRI contrast mechanisms, i.e., between 0.5 and 5% signal change of the raw images acquired in VASO, BOLD, or ASL fMRI (Huber et al., 2019). Both positive and negative signal changes resulting from these anatomical deformations will interfere with the true fMRI signals of interest that are more directly coupled to neuronal activity reflecting excitatory or inhibitory neural processes. Further, the spatial precision of fMRI may be significantly impacted. The peak of these erroneous signal changes can be located some centimeters away from the origin of the focally activated capillary bed embedding the neuronal activity, i.e., from the base to the crown of the exemplary gyral structure modeled in our simulations.

Importantly, spurious signal changes due to activity-induced tissue displacements are strongest in fMRI methods based on images with high contrast between neighboring tissues, i.e., with effect sizes considerably decreasing from VASO to BOLD and proton density contrast (e.g., ASL). Therefore, the relative fMRI sensitivity of ASL compared to VASO might be empirically underestimated – especially at high-resolution fMRI. Further, it should be noted that fMRI variants nulling signal of specific tissue compartments like VASO or comparable methods (Shen et al., 2009) can result in erroneously high relative signal changes in those voxels (division by zero) and should therefore be interpreted with highest care. For example, cortical layers with T1 values matching T1b will be nulled in VASO fMRI and therefore create singularities in the relative signal change map.

As expected, the displacement fMRI signal change increases with the image resolution due to partial volume effect, i.e., following the ratio of the tissue displacement D and the voxel size V. Generally, the relative signal change due to tissue displacement can be approximated by ΔS/S = C × D/V with the relative tissue contrast C. In the example of idealized proton density contrast (with the lowest artifact strength compared to BOLD or VASO), the maximum tissue contrast C is found on the pial surface (creating a sharp discontinuity of proton density, i.e., 80% for GM and 100% for CSF). Assuming a tissue displacement of only 0.02 mm (∼double of a capillary diameter) resulting from 2% focal volume and stiffness changes (the lower limit of the explored range herein) predicts signal changes of 2, 1, and 0.5% with isotropic voxel sizes of 0.25, 0.5, and 1.0 mm, which is in good agreement with the corresponding signal change maps of our fMRI simulations.

For simplicity we assumed a linear increase of T2 values across cortical layers going from the white interface to the pial surface. Different T2 distributions across layers would lead to slightly different results. A more complex laminar T2 distribution, e.g., the one observed by McColgan et al. (2021) with a rapid drop of T2 in the outermost layer likely due to superficial veins, leads to a contrast increase between CSF and superficial GM. This would further slightly amplify the effects predicted by our simulations, especially in fMRI sequences employing longer echo times like BOLD fMRI.

Generally, resolving cortical layers and detecting the predicted tissue shift artifacts requires very high spatial resolution enabled by higher magnetic field strength and the associated increase of signal-to-noise ratio. However, at a given resolution the tissue shift effect is not expected to depend strongly on field strengths (in the range commonly available in human or animal MRI). Independent of the T1 relaxation of brain tissue at different field strengths, the GM-CSF contrast that dominates the simulated effect will peak for VASO sequences suppressing the signal of GM tissue. In BOLD fMRI, a potential decrease of GM-CSF contrast due to slower T2 relaxation at lower field will generally be compensated by the choice of TE = T2 of GM which provides optimal BOLD fMRI sensitivity. However, the significant decrease of R2 of pial veins at lower field strengths (Gati et al., 1997) will slightly dampen the tissue shift effect for BOLD fMRI at lower field. Independent of field strength, the tissue shift effect is minimal for ASL sequences at short TE due to the relatively small proton density difference between CSF and brain tissue. In line with the presumed robustness of the tissue shift effect across static magnetic field strengths it might contribute to functional signal changes observed in multi-echo BOLD fMRI at very short echo times independently from the main magnetic field strength (Markuerkiaga et al., 2021).

The conservative assumptions used in our biomechanical simulations tend to underestimate activity-induced tissue displacements. Based on our analyses, different choices of the relevant mechanical parameters translate linearly (within the ranges explored) into tissue displacements and the associated erroneous fMRI signal changes, respectively. Representative parameters are the already discussed baseline stiffness of cortical tissue (taken from ex-vivo data) and the GM stiffness increase of 2% taken from hypercapnia MRE experiments acquired with 2.0 mm isotropic voxels (Hetzer et al., 2018b) barely matching the thickness of cortical GM (Fischl and Dale, 2000). Assuming that the stiffness increase following vasodilation is driven by the highly vascularized layers only, i.e., 1/3 of the voxel volume, a three times higher stiffness change in the activated layer segment might be more realistic (partial volume effect) translating into three times higher signal changes due to the increased tissue displacement.

Further, higher erroneous displacement fMRI signal changes can be expected following our observation that the displacement scales with the activated volume. Our choice of the thin layer 1 to be activated in our fMRI simulation falls within the lower limit of displacements of gyral tissue following neuronal stimulation. For example, we found a factor of 3 higher displacement when activating one central layer instead of layer 1. As neuronal activity is usually distributed across several layers in most fMRI paradigms, the activated tissue volume will probably be even higher while the highly vascularized and perfused central layers are expected to dominate local tissue displacement effects (Shen et al., 2008; Jung et al., 2021). We found that the relative displacement effect (normalized by the volume of the activated layer segment) slightly decreases from outer to inner layers. It should be noted that the influence of layer thickness variability across cortical GM on the tissue shift effect dominates the influence of the baseline stiffness variability (see also Figure 6) of roughly 10% across layers (Iwashita et al., 2014).

Furthermore, depending on the gyral geometry and the position of the activated layer segment therein, the resulting forces can also be effectively absorbed by the surrounding tissue without causing significant deformations. Such favorable neuroanatomical configurations, i.e., in functional areas with low gyrification or flat and thin gyri (Alemán-Gómez et al., 2013), might constitute a precondition for current state-of-the-art fMRI studies successfully showing high layer-specificity (Huber et al., 2018; Yu et al., 2019). Furthermore, with current layer fMRI techniques, it cannot be excluded that other noise sources (i.e., motion resulting from breathing and cardiac pulsations, or manual steps during layer segmentation) dominate and obscure the small functional tissue shift effects predicted in our work, especially as motion correction or pseudo-automatic layer segmentation are commonly focused to specific areas of interest (Huber et al., 2020).

Assuming that practical solutions are found for removing the harmonic non-rigid motion of brain tissue during the cardiac cycle (Terem et al., 2018; Adams et al., 2019), e.g., by triggering the MR acquisition, our simulations predict that future technical developments of high-resolution fMRI might have to consider local cortical deformations of tissue mechanically interconnected to the activated area. For example, automated dynamic layer segmentation methods capable of dissociating local deformations from intensity changes within activated layers might be needed to correct for local tissue shifts in order to detect fMRI signal changes more directly coupled to neuronal activity. Further, our simulations suggest that such approaches should also actively address Gibb’s ringing (Kellner et al., 2016) (creating alternating signal orthogonally to sharp image boundaries similar to layers) which will be amplified by the subtractive nature of fMRI in the context of tissue displacements.

We are not aware of any results in published layer fMRI studies directly corroborating the predicted gyral deformation effects, but this may be explained by them focusing on rather small regions of interest and applying artifact suppression methods. To address this question, a systematic meta-analysis of mesoscopic fMRI data may be conducted to detect potential inverted signal changes on the opposite sites of activated gyri as predicted by our simulations.

As tissue displacement effects are expected to unfold within the vasodilatory reaction time of some seconds following neuronal activity, fMRI with high spatiotemporal resolution will gain from addressing tissue displacement effects to better dissociate early neuronal-driven responses (i.e., rapid changes of blood relaxation time reflecting oxygen extraction) and the delayed vascular-driven responses (Jung et al., 2021). Models of the laminar BOLD response (Havlicek and Uludaǧ, 2020; Uludag and Havlicek, 2021) incorporating the effects of gyral displacements may not only sharpen the spatial precision of laminar fMRI (Fracasso et al., 2021) but also help interpreting fMRI signal changes associated with non-vascular mechanisms (Le Bihan et al., 2006; Roth and Basser, 2009; Patz et al., 2019).


Our simulation study delivers new insights into blood flow-induced local displacement of human brain tissue in vivo. However, there are some limitations: In our current proof-of-principle study, we restrict our simulation to a 2-dimensional model of one representative gyrus. Further research implementing a computationally demanding full 3D model of the brain with different activated gyri may shed further light on how our predicted effect sizes depend on the anatomical variability. The deformative effects predicted by our quasi-static model might be dampened or temporally shifted relative to the hemodynamic response depending on the unknown underlying CSF redistribution mechanism (Jin and Kim, 2010). As long as the physical origin and the laminar extension of the observed volume fraction change due to neuronal activity is not resolved, a computationally highly demanding dynamic viscoelastic model of the brain would have to include hydrodynamics and a high number of unknown parameters, i.e., friction, pressure differentials, physical volume exchange rates within tissue, or the dynamic vessel dilation process across layers, etc. The maps used by Jin et al. acquired at mesoscopic resolution would be a perfect candidate for future empirical studies exploring the sources of CSF redistributions caused by neuronal activity. Furthermore, the baseline tissue stiffness values were taken from the ex-vivo studies using indentation methods that might not exactly reflect in vivo values because brain properties change post-mortem (Bertalan et al., 2020). However, current in vivo MRE is limited in detecting cortical stiffness due to boundary effects on propagating shear waves (Hirsch et al., 2017). Moreover, we have neglected the impact of biomechanics of CSF in our simulation as well as the pia mater due to its thinness which is a few micrometers and much thinner than a cortical layer. Including the mechanical properties of CSF might slightly change the displacement field of cortical tissue. More interestingly, the CSF redistribution flow patterns resulting from gyral displacements could induce flow artifacts in fMRI signals (Duyn et al., 1994; Fultz et al., 2019) acquired near the pial surface where the highest CSF flow velocities are expected. In a future study we will include the hydrodynamics of this CSF redistribution in our simulations to analyze associated artifactual signal changes in high-resolution fMRI.

The study focused on biomechanics and did not address the myriad factors impacting fMRI measurements. Established BOLD-fMRI models show the importance of the interaction of blood flow, volume and oxygenation. They are increasingly refined to include the intricacies of the vascular system of the cortex including draining vein and orientation effects (Viessmann et al., 2019). Also the regulation of the blood flow remains an intense area of research on coupling of blood flow and different types of neuroelectric activity (Logothetis et al., 2001; Havlicek and Uludaǧ, 2020).


In this biomechanical simulation study, we predict that increasing the resolution of fMRI to resolve single cortical layers is hindered by undesired small local deformations of cortical tissue mechanically interconnected to the activated area. Following the vasodilatory response, volume and stiffness changes of the activated layer lead to displacements realistically varying from 10 to 100 μm in gyral structures. The artifacts and erroneous changes in the predicted fMRI signal reached up to 5% – a similar magnitude as typical neurovascular responses. The artifactual signal changes increase with decreasing voxel size and are strongest in fMRI methods based on images with high contrast between neighboring tissues, i.e., with effect sizes considerably decreasing from VASO to BOLD and ASL fMRI. Although our simulations require experimental validation and would benefit from biomechanical properties measured in vivo, they suggest that focal biomechanical changes may relevantly affect ultra-high resolution fMRI experiments and deserve further consideration.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

Author Contributions

MZ and SHe contributed to the conception and design of the study, performed the simulations and the statistical data analysis, and wrote the first draft of the manuscript. MZ, SHe, NS, CJ, IS, SHi, and NW contributed substantially to interpretation of the data, revised the manuscript critically for intellectual content and have approved the submitted version.


The research leading to these results has received funding from German Research Foundation (Exc 257 #276880906 and #39052203) and the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007-2013) / ERC grant agreement no. 616905. NW has received funding from the European Union’s Horizon 2020 Research and Innovation Programme under the grant agreement no. 681094 and the BMBF (01EW1711A&B) in the framework of ERA-NET NEURON.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary Material

The Supplementary Material for this article can be found online at:


  1. ^


Adams, A. L., Kuijf, H. J., Viergever, M. A., Luijten, P. R., and Zwanenburg, J. J. M. (2019). Quantifying cardiac-induced brain tissue expansion using DENSE. NMR Biomed. 32:e4050. doi: 10.1002/nbm.4050

PubMed Abstract | CrossRef Full Text | Google Scholar

Alemán-Gómez, Y., Janssen, J., Schnack, H., Balaban, E., Pina-Camacho, L., Alfaro-Almagro, F., et al. (2013). The human cerebral cortex flattens during adolescence. J. Neurosci. 33, 15004–15010. doi: 10.1523/JNEUROSCI.1459-13.2013

PubMed Abstract | CrossRef Full Text | Google Scholar

Amunts, K., Lepage, C., Borgeat, L., Mohlberg, H., Dickscheid, T., Rousseau, M. -É, et al. (2013). BigBrain: an ultrahigh-resolution 3d human brain model. Science 340, 1472–1475. doi: 10.1126/science.1235381

PubMed Abstract | CrossRef Full Text | Google Scholar

Armstrong, E., Schleicher, A., Omran, H., Curtis, M., and Zilles, K. (1995). The ontogeny of human gyrification. Cereb. Cortex 5, 56–63.

Google Scholar

Bayly, P. V., Taber, L. A., and Kroenke, C. D. (2014). Mechanical forces in cerebral cortical folding: a review of measurements and models. J. Mech. Behav. Biomed. Mater. 29, 568–581. doi: 10.1016/j.jmbbm.2013.02.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Bayly, P., Okamoto, R., Xu, G., Shi, Y., and Taber, L. (2013). A cortical folding model incorporating stress-dependent growth explains gyral wavelengths and stress patterns in the developing brain. Phys. Biol. 10:016005.

Google Scholar

Bertalan, G., Klein, C., Schreyer, S., Steiner, B., Kreft, B., Tzschätzsch, H., et al. (2020). Biomechanical properties of the hypoxic and dying brain quantified by magnetic resonance elastography. Acta Biomater. 101, 395–402. doi: 10.1016/j.actbio.2019.11.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Brain Initiative (2021). Dissemination of Non-Invasive Imaging Technologies Workshop (David Feinberg). Available online at: (accessed February 23, 2021).

Google Scholar

Budday, S., Nay, R., de Rooij, R., Steinmann, P., Wyrobek, T., Ovaert, T. C., et al. (2015). Mechanical properties of gray and white matter brain tissue by indentation. J. Mech. Behav. Biomed. Mater. 46, 318–330. doi: 10.1016/j.jmbbm.2015.02.024

PubMed Abstract | CrossRef Full Text | Google Scholar

Budday, S., Ovaert, T. C., Holzapfel, G. A., Steinmann, P., and Kuhl, E. (2020). Fifty shades of brain: a review on the mechanical testing and modeling of brain tissue. Arch. Comput. Methods Eng. 27, 1187–1230. doi: 10.1007/s11831-019-09352-w

CrossRef Full Text | Google Scholar

Budday, S., Sommer, G., Birkl, C., Langkammer, C., Haybaeck, J., Kohnert, J., et al. (2016). Mechanical characterization of human brain tissue. Acta Biomater. 48, 319–340. doi: 10.1016/j.actbio.2016.10.036

PubMed Abstract | CrossRef Full Text | Google Scholar

Budday, S., Steinmann, P., and Kuhl, E. (2014). The role of mechanics during brain development. J. Mech. Phys. Solids 72, 75–92. doi: 10.1016/j.jmps.2014.07.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Caan, M. W. A., Bazin, P.-L., Marques, J. P., de Hollander, G., Dumoulin, S. O., Zwaag, W., et al. (2019). MP2RAGEME: T1, T2, and QSM mapping in one sequence at 7 tesla. Hum. Brain Mapp. 40, 1786–1798. doi: 10.1002/hbm.24490

PubMed Abstract | CrossRef Full Text | Google Scholar

Carey, D., Caprini, F., Allen, M., Lutti, A., Weiskopf, N., Rees, G., et al. (2018). Quantitative MRI provides markers of intra-, inter-regional, and age-related differences in young adult cortical microstructure. Neuroimage 182, 429–440. doi: 10.1016/j.neuroimage.2017.11.066

PubMed Abstract | CrossRef Full Text | Google Scholar

Deichmann, R., Josephs, O., Hutton, C., Corfield, D. R., and Turner, R. (2002). Compensation of susceptibility-induced BOLD sensitivity losses in echo-planar fMRI imaging. NeuroImage 15, 120–135. doi: 10.1006/nimg.2001.0985

PubMed Abstract | CrossRef Full Text | Google Scholar

Draganski, B., Gaser, C., Busch, V., Schuierer, G., Bogdahn, U., and May, A. (2004). Neuroplasticity: changes in grey matter induced by training. Nature 427, 311–312. doi: 10.1038/427311a

PubMed Abstract | CrossRef Full Text | Google Scholar

Duyn, J. H., Moonen, C. T., van Yperen, G. H., de Boer, R. W., and Luyten, P. R. (1994). Inflow versus deoxyhemoglobin effects in BOLD functional MRI using gradient echoes at 1.5 T. NMR Biomed. 7, 83–88. doi: 10.1002/nbm.1940070113

PubMed Abstract | CrossRef Full Text | Google Scholar

Fåhræus, R., and Lindqvist, T. (1931). The viscosity of the blood in narrow capillary tubes. Am. J. Physiol. Leg. Content 96, 562–568. doi: 10.1152/ajplegacy.1931.96.3.562

CrossRef Full Text | Google Scholar

Fischl, B., and Dale, A. M. (2000). Measuring the thickness of the human cerebral cortex from magnetic resonance images. Proc. Natl. Acad. Sci. U.S.A. 97, 11050–11055. doi: 10.1073/pnas.200033797

PubMed Abstract | CrossRef Full Text | Google Scholar

Fracasso, A., Dumoulin, S. O., and Petridou, N. (2021). Point-spread function of the BOLD response across columns and cortical depth in human extra-striate cortex. Prog. Neurobiol. 202:102034. doi: 10.1016/j.pneurobio.2021.102034

PubMed Abstract | CrossRef Full Text | Google Scholar

Fracasso, A., Luijten, P. R., Dumoulin, S. O., and Petridou, N. (2018). Laminar imaging of positive and negative BOLD in human visual cortex at 7T. NeuroImage 164, 100–111. doi: 10.1016/j.neuroimage.2017.02.038

PubMed Abstract | CrossRef Full Text | Google Scholar

Franceschini, G., Bigoni, D., Regitnig, P., and Holzapfel, G. A. (2006). Brain tissue deforms similarly to filled elastomers and follows consolidation theory. J. Mech. Phys. Solids 54, 2592–2620. doi: 10.1016/j.jmps.2006.05.004

CrossRef Full Text | Google Scholar

Fukunaga, M., Li, T.-Q., van Gelderen, P., de Zwart, J. A., Shmueli, K., Yao, B., et al. (2010). Layer-specific variation of iron content in cerebral cortex as a source of MRI contrast. Proc. Natl. Acad. Sci. U.S.A. 107, 3834–3839. doi: 10.1073/pnas.0911177107

PubMed Abstract | CrossRef Full Text | Google Scholar

Fultz, N. E., Bonmassar, G., Setsompop, K., Stickgold, R. A., Rosen, B. R., Polimeni, J. R., et al. (2019). Coupled electrophysiological, hemodynamic, and cerebrospinal fluid oscillations in human sleep. Science 366, 628–631. doi: 10.1126/science.aax5440

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganpule, S., Daphalapurkar, N. P., Cetingul, M. P., and Ramesh, K. T. (2018). Effect of bulk modulus on deformation of the brain under rotational accelerations. Shock Waves 28, 127–139.

Google Scholar

Garcia, K. E., Kroenke, C. D., and Bayly, P. V. (2018a). Mechanics of cortical folding: stress, growth and stability. Philos. Trans. R. Soc. Lond. B Biol. Sci. 373:20170321.

Google Scholar

Garcia, K. E., Robinson, E. C., Alexopoulos, D., Dierker, D. L., Glasser, M. F., Coalson, T. S., et al. (2018b). Dynamic patterns of cortical expansion during folding of the preterm human brain. Proc. Natl. Acad. Sci. U.S.A. 115, 3156–3161. doi: 10.1073/pnas.1715451115

PubMed Abstract | CrossRef Full Text | Google Scholar

Gati, J. S., Menon, R. S., Ugurbil, K., and Rutt, B. K. (1997). Experimental determination of the BOLD field strength dependence in vessels and tissue. Magn. Reson. Med. 38, 296–302. doi: 10.1002/mrm.1910380220

PubMed Abstract | CrossRef Full Text | Google Scholar

Goense, J., Bohraus, Y., and Logothetis, N. K. (2016). fMRI at high spatial resolution: implications for BOLD-models. Front. Comput. Neurosci. 10:66. doi: 10.3389/fncom.2016.00066

PubMed Abstract | CrossRef Full Text | Google Scholar

Hardan, A. Y., Jou, R. J., Keshavan, M. S., Varma, R., and Minshew, N. J. (2004). Increased frontal cortical folding in autism: a preliminary MRI study. Psychiatry Res. 131, 263–268. doi: 10.1016/j.pscychresns.2004.06.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Havlicek, M., and Uludaǧ, K. (2020). A dynamical model of the laminar BOLD response. NeuroImage 204:116209. doi: 10.1016/j.neuroimage.2019.116209

PubMed Abstract | CrossRef Full Text | Google Scholar

Hetzer, S., Bormann, K., Dittmann, F., Hirsch, S., Braun, J., and Sack, I. (2018a). “Hypercapnia-induced vasodilatation increases brain stiffness,” in Proceedings of the International Society for Magnetic Resonance in Medicine (ISMRM), Paris. Talk #705.

Google Scholar

Hetzer, S., Dittmann, F., Bormann, K., Hirsch, S., Lipp, A., Wang, D. J. J., et al. (2018b). Hypercapnia increases brain viscoelasticity. J. Cereb. Blood Flow Metab. 39, 2445–2455. doi: 10.1177/0271678X18799241

PubMed Abstract | CrossRef Full Text | Google Scholar

Hetzer, S., Mildner, T., and Möller, H. E. (2011). A modified EPI sequence for high-resolution imaging at ultra-short echo time. Magn. Reson. Med. 65, 165–175. doi: 10.1002/mrm.22610

PubMed Abstract | CrossRef Full Text | Google Scholar

Hirsch, S., Braun, J., and Sack, I. (2017). Magnetic Resonance Elastography: Physical Background and Medical Applications. Weinheim: Wiley-VCH Verlag GmbH & Co. KGaA.

Google Scholar

Holzapfel, G. A. (2000). Nonlinear Solid Mechanics: A Continuum Approach for Engineering | Wiley. Available online at: (accessed February 2, 2021).

Google Scholar

Hua, J., Qin, Q., Pekar, J. J., and van Zijl, P. C. M. (2011). Measurement of absolute arterial cerebral blood volume in human brain without using a contrast agent. NMR Biomed. 24, 1313–1325. doi: 10.1002/nbm.1693

PubMed Abstract | CrossRef Full Text | Google Scholar

Huber, L., Finn, E. S., Chai, Y., Goebel, R., Stirnberg, R., Stöcker, T., et al. (2020). Layer-dependent functional connectivity methods. Prog. Neurobiol. 101835. doi: 10.1016/j.pneurobio.2020.101835

PubMed Abstract | CrossRef Full Text | Google Scholar

Huber, L., Handwerker, D. A., Jangraw, D. C., Chen, G., Hall, A., Stüber, C., et al. (2017). High-Resolution CBV-fMRI allows mapping of laminar activity and connectivity of cortical input and output in human M1. Neuron 96, 1253–1263.e7. doi: 10.1016/j.neuron.2017.11.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Huber, L., Ivanov, D., Handwerker, D. A., Marrett, S., Guidi, M., Uludað, K., et al. (2018). Techniques for blood volume fMRI with VASO: from low-resolution mapping towards sub-millimeter layer-dependent applications. NeuroImage 164, 131–143. doi: 10.1016/j.neuroimage.2016.11.039

PubMed Abstract | CrossRef Full Text | Google Scholar

Huber, L., Uludað, K., and Möller, H. E. (2019). Non-BOLD contrast for laminar fMRI in humans: CBF, CBV, and CMRO2. NeuroImage 197, 742–760. doi: 10.1016/j.neuroimage.2017.07.041

PubMed Abstract | CrossRef Full Text | Google Scholar

Huesmann, G. R., Schwarb, H., Smith, D. R., Pohlig, R. T., Anderson, A. T., McGarry, M. D. J., et al. (2020). Hippocampal stiffness in mesial temporal lobe epilepsy measured with MR elastography: preliminary comparison with healthy participants. NeuroImage Clin. 27:102313. doi: 10.1016/j.nicl.2020.102313

PubMed Abstract | CrossRef Full Text | Google Scholar

Iwashita, M., Kataoka, N., Toida, K., and Kosodo, Y. (2014). Systematic profiling of spatiotemporal tissue and cellular stiffness in the developing brain. Dev. Camb. Engl. 141, 3793–3798. doi: 10.1242/dev.109637

PubMed Abstract | CrossRef Full Text | Google Scholar

Jin, T., and Kim, S.-G. (2010). Change of the cerebrospinal fluid volume during brain activation investigated by T1ρ-weighted fMRI. NeuroImage 51, 1378–1383. doi: 10.1016/j.neuroimage.2010.03.047

PubMed Abstract | CrossRef Full Text | Google Scholar

Jung, W. B., Im, G. H., Jiang, H., and Kim, S.-G. (2021). Early fMRI responses to somatosensory and optogenetic stimulation reflect neural information flow. Proc. Natl. Acad. Sci. U.S.A. 118:e2023265118. doi: 10.1073/pnas.2023265118

PubMed Abstract | CrossRef Full Text | Google Scholar

Kellner, E., Dhital, B., Kiselev, V. G., and Reisert, M. (2016). Gibbs-ringing artifact removal based on local subvoxel-shifts. Magn. Reson. Med. 76, 1574–1581. doi: 10.1002/mrm.26054

PubMed Abstract | CrossRef Full Text | Google Scholar

Koser, D. E., Thompson, A. J., Foster, S. K., Dwivedy, A., Pillai, E. K., Sheridan, G. K., et al. (2016). Mechanosensing is critical for axon growth in the developing brain. Nat. Neurosci. 19, 1592–1598. doi: 10.1038/nn.4394

PubMed Abstract | CrossRef Full Text | Google Scholar

Laffey, J. G., and Kavanagh, B. P. (2002). Hypocapnia. N. Engl. J. Med. 347, 43–53. doi: 10.1056/NEJMra012457

PubMed Abstract | CrossRef Full Text | Google Scholar

Le Bihan, D., Urayama, S., Aso, T., Hanakawa, T., and Fukuyama, H. (2006). Direct and fast detection of neuronal activation in the human brain with diffusion MRI. Proc. Natl. Acad. Sci. U.S.A. 103, 8263–8268. doi: 10.1073/pnas.0600644103

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, G., Wang, L., Shi, F., Lyall, A. E., Ahn, M., Peng, Z., et al. (2016). Cortical thickness and surface area in neonates at high risk for schizophrenia. Brain Struct. Funct. 221, 447–461. doi: 10.1007/s00429-014-0917-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Lin, J. J., Salamon, N., Lee, A. D., Dutton, R. A., Geaga, J. A., Hayashi, K. M., et al. (2007). Reduced neocortical thickness and complexity mapped in mesial temporal lobe epilepsy with hippocampal sclerosis. Cereb. Cortex 17, 2007–2018. doi: 10.1093/cercor/bhl109

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, F., Velikina, J. V., Block, W. F., Kijowski, R., and Samsonov, A. A. (2017). Fast realistic MRI simulations based on generalized multi-pool exchange tissue model. IEEE Trans. Med. Imaging 36, 527–537. doi: 10.1109/TMI.2016.2620961

PubMed Abstract | CrossRef Full Text | Google Scholar

Logothetis, N. K., Pauls, J., Augath, M., Trinath, T., and Oeltermann, A. (2001). Neurophysiological investigation of the basis of the fMRI signal. Nature 412, 150–157. doi: 10.1038/35084005

PubMed Abstract | CrossRef Full Text | Google Scholar

Lu, H., Golay, X., Pekar, J. J., and Van Zijl, P. C. M. (2003). Functional magnetic resonance imaging based on changes in vascular space occupancy. Magn. Reson. Med. 50, 263–274. doi: 10.1002/mrm.10519

PubMed Abstract | CrossRef Full Text | Google Scholar

Maguire, E. A., Gadian, D. G., Johnsrude, I. S., Good, C. D., Ashburner, J., Frackowiak, R. S., et al. (2000). Navigation-related structural change in the hippocampi of taxi drivers. Proc. Natl. Acad. Sci. U.S.A. 97, 4398–4403. doi: 10.1073/pnas.070039597

PubMed Abstract | CrossRef Full Text | Google Scholar

Markuerkiaga, I., Marques, J. P., Bains, L. J., and Norris, D. G. (2021). An in-vivo study of BOLD laminar responses as a function of echo time and static magnetic field strength. Sci. Rep. 11:1862. doi: 10.1038/s41598-021-81249-w

PubMed Abstract | CrossRef Full Text | Google Scholar

McColgan, P., Helbling, S., Vaculèiaková, L., Pine, K., Wagstyl, K., Attar, F. M., et al. (2021). Relating quantitative 7T MRI across cortical depths to cytoarchitectonics, gene expression and connectomics. Hum. Brain Mapp. doi: 10.1002/hbm.25595

PubMed Abstract | CrossRef Full Text | Google Scholar

Mcintosh, R. L., and Anderson, V. (2010). A comprehensive tissue properties database provided for the thermal assessment of a human at rest. Biophys. Rev. Lett. 05, 129–151. doi: 10.1142/S1793048010001184

CrossRef Full Text | Google Scholar

Murphy, M. C., Huston, J., Jack, C. R., Glaser, K. J., Manduca, A., Felmlee, J. P., et al. (2011). Decreased brain stiffness in Alzheimer’s disease determined by magnetic resonance elastography. J. Magn. Reson. Imaging JMRI 34, 494–498. doi: 10.1002/jmri.22707

PubMed Abstract | CrossRef Full Text | Google Scholar

Muthupillai, R., Lomas, D. J., Rossman, P. J., Greenleaf, J. F., Manduca, A., and Ehman, R. L. (1995). Magnetic resonance elastography by direct visualization of propagating acoustic strain waves. Science 269, 1854–1857. doi: 10.1126/science.7569924

PubMed Abstract | CrossRef Full Text | Google Scholar

Oliveri, H., Franze, K., and Goriely, A. (2021). A theory for durotactic axon guidance. Phys. Rev. Lett. 126:118101. doi: 10.1103/PhysRevLett.126.118101

PubMed Abstract | CrossRef Full Text | Google Scholar

Parker, K. J. (2017). Are rapid changes in brain elasticity possible? Phys. Med. Biol. 62, 7425–7439. doi: 10.1088/1361-6560/aa8380

PubMed Abstract | CrossRef Full Text | Google Scholar

Patz, S., Fovargue, D., Schregel, K., Nazari, N., Palotai, M., Barbone, P. E., et al. (2019). Imaging localized neuronal activity at fast time scales through biomechanics. Sci. Adv. 5:eaav3816. doi: 10.1126/sciadv.aav3816

PubMed Abstract | CrossRef Full Text | Google Scholar

Piechnik, S. K., Evans, J., Bary, L. H., Wise, R. G., and Jezzard, P. (2009). Functional changes in CSF volume estimated using measurement of water T2 relaxation. Magn. Reson. Med. 61, 579–586. doi: 10.1002/mrm.21897

PubMed Abstract | CrossRef Full Text | Google Scholar

Polimeni, J. R., Renvall, V., Zaretskaya, N., and Fischl, B. (2018). Analysis strategies for high-resolution UHF-fMRI data. NeuroImage 168, 296–320. doi: 10.1016/j.neuroimage.2017.04.053

PubMed Abstract | CrossRef Full Text | Google Scholar

Rodriguez, E. K., Hoger, A., and McCulloch, A. D. (1994). Stress-dependent finite growth in soft elastic tissues. J. Biomech. 27, 455–467. doi: 10.1016/0021-9290(94)90021-3

CrossRef Full Text | Google Scholar

Rooney, W. D., Johnson, G., Li, X., Cohen, E. R., Kim, S.-G., Ugurbil, K., et al. (2007). Magnetic field and tissue dependencies of human brain longitudinal 1H2O relaxation in vivo. Magn. Reson. Med. 57, 308–318. doi: 10.1002/mrm.21122

PubMed Abstract | CrossRef Full Text | Google Scholar

Roth, B. J., and Basser, P. J. (2009). Mechanical model of neural tissue displacement during Lorentz effect imaging. Magn. Reson. Med. 61, 59–64. doi: 10.1002/mrm.21772

PubMed Abstract | CrossRef Full Text | Google Scholar

Shen, Q., Ren, H., and Duong, T. Q. (2008). CBF, BOLD, CBV, and CMRO(2) fMRI signal temporal dynamics at 500-msec resolution. J. Magn. Reson. Imaging JMRI 27, 599–606. doi: 10.1002/jmri.21203

PubMed Abstract | CrossRef Full Text | Google Scholar

Shen, Y., Kauppinen, R. A., Vidyasagar, R., and Golay, X. (2009). A functional magnetic resonance imaging technique based on nulling extravascular gray matter signal. J. Cereb. Blood Flow Metab. Off. J. Int. Soc. Cereb. Blood Flow Metab. 29, 144–156. doi: 10.1038/jcbfm.2008.96

PubMed Abstract | CrossRef Full Text | Google Scholar

Streitberger, K.-J., Lilaj, L., Schrank, F., Braun, J., Hoffmann, K.-T., Reiss-Zimmermann, M., et al. (2020). How tissue fluidity influences brain tumor progression. Proc. Natl. Acad. Sci. U.S.A. 117, 128–134. doi: 10.1073/pnas.1913511116

PubMed Abstract | CrossRef Full Text | Google Scholar

Streitberger, K.-J., Wiener, E., Hoffmann, J., Freimann, F. B., Klatt, D., Braun, J., et al. (2011). In vivo viscoelastic properties of the brain in normal pressure hydrocephalus. NMR Biomed. 24, 385–392. doi: 10.1002/nbm.1602

PubMed Abstract | CrossRef Full Text | Google Scholar

Tallinen, T., Chung, J. Y., Rousseau, F., Girard, N., Lefèvre, J., and Mahadevan, L. (2016). On the growth and form of cortical convolutions. Nat. Phys. 12, 588–593. doi: 10.1038/nphys3632

CrossRef Full Text | Google Scholar

Terem, I., Ni, W. W., Goubran, M., Rahimi, M. S., Zaharchuk, G., Yeom, K. W., et al. (2018). Revealing sub-voxel motions of brain tissue using phase-based amplified MRI (aMRI). Magn. Reson. Med. 80, 2549–2559. doi: 10.1002/mrm.27236

PubMed Abstract | CrossRef Full Text | Google Scholar

Tofts, P. (2010). Quantitative MRI of the Brain: Measuring Changes Caused by Disease. Hoboken, N.J: Wiley.

Google Scholar

Trampel, R., Bazin, P.-L., Pine, K., and Weiskopf, N. (2019). In-vivo magnetic resonance imaging (MRI) of laminae in the human cortex. NeuroImage 197, 707–715. doi: 10.1016/j.neuroimage.2017.09.037

PubMed Abstract | CrossRef Full Text | Google Scholar

Uludag, K., and Havlicek, M. (2021). Determining laminar neuronal activity from BOLD fMRI using a generative model. Prog. Neurobiol. 102055. doi: 10.1016/j.pneurobio.2021.102055

PubMed Abstract | CrossRef Full Text | Google Scholar

Viessmann, O., Scheffler, K., Bianciardi, M., Wald, L. L., and Polimeni, J. R. (2019). Dependence of resting-state fMRI fluctuation amplitudes on cerebral cortical orientation relative to the direction of B0 and anatomical axes. NeuroImage 196, 337–350. doi: 10.1016/j.neuroimage.2019.04.036

PubMed Abstract | CrossRef Full Text | Google Scholar

Wagstyl, K., Larocque, S., Cucurull, G., Lepage, C., Cohen, J. P., Bludau, S., et al. (2020). BigBrain 3D atlas of cortical layers: cortical and laminar thickness gradients diverge in sensory and motor cortices. PLoS Biol. 18:e3000678. doi: 10.1371/journal.pbio.3000678

PubMed Abstract | CrossRef Full Text | Google Scholar

Watanabe, T., Sasaki, Y., Shibata, K., and Kawato, M. (2017). Advances in fMRI real-time neurofeedback. Trends Cogn. Sci. 21, 997–1010. doi: 10.1016/j.tics.2017.09.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Weiskopf, N., Edwards, L., Helms, G., Mohammadi, S., and Kirilina, E. (2021). Quantitative magnetic resonance imaging of brain anatomy and in-vivo histology. Nat. Rev. Phys. 3, 570–588.

Google Scholar

Wuerfel, J., Paul, F., Beierbach, B., Hamhaber, U., Klatt, D., Papazoglou, S., et al. (2010). MR-elastography reveals degradation of tissue integrity in multiple sclerosis. NeuroImage 49, 2520–2525. doi: 10.1016/j.neuroimage.2009.06.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, G., Bayly, P. V., and Taber, L. A. (2009). Residual stress in the adult mouse brain. Biomech. Model. Mechanobiol. 8, 253–262. doi: 10.1007/s10237-008-0131-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, G., Knutsen, A. K., Dikranian, K., Kroenke, C. D., Bayly, P. V., and Taber, L. A. (2010). Axons pull on the brain, but tension does not drive cortical folding. J. Biomech. Eng. 132:071013. doi: 10.1115/1.4001683

CrossRef Full Text | Google Scholar

Yin, Z., Romano, A. J., Manduca, A., Ehman, R. L., and Huston, J. (2018). Stiffness and beyond: what MR elastography can tell us about brain structure and function under physiologic and pathologic conditions. Top. Magn. Reson. Imaging TMRI 27, 305–318. doi: 10.1097/RMR.0000000000000178

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, Y., Huber, L., Yang, J., Jangraw, D. C., Handwerker, D. A., Molfese, P. J., et al. (2019). Layer-specific activation of sensory input and predictive feedback in the human primary somatosensory cortex. Sci. Adv. 5:eaav9053. doi: 10.1126/sciadv.aav9053

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, X., Petersen, E. T., Ghariq, E., De Vis, J. B., Webb, A. G., Teeuwisse, W. M., et al. (2013). In vivo blood T(1) measurements at 1.5 T, 3 T, and 7 T. Magn. Reson. Med. 70, 1082–1086. doi: 10.1002/mrm.24550

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: biophysical modeling, simulation, tissue mechanics, blood flow, deformation, BOLD, VASO, fMRI

Citation: Zoraghi M, Scherf N, Jaeger C, Sack I, Hirsch S, Hetzer S and Weiskopf N (2021) Simulating Local Deformations in the Human Cortex Due to Blood Flow-Induced Changes in Mechanical Tissue Properties: Impact on Functional Magnetic Resonance Imaging. Front. Neurosci. 15:722366. doi: 10.3389/fnins.2021.722366

Received: 08 June 2021; Accepted: 23 August 2021;
Published: 21 September 2021.

Edited by:

Mark Wagshul, Albert Einstein College of Medicine, United States

Reviewed by:

Tao Jin, University of Pittsburgh, United States
Martin J. McKeown, University of British Columbia, Canada
Silvia Budday, University of Erlangen Nuremberg, Germany

Copyright © 2021 Zoraghi, Scherf, Jaeger, Sack, Hirsch, Hetzer and Weiskopf. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Mahsa Zoraghi,

These authors share senior authorship