Full length article| Volume 15, ISSUE 3, P654-663, May 01, 2022

# The effect of meninges on the electric fields in TES and TMS. Numerical modeling with adaptive mesh refinement

Open AccessPublished:April 17, 2022

## Highlights

• The impact of the meningeal layers on the electric field in TES and TMS is investigated.
• The boundary element fast multipole method with adaptive mesh refinement is used.
• TES fields in the cortex increase by 30% when the meninges are included.
• TMS fields decrease by 5% on average when the meninges are included.
• The solver, the mesh generator, and an example data set is made publicly available.

## Abstract

### Background

When modeling transcranial electrical stimulation (TES) and transcranial magnetic stimulation (TMS) in the brain, the meninges – dura, arachnoid, and pia mater – are often neglected due to high computational costs.

### Objective

We investigate the impact of the meningeal layers on the cortical electric field in TES and TMS while considering the headreco segmentation as the base model.

### Method

We use T1/T2 MRI data from 16 subjects and apply the boundary element fast multipole method with adaptive mesh refinement, which enables us to accurately solve this problem and establish method convergence at reasonable computational cost. We compare electric fields in the presence and absence of various meninges for two brain areas ($M1HAND$ and $DLPFC$) and for several distinct TES and TMS setups.

### Results

Maximum electric fields in the cortex for focal TES consistently increase by approximately 30% on average when the meninges are present in the CSF volume. Their effect on the maximum field can be emulated by reducing the CSF conductivity from 1.65 S/m to approximately 0.85 S/m. In stark contrast to that, the TMS electric fields in the cortex are only weakly affected by the meningeal layers and slightly (∼6%) decrease on average when the meninges are included.

### Conclusion

Our results quantify the influence of the meninges on the cortical TES and TMS electric fields. Both focal TES and TMS results are very consistent. The focal TES results are also in a good agreement with a prior relevant study. The solver and the mesh generator for the meningeal layers (compatible with SimNIBS) are available online.

## 1. Introduction

Mapping studies that relate electric field strengths to behavioral effects are particularly dependent on modeling accuracy [
• Weise K.
• Numssen O.
• Thielscher A.
• Hartwigsen G.
• Knösche T.R.
A novel approach to localize cortical TMS effects.
,
• Numssen O.
• Zier A.L.
• Thielscher A.
• Hartwigsen G.
• Knösche T.R.
• Weise K.
Efficient high-resolution TMS mapping of the human motor cortex by nonlinear regression.
]. Such studies also contribute to the dosing problem with individual head models [
• Peterchev A.V.
• Wagner T.A.
• Miranda P.C.
• Nitsche M.A.
• Paulus W.
• Lisanby S.H.
• Pascual-Leone A.
• Bikson M.
Fundamentals of transcranial electric and magnetic stimulation dose: definition, selection, and reporting practices.
]. Including multiple tissue classes into the modeling is beneficial for accurate electric field simulations [
• Puonti O.
• Van Leemput K.
• Saturnino G.B.
• Siebner H.R.
• Thielscher A.
Accurate and robust whole-head segmentation from magnetic resonance images for individualized head modeling.
]. Current head modeling pipelines – mri2mesh [
• Windhoff M.
• Opitz A.
• Thielscher A.
Electric field calculations in brain stimulation based on finite elements: an optimized processing pipeline for the generation and usage of accurate individual head models.
• Puonti O.
• Iglesias J.E.
• Van Leemput K.
Fast and sequence-adaptive whole-brain segmentation using parametric Bayesian modeling.
] or ROAST [
• Huang Y.
• Datta A.
• Bikson M.
• Parra L.C.
Realistic volumetric-approach to simulate transcranial electric stimulation-ROAST-a fully automated open-source pipeline.
] – segment the head into skin, skull, cerebrospinal fluid (CSF), gray matter (GM), and white matter (WM). The meninges (dura, arachnoid, and pia mater) cannot be segmented using the existing tools.
Previous work [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] has suggested that the meninges are included within the volume conventionally modeled as homogenous, highly-conductive CSF. However, they possess much lower electrical conductivities, thus representing additional layers to pass before reaching the brain [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
]. Previous modeling and experimental work [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] has shown that the meningeal layers may significantly (by 16–60%) increase TES induced cortical fields as compared to the standard modeling. To emulate the effect of meninges, it has been suggested to replace the conventional CSF conductivity of 1.65 S/m by 0.85 S/m [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
]. This result illustrates the potential impact of the meningeal layers on the electric field estimates. At the same time, this result was based on a spherical head model, which simplifies the electric field distribution within the brain [
• Opitz A.
• Windhoff M.
• Heidemann R.M.
• Turner R.
• Thielscher A.
How the brain tissue shapes the electric field induced by transcranial magnetic stimulation.
].
For this reason, the present study aims to independently validate and generalize this observation and justify the suggested CSF conductivity approximation considering realistic head models. We perform numerical modeling for 16 subjects drawn from the Human Connectome Project (HPC) database [
• Van Essen D.C.
• Ugurbil K.
• Auerbach E.
• Barch D.
• Behrens T.E.
• Bucholz R.
• Chang A.
• Chen L.
• Corbetta M.
• Curtiss S.W.
• Della Penna S.
• Feinberg D.
• Glasser M.F.
• Harel N.
• Heath A.C.
• Larson-Prior L.
• Marcus D.
• Michalareas G.
• Moeller S.
• Oostenveld R.
• Petersen S.E.
• Prior F.
• Schlaggar B.L.
• Smith S.M.
• Snyder A.Z.
• Xu J.
• Yacoub E.
The Human Connectome Project: a data acquisition perspective.
] with 0.7 mm isotropic resolution. Along with TES, we also include transcranial magnetic stimulation (TMS) into our analysis, thus extending the scope of the prior work [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
].
The major methodological difficulty is in the accurate modeling of thin meningeal layers with average thicknesses of 1.11 mm, 0.2 mm, and 0.1 mm for dura, arachnoid, and pia mater, respectively [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
]. These layers are located very close to the neighboring brain compartments. Irrespective of which numerical method is used, the geometrical element sizes should be smaller than or equal to the separation distances between individual tissue interfaces for accurate computation. This could lead to a large computational burden. Modeling thin layers is computationally expensive when using the finite element method (FEM) due to the requirement of an adequate volumetric mesh [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
].
To overcome this issue, we used the boundary element fast multipole method (BEM-FMM) [
• Makarov S.N.
• Wartman W.A.
• Nguyen B.T.
• Noetscher G.
• Ahveninen J.P.
• Fujimoto K.
• Weise K.
• Nummenmaa A.
Boundary element fast multipole method for modeling electrical brain stimulation with voltage and current electrodes.
], which only requires the discretization of the surfaces. Furthermore, we extended the existing BEM-FMM algorithm by an automated adaptive mesh refinement mechanism introduced and described in detail in Appendix A. With this extension, we were able to determine accurate and self-converged results for the electric field with reasonable computational cost. Our method could be applicable to any variable-thickness meningeal geometry, and also to meninges in the form of thin separate islands.

## 2. Materials and methods

MRI T1/T2 data for 16 healthy HCP subjects [
• Van Essen D.C.
• Ugurbil K.
• Auerbach E.
• Barch D.
• Behrens T.E.
• Bucholz R.
• Chang A.
• Chen L.
• Corbetta M.
• Curtiss S.W.
• Della Penna S.
• Feinberg D.
• Glasser M.F.
• Harel N.
• Heath A.C.
• Larson-Prior L.
• Marcus D.
• Michalareas G.
• Moeller S.
• Oostenveld R.
• Petersen S.E.
• Prior F.
• Schlaggar B.L.
• Smith S.M.
• Snyder A.Z.
• Xu J.
• Yacoub E.
The Human Connectome Project: a data acquisition perspective.
] with isotropic voxel resolution of 0.7 mm used previously [
• Makarov S.N.
• Wartman W.A.
• Noetscher G.M.
• Fujimoto K.
• Zaidi T.
• Burnham E.H.
• Daneshzand M.
• Nummenmaa A.
Degree of improving TMS focality through a geometrically stable solution of an inverse TMS problem.
] have been segmented using the automated segmentation routine headreco (based on SPM/CAT [
Structural Brain Mapping Group
Computational anatomy toolbox (CAT).
]), available in the SimNIBS 3.1 software. This segmentation is termed as the base (non-meningeal) segmentation; it is shown in Fig. 1a–d for one subject. Fig. 1a–d also show that the dura mater is included into the CSF volume in motor, premotor, and somatosensory areas. More details are provided in Appendix D.

### 2.2 Assembly of brain meninges

Given the current resolution limits of the MRI devices, the available segmentation routines [
• Fischl B.
FreeSurfer.
,
Structural Brain Mapping Group
Computational anatomy toolbox (CAT).
] are not capable of predicting the meninges’ exact locations. Therefore, we created the meningeal layers by conjecturing their shape based on average anatomical data [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] and following a fully automated routine; the final result is illustrated in Fig. 1e and f. This routine creates meninges from constrained deformations of the base segmentation and is included in the example dataset of interest (DOI) [
• TES & TMS BEM-FMM Solver with Adaptive Mesh Refinement
Connectome subject 120111. MATLAB windows/linux code and user's manual. β-Version january 2022.
]. According to Fig. 1c and d, headreco puts the dura mater into the CSF volume (cf. also [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
]). We therefore proceed as follows:
• 1.
First, the original CSF shell is designated as dura. Next, the new CSF shell is created by moving the original one inwards along the normal direction to account for dura and arachnoid mater thicknesses: 1.11 mm for dura +0.2 mm for arachnoid. If intersections occur with the GM surface (later pia mater surface), the dura and arachnoid thicknesses are spatially smoothly reduced to provide at least 0.1 mm separation in the normal direction between the new CSF and pia mater – cf. insets in Fig. 1e and f.
• 2.
The original GM shell is designated as pia mater. The new GM shell is created by moving the original shell inwards along the normal direction to account for a pia mater thickness of 0.1 mm. If there is not enough space (rare), then half of the distance between GM and WM is chosen as the pia thickness.
Here, we used integral data given in Ref. [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] for the meningeal layer thicknesses, which were estimated from literature [
• Bashkatov A.N.
• Genina E.A.
• Sinichkin Y.P.
• Kochubey V.I.
• Lakodina N.A.
• Tuchin V.V.
Glucose and mannitol diffusion in human dura mater.
,
• Kuchiwaki H.
• Inao S.
• Ishii N.
• Ogura Y.
• Gu S.P.
Human dural thickness measured by ultrasonographic method: reflection of intracranial pressure.
,
• Saboori P.
Histology and morphology of the brain subarachnoid trabeculae.
,
• Fournier M.
• Combès B.
• Roberts N.
• Braga J.
• Prima S.
Mapping the distance between the brain and the inner surface of the skull and their global asymmetries Medical Imaging 2011: image Processing Medical Imaging 2011.
] noting that references vary due to the heterogeneous geometry of the layers themselves. All three meninges of the cortex fold and descend deep down into the longitudinal fissure, where we have followed the original base segmentation. Since the selected stimulation sites ($M1HAND$ and $DLPFC$) are quite distant from this feature, we assume that this has only a minor influence on the results.

### 2.3 Other modeling sets

In addition to these full models, several additional head models were constructed. One set of models includes only the dura mater combined with arachnoid mater because it is the thickest membrane, and we expect it to have the largest effect on the electric field. In another set of models, only the thin pia mater was included. In the last set of models, the dura and arachnoid have been gradually moved in the normal direction (progressively encroaching on the skull volume). This was done to study the field sensitivity with respect to membrane position and to make a statement regarding required resolution of future imaging methods.

### 2.4 Stimulation setups and ROI definition

For focal TES, two ring electrode configurations [
• Datta A.
• Bansal V.
• Diaz J.
• Patel J.
• Reato D.
• Bikson M.
Gyri-precise head model of transcranial direct current stimulation: improved spatial focality using a ring electrode versus conventional rectangular pad.
] labeled as Montage 1 and Montage 2 in Fig. 2 have been used. The center electrodes were placed over two targets: the approximate center of the $M1HAND$ area of the left hemisphere and the approximate center of the dorsolateral left prefrontal cortex ($DLPFC$) at EEG electrode F3 [
• Jurcak V.
• Tsuzuki D.
• Dan I.
10/20, 10/10, and 10/5 systems revisited: their validity as relative head-surface-based positioning systems.
,
• Herwig U.
• Satrapi P.
• Schönfeldt-Lecuona C.
Using the international 10-20 EEG system for positioning of transcranial magnetic stimulation.
]. The other four electrodes were positioned in the anterior-posterior and median-lateral directions, at distances of 27.5 mm (Montage 1) or 40 mm (Montage 2) from the center electrode. The electrode radii were selected as 5 mm or 10 mm, respectively. The center electrode was set at -1V, and all other electrodes were set at +1V. The results were linearly scaled to achieve a total injected current of 1 mA. Fig. 2a–d shows these four cases.
For conventional TES, Montage A in Fig. 2e targets the $M1HAND$ area of the left hemisphere with the second electrode placed at F4. The square electrode size is 2 × 2 cm. Montage B in Fig. 2f is the same as A, except the electrode size is 4 × 4 cm. In Montage C, both electrodes are additionally shifted posteriorly by 2 cm to place the maximum ROI field in the $M1HAND$ area. Montage D is the same as C, except an ipsilateral configuration was assumed for the second electrode instead of a contralateral configuration.
For TMS computations, we placed a large Magventure C-B60 figure-of-eight coil or a small Cool-B35 coil over the same $M1HAND$ or $DLPFC$ areas, parallel to the skin, and with an orientation of 45° towards the fissura longitudinalis as shown in Fig. 2i-l. The coil-skin distances along the centerline were 10 mm accounting for the housing of the coil.
Conductivity values of the meninges follow Ref. [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] and were chosen as: dura - 0.1 S/m; arachnoid - 0.125 S/m; pia - 0.15 S/m. Other conductivity values are the default conductivities of SimNIBS [
• Wagner T.A.
• Zahn M.
• Grodzinsky A.J.
• Pascual-Leone A.
Three-dimensional head model simulation of transcranial magnetic stimulation.
,
• Thielscher A.
• Opitz A.
• Windhoff M.
Impact of the gyral geometry on the electric field induced by transcranial magnetic stimulation.
]: skin - 0.465 S/m; skull - 0.01 S/m; CSF - 1.654 S/m; GM - 0.275 S/m; WM - 0.126 S/m. The region of interest (ROI, ∼4000 observation points) is located on the mid-surface between the inner and outer GM surfaces within a sphere (D = 40 mm) centered on a target at the center of the $M1HAND$ region or the $DLPFC$ region.

### 2.5 Modeling engine

The boundary element fast multipole method (BEM-FMM) [
• Makarov S.N.
• Noetscher G.M.
• Raij T.
• Nummenmaa A.
A quasi-static boundary element approach with fast multipole acceleration for high-resolution bioelectromagnetic models.
] was used to determine the electric field distributions for both TES and TMS. It operates in terms of the surface charge density residing on interfaces and electrode surfaces. BEM-FMM has been shown to outperform the standard first-order finite element method for TMS modeling [
• Gomez L.J.
• Dannhauer M.
• Koponen L.M.
• Peterchev A.V.
Conditions for numerically accurate TMS electric field simulation.
] when isotropic and constant-conductivity compartments are used. It has been applied to TMS modeling [
• Makarov S.N.
• Wartman W.A.
• Noetscher G.M.
• Fujimoto K.
• Zaidi T.
• Burnham E.H.
• Daneshzand M.
• Nummenmaa A.
Degree of improving TMS focality through a geometrically stable solution of an inverse TMS problem.
,
• Makarov S.N.
• Wartman W.A.
• Daneshzand M.
• Fujimoto K.
• Raij T.
• Nummenmaa A.
A software toolkit for TMS electric-field modeling with boundary element fast multipole method: an efficient MATLAB implementation.
,
• Daneshzand M.
• Makarov S.N.
• de Lara L.I.N.
• Guerin B.
• McNab J.
• Rosen B.R.
• Hämäläinen M.S.
• Raij T.
• Nummenmaa A.
Rapid computation of TMS-induced E-fields using a dipole-based magnetic stimulation profile approach.
], to EEG/MEG modeling [
• Makarov S.N.
• Hämäläinen M.
• Noetscher G.M.
• Ahveninen J.
• Nummenmaa A.
Boundary element fast multipole method for enhanced modeling of neurophysiological recordings.
], and most recently to TES modeling [
• Makarov S.N.
• Wartman W.A.
• Nguyen B.T.
• Noetscher G.
• Ahveninen J.P.
• Fujimoto K.
• Weise K.
• Nummenmaa A.
Boundary element fast multipole method for modeling electrical brain stimulation with voltage and current electrodes.
]. Its advantage is an ability to efficiently model a complicated multi-shell topology.

### 2.6 Reference solutions

Although BEM-FMM can handle extremely closely spaced shells [
• Makarov S.N.
• Wartman W.A.
• Nguyen B.T.
• Noetscher G.
• Ahveninen J.P.
• Fujimoto K.
• Weise K.
• Nummenmaa A.
Boundary element fast multipole method for modeling electrical brain stimulation with voltage and current electrodes.
], many analytical neighbor potential integrals per facet must be stored in RAM [
• Makarov S.N.
• Wartman W.A.
• Nguyen B.T.
• Noetscher G.
• Ahveninen J.P.
• Fujimoto K.
• Weise K.
• Nummenmaa A.
Boundary element fast multipole method for modeling electrical brain stimulation with voltage and current electrodes.
]. This makes the solution less practical. The headreco shells are separated by less than the average edge length of 1.4 mm at many places. The meningeal interfaces make the problem much worse. Therefore, two “ground truth” models of a very high computational cost were constructed.
The first such model – Ref#1 – was constructed from the original headreco meshes using two levels of a uniform edge-based congruent triangle subdivision (1:4). This results in 17 million facets per head with an average edge length of 0.36 mm, and the same average triangle quality. This model provides an accurate solution without meninges.
The second model – Ref#2 – is constructed from the original headreco meshes using three levels of the uniform edge-based congruent subdivision (1:4) for dura mater, arachnoid mater, CSF, pia mater, gray matter, and two levels otherwise. This results in 70 million facets per head with the average edge length of 0.18 mm for the critical compartments. This model arguably provides an accurate solution with the meninges.

Solving the large models constructed above was possible but very time consuming. Therefore, we propose an alternative solution, which is based on an adaptive mesh refinement algorithm illustrated in Fig. 3.
This means that the mesh is refined automatically and selectively in those regions where it is required, by subdividing triangular facets in question into smaller facets and recomputing the solution. The subdivision process is controlled by a cost function. Appendix A specifies the algorithm, its implementation, and associated convergence. A final refined head model resulting in a converged solution for TES is illustrated in Fig. 3.

### 2.8 Software and data availability

The BEM-FMM solver with the adaptive mesh refinement as well as with the automated mesh generator for the meninges with an example are available online [
• TES & TMS BEM-FMM Solver with Adaptive Mesh Refinement
Connectome subject 120111. MATLAB windows/linux code and user's manual. β-Version january 2022.
].

## 3. Results

### 3.1 Modeling results when meninges are included into CSF volume

For the 16 subject models considered, the ground truth results for the uniformly refined models and the modeling results obtained with the adaptive mesh refinement are nearly identical. Their equivalence, parameters of the numerical algorithm, and algorithm's performance are reported in detail in Appendix A. For this reason, the electric field differences are reported here using the adaptive mesh refinement method.
Table 1 summarizes the ROI field differences between models with and without the meninges, respectively, for focal TES. Table 2 summarizes similar results for conventional TES. The corresponding TMS results are given in Table 3. All three tables report relative signed differences in the maximum ROI electric field magnitude (the 98th percentile), in the mean ROI electric field magnitude, and the average relative 2-norm difference for the vector electric field, for all 4000 ROI sampling points. All these results have been averaged over 16 subjects. All three tables report data for four different stimulation combinations and three different modeling sets (dura and arachnoid only, pia only, all three meningeal layers). The reported differences are:
$Equation 1.$
(1)

where the 98th percentile maximum of the electric field magnitude, $|E|$, is found over the entire ROI;
$Equation 2.$
(2)

where the mean of the electric field magnitude, $|E|$, is computed over the entire ROI;
$Equation 3.$
(3)

where the 2-norm of the full vector field $E$ is computed over the entire ROI. Additional statistical data (STDs and p-values of the paired t-test) are given in Appendix B.
Table 1Averaged relative differences over 16 subjects in the ROI electric field for models with and without the meninges: TES modeling with focal montages. Three cases – when only dura and arachnoid mater are included, when only pia mater is included, or when all meninges are included – have been computed and compared with the base non-meningeal solutions. Relevant statistical data (standard deviation and p-values of the paired t-test) are given in Appendix B.
Model

Includes
Stimulation setupSigned $Diffmax$, %Signed $Diffmean$, %$Diff2norm$, %
Dura and Arachnoid onlyTarget M1, Montage 1+40+2733
Target DLPFC, Montage 1+36+2730
Target M1, Montage 2+40+2733
Target DLPFC, Montage 2+35+2630
Pia onlyTarget M1, Montage 1−6−66
Target DLPFC, Montage 1−7−77
Target M1, Montage 2−8−77
Target DLPFC, Montage 2−8−88
All meningesTarget M1, Montage 1+31+1624
Target DLPFC, Montage 1+28+1520
Target M1, Montage 2+23+1217
Target DLPFC, Montage 2+27+1321
Table 2Averaged relative differences over 16 subjects in the ROI electric field for models with and without the meninges: TES modeling with conventional montages. Three cases – when only dura and arachnoid mater are included, when only pia mater is included, or when all meninges are included – have been computed and compared with the base non-meningeal solutions. Relevant statistical data (standard deviation and p-values of the paired t-test) are given in Appendix B.
Model

Includes
Stimulation setupSigned $Diffmax$, %Signed $Diffmean$, %$Diff2norm$, %
Dura and Arachnoid onlyTarget M1, Montage A+26+1820
Target M1, Montage B+22+1618
Target M1, Montage C+22+1416
Target M1, Montage D+19+1415
Pia onlyTarget M1, Montage A−7−79
Target M1, Montage B−9−810
Target M1, Montage C−10−1011
Target M1, Montage D−7−78
All meningesTarget M1, Montage A+15+615
Target M1, Montage B+6−112
Target M1, Montage C+5−48
Target M1, Montage D+10+57
Table 3Averaged relative differences over 16 subjects in the ROI electric field for models with and without the meninges: TMS modeling. Three cases – when only dura and arachnoid mater are included, when only pia mater is included, or when all meninges are included – have been computed and compared with the base non-meningeal solutions. Relevant statistical data (standard deviation and p-values of the paired t-test) are given in Appendix B.
Model

Includes
Stimulation setupSigned $Diffmax$, %Signed $Diffmean$, %$Diff2norm$, %
Dura and Arachnoid onlyTarget M1, C-B60 coil+1+15
Target DLPFC, C-B60 coil−2−0.14
Target M1, Cool-B35 coil+2+15
Target DLPFC, Cool-B35 coil−2−0.24
Pia onlyTarget M1, C-B60 coil−5−23
Target DLPFC, C-B60 coil−5−23
Target M1, Cool-B35 coil−4−16
Target DLPFC, Cool-B35 coil−4−23
All meningesTarget M1, C-B60 coil−5−16
Target DLPFC, C-B60 coil−8−36
Target M1, Cool-B35 coil−4−16
Target DLPFC, Cool-B35 coil−7−35
Along with Eq. (3), we also computed the topographical difference or the relative difference measure (RDM). Its average value approaches 60% of the relative 2-norm difference (3) on average and is therefore not reported separately.
Additionally, Fig. 4 shows the electric field magnitude in the ROI with and without the meninges for one subject (HCP subject 120111) for focal TES (Montage 1), conventional TES (Montage C) and TMS (with the C-B60 coil), all targeting the $M1HAND$ area.
According to Table 3 (and Fig. 4), the effect of the meninges on the TMS cortical fields generally appears to be small. This is in stark contrast to the focal TES cortical fields from Table 1 (and Fig. 4) where this effect is quite large. Therefore, the membrane emulation results as well as the sensitivity results versus membrane positions and conductivities will be reported below only for the focal TES case.

### 3.2 Variation of CSF conductivity to emulate the meningeal effect (focal TES)

Jiang et al. (2020) [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] suggested to decrease the conventional CSF conductivity of 1.65 S/m to 0.85 S/m to emulate the effect of the meninges without including them into the base model explicitly. Fig. 5 shows the average effects of CSF conductivity variations from 0.75 S/m to 1.65 S/m in steps of 0.05 S/m, in the base model for the focal TES Montage 1 targeting the $M1HAND$ area with respect to all three difference measures reported previously in Table 1. The first conductivity domain marked light yellow in Fig. 5 corresponds to the value of 0.85 S/m. In this domain, the relative signed difference in the maximum ROI field crosses zero and does become quite small while its standard deviation does not exceed 10%. For all focal TES montages from Table 1, the $Diffmax$ values for the full set of membranes (four last rows of Table 1) reduce from +31, +28, +23, +27% to +3, +1, +3, +1% that is, they decrease by a factor of ∼10 when a CSF conductivity of 0.85 S/m is used! However, the 2-norm difference and the mean-field difference are still significant at 0.85 S/m. The second domain in Fig. 5 corresponds to the value of 1.25 S/m. In this region, the 2-norm difference and the mean-field difference reach their absolute minima. The $Diffmean$ values for the full set of membranes reduce from +16, +15, +12, +13% to +1, +2, −3, −1%; that is, they again decrease by a nearly a factor of ∼10. However, the relative signed difference in the maximum ROI field is relatively large.

### 3.3 Variation of the meningeal conductivity (focal TES)

Appendix C considers an extreme case where both dura and arachnoid mater have the conductivity of skull. This might happen if a hypothetical non-meningeal segmentation would move the skull boundary down to “subtract” the meninges from the CSF volume.

### 3.4 Influence of segmentation errors of dura and arachnoid mater (focal TES)

The previous sections have considered the most likely case, where the meninges have been included into the homogeneous CSF volume [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] of the original “non-meningeal” segmentation using headreco. In order to investigate the effects of possible segmentation errors, we gradually moved the dura and arachnoid mater in the direction of the skull. They would then be located, partially or fully, in the skull volume of the original segmentation. Table 4 reports the average relative ROI field change (i.e. sensitivity) per 0.1 mm of such a complex movement computed numerically. Table 5 reports the results when the dura and arachnoid mater are fully included into the skull volume. This table compares the meningeal and non-meningeal (base segmentation) solutions, similarly to Table 1.
Table 4Average relative unsigned sensitivity over 16 subjects of the ROI electric field for focal TES per a 0.1 mm meningeal position movement when the dura and arachnoid mater are gradually moved upwards, i.e. from the CSF volume into the skull volume. The pia mater stays the same. The maximum and mean fields decrease when moving upwards and approximately approach the result of the base segmentation without meninges when the meninges become fully included into the skull volume instead of the CSF volume (cf. Table 5 below).
Stimulation setupSensitivity for $max98%(|E|)$Sensitivity for $mean(|E|)$Sensitivity for $E2norm$
Target M1, Montage 13%1%2%
Target DLPFC, Montage 13%1%1%
Target M1, Montage 23%1%2%
Target DLPFC, Montage 22%1%2%
Table 5Averaged signed relative differences over 16 subjects in the ROI electric field for focal TES models with and without the meninges. The model without the meninges is the base segmentation. The model with the meninges now fully includes them into the skull volume instead of the CSF volume.
Stimulation setup$Diffmax$, %$Diffmean$, %$Diff2norm$, %
Target M1, Montage 1+6+47
Target DLPFC, Montage 1+5+36
Target M1, Montage 2+1+25
Target DLPFC, Montage 2+1+25

## 4. Discussion

### 4.1 Focal TES results: agreement with the prior study

The present modeling results of focal TES montages agree with the prior study [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] in that they also predict a significant average increase in the cortical electric field (an average increase of approximately 30% in the maximum field strength) caused by the meningeal layers. Next, we confirm that almost exactly the same value of 0.85 S/m for the CSF conductivity [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] is necessary to emulate the effect of the meninges when the maximum electric field strength within the ROI is concerned. However, in order to decrease the average electric field difference in the ROI, a value of 1.25 S/m should be considered. These two results are consistent; they are confirmed for all four focal TES stimulation setups from Table 1.

### 4.2 Focal TES results: new findings

Surprisingly, the very thin pia mater alone seems to have a rather significant effect on TES intracortical fields. According to Table 1, these fields decrease by approximately 7% on average for all focal TES stimulation setups when only the pia mater is included. We explain this fact by the appearance of a thin layer of electric dipoles (a double layer) across the thin and low-conducting pia mater. Although the field of dipoles is small at larger distances, their near field is significant and could therefore affect the very nearby ROI.
On the other hand, the dura and arachnoid mater have a strong opposite effect on the intracortical fields since they effectively reduce the field-blocking CSF volume and make the transition between bone and CSF more “smooth”. As a result, an increase of 37% in the field strength is observed in Table 1 on average.
Together, these two effects partially cancel each other, but the effects of dura and arachnoid mater still dominate. A total increase of around 30% in the field strength is therefore observed in Table 1 on average when all three meningeal layers are present in the otherwise-homogeneous CSF volume. These three observations are consistent; they are confirmed for all four TES stimulation setups from Table 1.

### 4.3 Focal TES results: critical influence of meningeal positions

It follows from Sections 3.1, 3.3 that the main question is whether the base non-meningeal segmentation (headreco in the present study) includes the dura and arachnoid mater into the CSF volume or into the skull volume. In the former case, large differences in the maximum ROI electric field on the order of 30% are expected when the meninges are extracted from the CSF volume (Table 1). In the latter case, an average difference in the maximum ROI electric field of only 3% is expected when the meninges are extracted from the skull volume (Table 5). The transition between the two cases is very sharp: the relative sensitivity of the maximum ROI field is 3% per 0.1 mm meningeal movement on average (Table 4). Other potential sources of electric field error are artifacts in the MR images, which may result in electric field differences of up to 30% [
• Puonti O.
• Van Leemput K.
• Saturnino G.B.
• Siebner H.R.
• Thielscher A.
Accurate and robust whole-head segmentation from magnetic resonance images for individualized head modeling.
].

### 4.4 Meningeal positions for base headreco segmentations

Ref. [
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
] and Fig. 1 of this study both suggest that the base segmentation subsumes the dura and arachnoid mater into the CSF volume. This may not always be the case. Appendix D reports headreco segmentations overlapped with the original T1/T2 images for several subjects of this study. A visual inspection (performed wherever possible) indicates that dura and arachnoid mater are subsumed into the CSF volume for the primary and supplementary motor cortexes, primary somatosensory cortex, and partially the posterior parietal cortex. However, for the occipital lobe and frontal lobes, the base headreco segmentation better follows the anatomical CSF boundary. This shows that this aspect should be given high attention in the development of future segmentation methods.

### 4.5 Results for conventional bipolar TES montages

For the conventional TES montages (Table 2), we observed that the pia mater has a similar effect on the electric field compared to focal TES (Table 1). On the other hand, the effect of the dura mater appears smaller than for focal TES. In combination, these two effects tend to cancel each other more strongly, depending on the individual geometry. As a result, a total increase of only 9% in the field strenth is observed in Table 2 on average, on the border of statistical significance (supplementary Table B2 of Appendix B). We hypothesize that this observation may be related to the larger spread of the electric field associated with conventional as opposed to focal TES montages.

### 4.6 TMS results

In case of TMS, the effect of the meninges is much less pronounced (Table 3). However, a small, yet statistically significant, decrease in the maximum ROI field of 6% on average is consistently observed. The shunting properties of CSF and its thickness have less effect here since the dominant E-field of the coil itself is mostly parallel to the CSF boundaries: because the normal component of the E-field is small, the induced charges are also small, so they only generate a weak secondary field that can barely influence the total electric field. Still, the normal component of the E-field might be substantial at the sulcal walls, where a thin blocking dipole layer appears at the pia mater and the gray matter surface. Hence, most of the effects in Table 3 arise from the thin pia mater because it directly borders the gray matter, which is very close to the ROI.

### 4.7 Solutions with adaptive mesh refinement

The use of adaptive mesh refinement is new and innovative in the context of brain stimulation modeling. It is also used in other engineering fields, for example in communications and defense applications. The BEM-FMM method used here solves the uniformly refined ground truth TES/TMS models of exemplary sizes (∼70 million surface facets) and with mesh resolutions of 0.18 mm in approximately 2 h using a multicore 2.8 GHz server with 512 Gbytes of RAM. The use of adaptive mesh refinement based on Eq. (4) of Appendix A allows us to reduce this time to approximately 10–20 min while utilizing less demanding computer resources and achieving better mesh resolution (on the order of 0.05 mm) in critical areas.

## 5. Conclusionsconclusion

The effect of thin meningeal layers on the electric field is a difficult matter. Thanks to BEM-FMM with adaptive mesh refinement, we were able to accurately model ∼1,000 different TES and TMS cases for 16 subjects with the average meningeal thicknesses. Our conclusions are:
• i.
Explicit inclusion of the three meningeal layers is important for focal TES computations if a base non-meningeal segmentation subsumes the meninges into the CSF volume (cf. Appendix D). It may result in differences (increases) in the maximum cortical electric field of 30% on average when focal TES is employed (Table 1).
• ii.
Explicit inclusion of the three meningeal layers is less important for focal TES computations if a base non-meningeal segmentation correctly predicts the CSF volume but combines the meninges with the skull volume (Appendix D). Differences (increases) in the maximum cortical electric fields of 3% on average (Table 5) may result.
• iii.
The transition between the two cases is sharp: the relative sensitivity of the maximum ROI E-field is 3% per 0.1 mm meningeal movement on average (Table 4).
• iv.
The pia mater alone may generate an 8% difference in the maximum cortical focal TES electric fields on average. All three meninges must be considered simultaneously since their effects partially cancel (Table 1, Table 2).
• v.
While all results for the focal TES (Table 1 and Table B1) are very consistent and statistically significant across subjects, the results for the conventional bipolar TES (Table 2 and Table B2) are less classifiable; they show marginal statistical significance of the overall meningeal effect.
• vi.
Inclusion of the three meningeal layers into TMS computations seems to be much less important. Still, a maximum cortical electric field difference (decrease) of approximately 6% on average may result (Table 3).
• vii.
The worst-case meningeal effect on focal TES (∼30% maximum field change) is comparable to that of a segmentation error of ∼1.2 mm in the eccentric CSF volumes. However, the meningeal effect cannot be entirely undone by simply modifying the surfaces of the base non-meningeal segmentation (Appendix C).
• viii.
The influence of the meninges varies greatly between the stimulation modalities. In the case of classical TES, we observed a stronger position dependence compared to focal TES and TMS. If precise quantification of the electric field magnitude is important (e.g., in the context of patient safety), we recommend individualized numerical simulations that also account for the meninges. Future studies may investigate the impact of uncertainties in the tissue conductivities on the electric field distributions [
• Saturnino G.B.
• Thielscher A.
• Knösche T.R.
• Weise K.
A principled approach to conductivity uncertainty analysis in electric field calculations.
,
• Weise K.
• Poßner L.
• Müller E.
• Gast R.
• Knösche T.R.
Pygpc: a sensitivity and uncertainty analysis toolbox for Python.
].

## CRediT authorship contribution statement

Konstantin Weise: Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Project administration, Supervision, Writing – original draft, Writing – review & editing. William A. Wartman: Data curation, Formal analysis, Investigation, Methodology, Software, Visualization, Writing – review & editing. Thomas R. Knösche: Conceptualization, Investigation, Project administration, Supervision, Writing – review & editing. Aapo R. Nummenmaa: Conceptualization, Funding acquisition, Project administration, Supervision, Writing – review & editing. Sergey N. Makarov: Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Software, Supervision, Visualization, Writing – original draft, Writing – review & editing.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Acknowledgements

Konstantin Weise: This work was partially supported by the German Science Foundation (DFG) (grant number WE 5985/1-2 ) and the NVIDIA Corporation (donation of a Titan Xp graphics card to KW).
Aapo Nummenmaa & Sergey Makarov: The authors are grateful to Dr. Leslie Greengard and Dr. Manas Rachh of Flatiron Inst. and Courant Inst. of Mathematical Sciences (NYC, USA) for their continuous algorithmic and conceptual FMM support. This work has been partially funded by P41EB030006 , a P41 Center for Mesoscale Modeling Resource Grant supported by the National Institute of Biomedical Imaging and Bioengineering ( NIBIB, National Institutes of Health h, NIH , USA , and by NIMH / NINDS / NIH grants 1R01MH111829 and 1R01MH128421 , National Institutes of Health , NIH , USA .

## Appendix A. Supplementary data

• Meninges segmentation location
• Dura mater in MRI

## References

• Weise K.
• Numssen O.
• Thielscher A.
• Hartwigsen G.
• Knösche T.R.
A novel approach to localize cortical TMS effects.
Neuroimage. 2020 Apr 1; 209116486https://doi.org/10.1016/j.neuroimage.2019.116486
• Numssen O.
• Zier A.L.
• Thielscher A.
• Hartwigsen G.
• Knösche T.R.
• Weise K.
Efficient high-resolution TMS mapping of the human motor cortex by nonlinear regression.
Neuroimage. 2021 Dec 15; 245 (Epub 2021 Oct 12. PMID: 34653612)118654https://doi.org/10.1016/j.neuroimage.2021.118654
• Peterchev A.V.
• Wagner T.A.
• Miranda P.C.
• Nitsche M.A.
• Paulus W.
• Lisanby S.H.
• Pascual-Leone A.
• Bikson M.
Fundamentals of transcranial electric and magnetic stimulation dose: definition, selection, and reporting practices.
Brain Stimul. 2012 Oct; 5: 435-453https://doi.org/10.1016/j.brs.2011.10.001
• Puonti O.
• Van Leemput K.
• Saturnino G.B.
• Siebner H.R.
• Thielscher A.
Accurate and robust whole-head segmentation from magnetic resonance images for individualized head modeling.
Neuroimage. 2020 Oct 1; 219117044https://doi.org/10.1016/j.neuroimage.2020.117044
• Windhoff M.
• Opitz A.
• Thielscher A.
Electric field calculations in brain stimulation based on finite elements: an optimized processing pipeline for the generation and usage of accurate individual head models.
Hum Brain Mapp. 2013 Apr; 34 (Epub 2011 Nov 23. PMID: 22109746; PMCID: PMC6870291): 923-935https://doi.org/10.1002/hbm.21479
• Puonti O.
• Iglesias J.E.
• Van Leemput K.
Fast and sequence-adaptive whole-brain segmentation using parametric Bayesian modeling.
Neuroimage. 2016; 143: 235-249https://doi.org/10.1016/j.neuroimage.2016.09.011
• Huang Y.
• Datta A.
• Bikson M.
• Parra L.C.
Realistic volumetric-approach to simulate transcranial electric stimulation-ROAST-a fully automated open-source pipeline.
J Neural Eng. 2019 Jul 30; 16056006https://doi.org/10.1088/1741-2552/ab208d
• Jiang J.
• Truong D.Q.
• Esmaeilpour Z.
• Huang Y.
• Bikson M.
Enhanced tES and tDCS computational models by meninges emulation.
J Neural Eng. 2020 Jan 14; 17016027https://doi.org/10.1088/1741-2552/ab549d
• Opitz A.
• Windhoff M.
• Heidemann R.M.
• Turner R.
• Thielscher A.
How the brain tissue shapes the electric field induced by transcranial magnetic stimulation.
Neuroimage. 2011 Oct 1; 58: 849-859https://doi.org/10.1016/j.neuroimage.2011.06.069
• Van Essen D.C.
• Ugurbil K.
• Auerbach E.
• Barch D.
• Behrens T.E.
• Bucholz R.
• Chang A.
• Chen L.
• Corbetta M.
• Curtiss S.W.
• Della Penna S.
• Feinberg D.
• Glasser M.F.
• Harel N.
• Heath A.C.
• Larson-Prior L.
• Marcus D.
• Michalareas G.
• Moeller S.
• Oostenveld R.
• Petersen S.E.
• Prior F.
• Schlaggar B.L.
• Smith S.M.
• Snyder A.Z.
• Xu J.
• Yacoub E.
The Human Connectome Project: a data acquisition perspective.
Neuroimage. 2012; 62 (Online (Jan. 2020):): 2222-2231
• Makarov S.N.
• Wartman W.A.
• Nguyen B.T.
• Noetscher G.
• Ahveninen J.P.
• Fujimoto K.
• Weise K.
• Nummenmaa A.
Boundary element fast multipole method for modeling electrical brain stimulation with voltage and current electrodes.
J Neural Eng. 2021 Jul 26; https://doi.org/10.1088/1741-2552/ac17d7
• Makarov S.N.
• Wartman W.A.
• Noetscher G.M.
• Fujimoto K.
• Zaidi T.
• Burnham E.H.
• Daneshzand M.
• Nummenmaa A.
Degree of improving TMS focality through a geometrically stable solution of an inverse TMS problem.
Neuroimage. 2021 Jul 28; 241: 118437https://doi.org/10.1016/j.neuroimage.2021.118437
• Fischl B.
FreeSurfer.
NeuroImage. 2012; 62 (PubMed PMID: 22248573): 774-781https://doi.org/10.1016/j.neuroimage.2012.01.021
• Structural Brain Mapping Group
Computational anatomy toolbox (CAT).
(Univ. of Jena, Germany. Accessed 04/05/21. Online:)
• TES & TMS BEM-FMM Solver with Adaptive Mesh Refinement
Connectome subject 120111. MATLAB windows/linux code and user's manual. β-Version january 2022.
(Online:)
• Bashkatov A.N.
• Genina E.A.
• Sinichkin Y.P.
• Kochubey V.I.
• Lakodina N.A.
• Tuchin V.V.
Glucose and mannitol diffusion in human dura mater.
Biophys J. 2003 Nov; 85: 3310-3318https://doi.org/10.1016/S0006-3495(03)74750-X
• Kuchiwaki H.
• Inao S.
• Ishii N.
• Ogura Y.
• Gu S.P.
Human dural thickness measured by ultrasonographic method: reflection of intracranial pressure.
J Ultrasound Med. 1997 Nov; 16: 725-730https://doi.org/10.7863/jum.1997.16.11.725
• Saboori P.
Histology and morphology of the brain subarachnoid trabeculae.
Anat Res Int. 2015; 2015: 279814https://doi.org/10.1155/2015/279814
• Fournier M.
• Combès B.
• Roberts N.
• Braga J.
• Prima S.
Mapping the distance between the brain and the inner surface of the skull and their global asymmetries Medical Imaging 2011: image Processing Medical Imaging 2011.
in: Image processing. vol. 7962. International Society for Optics and Photonics, 201179620Y
• Datta A.
• Bansal V.
• Diaz J.
• Patel J.
• Reato D.
• Bikson M.
Gyri-precise head model of transcranial direct current stimulation: improved spatial focality using a ring electrode versus conventional rectangular pad.
Brain Stimul. 2009; 2 (e1): 201-207https://doi.org/10.1016/j.brs.2009.03.005
• Jurcak V.
• Tsuzuki D.
• Dan I.
10/20, 10/10, and 10/5 systems revisited: their validity as relative head-surface-based positioning systems.
Neuroimage. 2007 Feb 15; 34: 1600-1611https://doi.org/10.1016/j.neuroimage.2006.09.024
• Herwig U.
• Satrapi P.
• Schönfeldt-Lecuona C.
Using the international 10-20 EEG system for positioning of transcranial magnetic stimulation.
Brain Topogr. 2003 Winter; 16: 95-99https://doi.org/10.1023/b:brat.0000006333.93597.9d
• Wagner T.A.
• Zahn M.
• Grodzinsky A.J.
• Pascual-Leone A.
Three-dimensional head model simulation of transcranial magnetic stimulation.
IEEE Trans Biomed Eng. 2004 Sep; 51: 1586-1598https://doi.org/10.1109/TBME.2004.827925
• Thielscher A.
• Opitz A.
• Windhoff M.
Impact of the gyral geometry on the electric field induced by transcranial magnetic stimulation.
Neuroimage. 2011 Jan 1; 54: 234-243https://doi.org/10.1016/j.neuroimage.2010.07.061
• Makarov S.N.
• Noetscher G.M.
• Raij T.
• Nummenmaa A.
A quasi-static boundary element approach with fast multipole acceleration for high-resolution bioelectromagnetic models.
IEEE Trans Biomed Eng. 2018 Mar 7; https://doi.org/10.1109/TBME.2018.2813261
• Gomez L.J.
• Dannhauer M.
• Koponen L.M.
• Peterchev A.V.
Conditions for numerically accurate TMS electric field simulation.
Brain Stimul. 2020 Jan-Feb; 13: 157-166https://doi.org/10.1016/j.brs.2019.09.015
• Makarov S.N.
• Wartman W.A.
• Daneshzand M.
• Fujimoto K.
• Raij T.
• Nummenmaa A.
A software toolkit for TMS electric-field modeling with boundary element fast multipole method: an efficient MATLAB implementation.
J Neural Eng. 2020 Apr 1; https://doi.org/10.1088/1741-2552/ab85b3
• Daneshzand M.
• Makarov S.N.
• de Lara L.I.N.
• Guerin B.
• McNab J.
• Rosen B.R.
• Hämäläinen M.S.
• Raij T.
• Nummenmaa A.
Rapid computation of TMS-induced E-fields using a dipole-based magnetic stimulation profile approach.
Neuroimage. 2021 Aug 15; 237118097https://doi.org/10.1016/j.neuroimage.2021.118097
• Makarov S.N.
• Hämäläinen M.
• Noetscher G.M.
• Ahveninen J.
• Nummenmaa A.
Boundary element fast multipole method for enhanced modeling of neurophysiological recordings.
IEEE Trans Biomed Eng. 2020 May 29; https://doi.org/10.1109/TBME.2020.2999271
• Saturnino G.B.
• Thielscher A.