Emergence of localized persistent weakly–evanescent cortical brain wave loops

An inhomogeneous anisotropic physical model of the brain cortex is presented that predicts the emergence of non–evanescent (weakly damped) wave–like modes propagating in the thin cortex layers transverse to both the mean neural fiber direction and to the cortex spatial gradient. Although the amplitude of these modes stays below the typically observed axon spiking potential, the lifetime of these modes may significantly exceed the spiking potential inverse decay constant. Full brain numerical simulations based on parameters extracted from diffusion and structural MRI confirm the existence and extended duration of these wave modes. Contrary to the commonly agreed paradigm that the neural fibers determine the pathways for signal propagation in the brain, the signal propagation due to the cortex wave modes in the highly folded areas will exhibit no apparent correlation with the fiber directions. The results are consistent with numerous recent experimental animal and human brain studies demonstrating the existence electrostatic field activity in the form of traveling waves (including studies where neuronal connections were severed) and with wave loop induced peaks observed in EEG spectra. The localization and persistence of these cortical wave modes has significant implications in particular for neuroimaging methods that detect electromagnetic physiological activity, such as EEG and MEG, and for the understanding of brain activity in general, including mechanisms of memory.

The majority of approaches to characterizing brain dynamical behavior are based on the assumption that signal propagation along well known anatomically defined pathways, such as major neural fiber bundles, tracts or groups of axons (down to a single axon connectivity) should be sufficient to deduce the dynamical characteristics of brain activity at different spatiotemporal scales. Experimentally, data on which this assumption is employed range from high temporal resolution neural oscillations detected at low spatial resolution by EEG/MEG [1] to high spatial resolution functional MRI resting state modes oscillation detected at low temporal resolution [2]. As a consequence, a great deal of research activity is directed at the construction of connectivity maps between different brain regions (e.g., the Human Connectome Project [3]), and using those maps to study dynamical network properties with the help of different models of signal communication through this network along structurally aligned pathways [4], i.e. by allowing input from different scales or introducing axon propagation delays [5].
However, recent detection [6] of cortical wave activity spatiotemporally organized into circular wave-like patterns on the cortical surface, spanning the area not directly related to any of the structurally aligned pathways, but nevertheless persistent over hours of sleep with millisecond temporal precision, presents a formidable challenge for network theories to explain such a remarkable * vit@ucsd.edu † lfrank@ucsd.edu synchronization across a multitude of different local networks. Additionally, many studies show evidence that electrostatic field activity in animal or human hippocampus (as well as cortex) are traveling waves [7] that can affect neuronal activity by modulating the firing rates [8] and may possibly play an important functional role in diverse brain structures [9]. And perhaps more importantly it has been experimentally shown [10] that periodic activity can self-propagate by endogenous electric fields even through a physical cut in vitro that destroys all mechanisms of neuron to neuron communication.
In this paper we investigate a more general physical wave mechanism that allows cortical surface wave propagation in the cross fiber directions due to the interplay between tissue inhomogeneity and anisotropy in the thin surface cortex layer. This new mechanism has been overlooked by previous models of brain wave characterization and thus is absent from current network pathway reconstruction and analysis approaches.
The main claim of this paper is that there is a simple and elegant physical mechanism behind the existence of these cross-fiber waves that can explain the emergence and persistence of wave loops and wave propagation along the highly folded cortical regions with a relatively slow damping. The lifetime of these wave-like cortex activity events can significantly exceed the decay time of the typical axon action potential spikes and thus can provide "memory-like" response in the cortical areas generated as a result of "along-the-axon" spiky activations. This new mechanism may provide an alternative approach for the integration of microscopic brain properties and for the development of a "physical" model for memory.
The paper derives cortex wave dispersion relation from well known relatively basic physical principles and provides illustrations of why and how cortical tissue inhomogeneity and anisotropy influence propagation magnitude, time-scales, and directions and supports extended and highly structured regions of existence in dissipative media using simple 1-and 2-D anisotropic models, as well as more realistic full brain model based on a set of parameters extracted from real diffusion and structural MRI acquisitions. The wave properties (frequency ranges, phase and group velocities, possible spectra) are then compared with real EEG wave acquisitions. Finally, examples of wave propagation are studied analytically and numerically using a simple idealized but informative spherical shell cortex model (i.e. thin inhomogeneous layer around a sphere with homogeneous anisotropic conducting medium) as well as a more realistic anisotropic inhomogeneous full brain model with actual cortical fold geometry that clearly shows the emergence of localized persistent wave loops or rotating wave patterns at various scales, including at similar scales of global rotation recently detected experimentally [6].
An important aspect of the cortical waves model we present is that it is based on relatively simple but physically motivated averaged electrostatic properties of human neuronal tissue within realistic data-derived brain tissue distributions, geometries, and anisotropy. While the source of these averaged tissue properties includes the extraordinarily complex network of neuronal fiber connections supporting a multitude of underlying cellular, subcellular and extracellular processes, we demonstrate that the inclusion of such details for the activation/excitation process is not necessary to produce coherent, stable, macroscopic cortical waves. Instead, a simple and elegant physical model for wave propagation in a thin dissipative inhomogeneous and anisotropic cortical layer of these averaged properties is sufficient to predict the emergence of coherent, localized, and persistent wave loop patterns in the brain.
We will start with the most general form of description of brain electromagnetic activity using Maxwell equations in a medium Using the electrostatic potential E = −∇φ, Ohm's law J = σ · E (where σ ≡ {σ ij } is an anisotropic conductivity tensor), a linear electrostatic property for brain tissue D = εE, assuming that the permittivity is a "good" function (i.e. it does not go to zero or infinity everywhere) and taking the change of variables ∂x → ε∂x , the charge continuity equation for the spatial-temporal evolution of the potential φ can be written in terms of a permittivity scaled conductivity tensor Σ = {σ ij /ε} as where we have included a possible external source (or forcing) term F. For brain fiber tissues the conductivity tensor Σ might have significantly larger values along the fiber direction than across them. The charge continuity without forcing i.e., (F = 0) can be written in tensor notation as where repeating indices denotes summation. Simple linear wave analysis, i.e. substitution of φ ∼ exp (−i(k · r − Ωt)), gives the following complex dispersion relation which is composed of the real and imaginary components: The condition for non-(or weak-) evanescence is that the oscillatory (i.e., imaginary) component of φ, characterized by the frequency ω, is much larger than the decaying (i.e., real) component, characterized by the damping γ: i.e. that the condition |γ/ω| 1 must be satisfied. This requirement is clearly not satisfied if reported average isotropic and homogeneous parameters are used to describe brain tissues. For typical low frequency ( 10Hz) white and gray matter conductivity and permittivity (i.e. from [11]) ε GM = 4.07 · 10 7 ε 0 , ε W M = 2.76 · 10 7 ε 0 , σ GM = 2.75 · 10 −2 S/m, σ W M = 2.77 · 10 −2 S/m, where ε 0 = 8.854187817·10 −12 F/m is the vacuum permittivity) the damping rate γ is in the range of 75-115 s -1 which would be expected to give strong wave damping.
In order to better understand the effects of brain tissue micro-and macro-structure on the manifestation of propagating brain waves, it is instructive to consider two idealized tissue models. In the first model ( Fig. 1a) all brain fibers are packed in a half space aligned in z direction and their number decreases in x direction in a relatively thin layer at the boundary. We assume that small cross fiber currents can be characterized by a small parameter and introduce the conductivity tensor as where υ ≡ υ(x). For the υ(x) dependence we will assume that the conductivity is changing only through a relatively narrow layer at the boundary (as illustrated in Fig. 1a) and the conductivity gradient is directed along x axis. Then we will look for a solution for the potential φ located in the thin layer of inhomogeneity (that is we substitute x → x) that depends on z and y only as φ = φ (z) + φ ⊥ (y).
The first equation (6) describes a potential along the fiber direction and is a damped oscillator equation that has a decaying solution. But the second equation (7) describes a potential perpendicular to the fiber direction and does not include a damping term, hence it describes a pure wave-like solution that propagates in the thin layer transverse to the main fiber direction. Thus although this wave-like solution φ ⊥ has a smaller amplitude than along the fiber action potential φ , it can nevertheless have a much longer lifetime.
We would like to stress that the equations (6) and (7) are here only for illustrative purposes to reiterate a relatively obvious but often overlooked consideration that under anisotropic inhomogeneous conditions some direction may happen to be better suited for wave propagation than the others. Those equations should not be considered as the most important part as the complete dispersion relation (3) will be used to study waves propagation later.
In order to account for geometric variations, we construct a slightly more complex two dimensional model that can be viewed as a very crude approximation of a cortical fold (Fig. 1b). The equation for the φ ⊥ (y) will again be of a wave type, similar to (7), with the addition of a z component of the conductivity gradient and with a similar wave-like solution that will include an additional term induced by the inhomogeneity in z.
The sole purpose of those examples is to provide illustrations that dissipative media with complex structure may show surface wave-like solutions. Surface waves at the boundary of various elastic media have been extensively studied and used in various areas of science (including acoustics, hydrodynamics, plasma physics, etc.) since the work of Lord Rayleigh [12]. The existence of surface waves at the dissipative medium boundary is also known [13].
In order to extend the above analysis to a more realistic model of brain tissue architecture we extracted volumetric structural brain parameters from high resolution anatomical MRI datasets as well as brain fiber anisotropy from diffusion weighted MRI datasets. All anatomical and diffusion MRI datasets are from the Human Connectome Project [3]. The details for all of the processing steps can be found in [14][15][16]. More refined procedures for constructing the conductivity tensor and anisotropy in different brain regions, cortical areas in particular, would clearly be beneficial and will be addressed in the future.
To provide an illustration of where the conductivity anisotropy and inhomogeneity can form appropriate conditions for cortex surface waves generation we created plots that can be used to characterize the ratio |γ/ω|. To construct the ratio we calculated two vectors, ∂ i Σ ij and Σ ij k i and compared their norms. For wave vector k we used a vector with the same direction as in ∂ i Σ ij vector and magnitude |k| = |∇Σ|/|Σ|, i.e. our intention is to compare norms of dissipative and wave-like terms at kh ≈ 1 where h is the cortical thickness. Fig. 2 shows plot of the ratio of |Σ ij k i | and |∂ i Σ ij | at two different depths inside the brain with dissipative regime in the inner cortex (left) and wave-like conditions in the outer cortex (right).
Isovolume maps for comparison of dissipation vs wave-like effects at different cortex layers. The left panel shows the isosurface from the cortex area that may be representative of white-gray matter interface. The right panel shows the isosurface that is located in the outer cortex area of gray matter. The color scheme uses shades of green to mark regions where dissipative term dominates, i.e. |Σijki| ≥ |∂iΣij|, shades of red where |Σijki| < |∂iΣij| ≤ 2|Σijki|, and shades of blue where the wave-like term is more than two times dominant, i.e. 2|Σijki| < |∂iΣij|. The inner cortex shown in the left panel clearly display prevalence of dissipation, whereas the outer cortex shown in the right panel allows for wave-like cortex activity in a majority of the locations.
Both anisotropy and inhomogeneity are important for the existence of the cortex surface waves. For example for inhomogeneous but isotropic tissue the conductivity tensor Σ ij will simply be Sδ ij , where S is a scalar inhomogeneous conductivity. Therefore both the phase velocity (v ph = ω/k) and the group velocity (v gr = ∂ω/∂k) will include terms ∇S and ∇S · k, meaning that those waves are not able to propagate normally to the local conductivity gradient ∇S. This restriction is absent in cortex areas when both anisotropy and inhomogeneity are present.
To provide some (possibly overly optimistic) estimates based on a typical human brain dimensions we can consider a simple spherical cortex shell model with a cortical layer of fixed thickness h ≈ 1.5-3mm spread over a hemisphere of radius R ≈ 75mm, with all parameters kept constant inside the hemisphere (for r < R) and changing as a function of radius r in a cortical layer. Even without taking into account the known strong anisotropy of neural tissue, these simple geometric considerations provides for the longest waves (with the smallest amount of damping) with |γ/ω| ∼ kh ∼ 0.02-0.04. Anisotropy will reduce this estimate even further, thus further strengthening the condition necessary to support stable waves. This simple spherical shell model can be used to illustrate the natural appearance of standing-type cortical waves, which can be easily understood from simple geometrical optics arguments. Using the dispersion relation (3) with only single component of the conductivity tensor Σ zz ≡ S(r), the wave frequency ω can be expressed as where we have neglected the term Sk 2 z because kh 1. Then from the geometrical optics ray equations we can get These equations will generate rays inside the spherical cortex shell showing wave propagation across both the fibers and the conductivity gradient in the cortex subregion where ω dω/dr < 0. For k z = k/ √ 2 and z = −ωr/(dω/dr) the wave path has the simplest form -the wave follows the same loop through the cortex over and over again. Different families of frequencies ω and wavevectors k will result in the appearance of cortical wave loops at different cortex locations.
In order to determine a possible energy distribution across different frequencies for these cortical waves some knowledge about the forcing term F in (1) is required. A rough estimate for this distribution in the form of a power spectra scaling can be carried out using some simple assumptions. Assuming that the forcing consists of spiking input localized at random locations and times, it can be described as a sum of delta functions, F = frequency spectrum in the Fourier domain. Then, from the last term in (2) we can estimate the exponent α = 2 in a power law scaling of the potential φ frequency spectrum (i.e. for The presence of the cortical wave loops described above can modify this ω −2 dependence and thus modify the spectrum in such a way that spectrum peaks are generated that correspond to these loop wave currents. From (8) we can estimate the range of frequencies where those loops can possibly be present. Taking 1/r ∼ 1/R, dS/dr ∼ Σ zz /h, z = −ωr/(dω/dr) ∼ √ Rh gives a frequency estimate as f = ω/2π ∼ Σ zz /(2πk √ 2Rh), hence for the largest and smallest wavenumbers defined by the smallest (1.5mm) and largest (75mm) loop radii, the frequency f spans the range 1.2-92Hz (v ph ≈ 0.002−7m/s). The large scale circular cortical waves of [6] (9-18Hz, 2-5m/s) are clearly inside this range. The above spherical shell wave propagation example (as well as wave simulations presented below) are not able to predict the exact values for wave amplitude as it requires calculation of the balance between excitation and dissipation of these waves at various frequency ranges and the simple estimates using amplitudes of spiking are not that particularly informative as they simply predict the values that are in the range of typically observed field potentials. Nevertheless, even without the amplitude estimation our analysis is able to predict the preferred directions of wave propagation and loop pattern formation based on geometrical and tissue properties. Moreover, this framework provides the mechanisms to incorporate a more complete quantitative description of axonal spiking, which is beyond the scope of this paper but currently under investigation in our lab.
Evidence for the existence of these wave loop induced spectral peaks is shown in Fig.3, which shows the spectral power of the EEG signal for six subjects [17], averaged over all sensors. The dashed lines outline the predicted f −2 in the lower (f 1.2Hz) and higher (f 92Hz) parts of the spectra. The dashed-dotted vertical lines denote the frequency range where the cortical loops may exist, agreeing very well with typically observed EEG excessive activity range, from low frequency delta (0.5-4Hz) to high frequency gamma (25-100Hz) bands.
The effect described theoretically above can be demonstrated through numerical simulation of wave propagation in a thin dissipative inhomogeneous and anisotropic cortical layer. As a starting point for numerical study we included Fig. 4, that shows spatial snapshots of the dynamical behavior of randomly generated wave trajectories (the movies are in [18]) in a regime with γ eff L loop /v gr ≤ 1 (where γ eff is an effective wave dissipation rate, i.e. a difference between average dissipation and spiky activations rates for a wave packet propagating with group velocity v gr along some characteristic loop of L loop length. Wave packets are simulated using the ray equations (9) where the general form of anisotropic dispersion relation (3) and (4) is used.
In the idealized spherical model some of the emergent persistent localized cortical wave patterns are precisely the simple loop pattern predicted by (10), despite the complex initial spatiotemporal pattern of the initial wave trajectory. As further predicted, the simulations informed by the real human data with the same cortex fold geometry as in Fig. 2 (middle and right panels) also produce stable loop structures, now embedded within the complex geometry of the cortical folds.
All wave simulations were initialized with wave packets of random parameters (frequency, wave number, location, etc), but clearly show an emergence of localized persistent closed loop patterns at different spatial and temporal scales, from scales as large as global the whole brain rotational wave activity experimentally detected in [6] to as small as the resolution used for the cortical layer thickness detection.
In conclusion, in this paper we have presented an inhomogeneous anisotropic physical model of wave propagation in the brain cortex. The model predicts that in addition to the well-known damped oscillator-like wave activity in brain fibers, there is another class of brain waves that are not directly related to major fibers, but instead propagate perpendicular to the fibers along the highly folded cortical regions in a weakly-evanescent manner that results in their persistence on time scales long compared to waves along brain fibers. Thus waves can potentially propagate in any direction, including the direction along the fibers or in the direction of the inhomogeneity gradient. However, the dissipation of the waves is smallest when they propagate cross fibers, therefore, on average the cross fiber direction of propagation should be seen more often. For the first time, we have obtained the dispersion relation for those surface cortex waves, and have shown, both analytically and numerically, a plausible argument for their existence. Through numerical analysis we have developed a procedure to generate an inhomogeneous and anisotropic distribution of conductivity tensors using anatomical and diffusion brain MRI data. While the detailed numerical studies of effects and importance of these waves and their possible biological role are out of the scope of this paper, we have presented preliminary results that suggest that the time of life for these wave-like cortex activity events may significantly exceed the decay time of the typical axon action potential spikes. Thus, they can provide an persistent neuronal response in the cortical areas generated as a result of "along-the-axon" spiky activations.
The interaction of the traveling waves with spiking activity is of course an interesting and important question. Our model supports the inclusion of spiking activity as a source term to the wave model (one example is sketched in appendix D. In particular this may allow the development of models of brain rhythm generation based on coupling with spiking sources. Exploration of the interplay of the traveling waves and spiking activity with this model will certainly be the focus of future work.
The ranges of parameters for the waves produced by our model are in agreement with presented by several studies [7] that present evidence that electrostatic field activity in several areas of animal and human brains are traveling waves that can affect neuronal activity by modulating the firing rates [8] and may possibly play an important functional role in diverse brain structures [9]. The natural self organization of these traveling waves into loop-like structures that our model produces agrees well with recently detected [6] cortical wave ac-tivity spatiotemporally organized into circular wave-like patterns on the cortical surface. Self-propagation of endogenous electric fields through a physical cut in vitro when all mechanisms of neuron to neuron communication has been destroyed [10] can potentially be attributed to these waves as well. We have also demonstrated that the peaks these wave loops would induce in EEG spectra are consistent with those typically observed EEG data.
Direct experimental results have shown that even despite the small amplitudes of the external field potentials relative to the threshold of a spiking neuron, external fields can play a substantial role in the spiking activity and "even very small and slowly changing fields that triggered V e changes under 0.2 mV led to phase locking of spikes to the external field and to a greatly enhanced spike-field synchrony" [19]. Therefore our models ability to predict regions of cortex where external wave activity can emerge and form a sustained loop pattern has the potential to be important for understanding where the neuron spiking synchronization will have better chances to be achieved, hence it can provide substantial input in understanding effects on neural information processing and plasticity.
The ability of this new physical model for the generation, propagation, and maintenance of brain waves may have significant implications for the analysis of electrophysiological brain recording and for current theories about human brain function. Furthermore, the dependence of these waves on brain geometry, such as cortical thickness, has potentially significant implications for understanding brain function in abnormal states, such as Alzheimer's Disease, where cortical thickness changes are evident, and the dependence on tissue status may be important in conditions such as Traumatic Brain Injury, where tissue damage may alter its anisotropic and inhomogeneous properties. We first provide details about the procedures used for generating inhomogeneous and anisotropic components of the permittivity scaled conductivity tensor Σ. Spherical shell cortex model is represented by fixed anisotropy tensor σ ij (r) ≡ σ ij scaled by a radial inhomogeneous density ρ(r) ≡ ρ(r), such that the total conductivity tensor Σ ij (r) is defined by where σ ij is a constant, anisotropic, positive semidefinite symmetric tensor and ρ(r) is piecewise continuous function of radius r (0 ≤ r ≤ 1) The spherical shell cortex wave simulations include three different anisotropic tensor σ ij choices, where (for ε = 0) σ (1) represents the conductivity tensor used in analytical solution of equations (10) of the paper, i.e. currents only in z direction, σ (2) represents different orientation, where currents are allowed in the direction of 45 • relative to all axis, and σ (3) represents more complicated current anisotropy, with crossing currents flowing in x and z direction, but with no currents in y direction. Parameters n, r 0 and r 1 used to control the thickness of the inhomogeneous "cortical" layer (r 0 =0.5, r 1 =0.9, α ∞ =500 and n =0,2,14 and 25 used for Figs. A1 to A4).

Appendix B: Cortical fold model
For cortical fold model the conductivity tensor Σ ij (r) is again defined as a product of inhomogeneous ρ(r) and anisotropic σ ij (r) parts (A1), but both parts are now functions of brain locations.

Inhomogeneity estimation
The inhomogeneous density function ρ(r) is estimated from high resolution anatomical MRI data by processing it with SWD [14] (skull stripping, field of view normalization, noise filtering) and then registering to MNI152 space [20,21] with SYM-REG [16]. The final 1mm 3 182x218x182 volume is used as the inhomogeneous density ρ(r).

Anisotropy estimation
For the anisotropy tensor σ ij (r) several different test cases were employed.

a. Fixed anisotropy orientation and value
The same σ (1) , σ (2) and σ (3) fixed (location independent) anisotropic tensors as in the spherical shell cortex model (A2) were assigned to every location in the brain. This is the simplest case that may allow to separate the effects of inhomogeneity and anisotropy on loop pattern formation.

b. Varying anisotropy orientation and fixed anisotropy value
Multiple diffusion direction and diffusion gradient strength MRI data were used to estimate diffusion tensor D ij [22].
The anisotropy orientation was estimated using eigenvector d (1) of the diffusion tensor D ij with the largest eigenvalue λ (1) d . The anisotropy tensor σ ij (r) was defined as where the value of anisotropy is constant across the volume (ε = 0, 0.01 and 0.1 were used for different test examples). R is a rotation matrix between directions (0,0,1) and d (1) ≡ (x, y, z) that can be expressed as

. Varying anisotropy orientation and value
Similarly to the previous case the diffusion tensor D ij was estimated for each voxel and eigenvector d (1) with the largest eigenvalue λ (1) d was used to define the anisotropy tensor major axis. The position dependent anisotropy value was defined by intoducing parameter ε as either λ (in both cases resulting in axisymmetric form of the conductivity tensor) and also using two different values λ The microstructure anisotropy at the level of a single cell and its extracellular vicinity may be significantly higher than detected by diffusion MRI estimates. Recent measurements of effective impedance that involve both intercellular and extracellular electrodes [23] show significantly lower conductivity (by orders of magnitude) than conductivity obtained by measurements between extracellular electrodes only. Although attributing this to lower than typically assumed conductivity (or higher impedance) of extracellular medium may be questionable [24], this clearly confirms high anisotropy for the effective conductivity that should be used for analysis of effects of intercellular and membrane sources (and these are the most important type of sources for generation of surface brain waves considered in this paper). This may be important for future generation and refinement of the macroscopic conductivity estimates from microstructure data.
We assumed here that anisotropic form of σ (1) ij with ε = 0.001, 0.01 or 0.1 can be used to describe different levels of intercellular conductivity as well as microstructure extracellular conductivity in the vicinity of cell membranes (where z-direction corresponds to the direction of the fiber). Using full brain tractography results [25] we generated the anisotropic part of the conductivity tensor σ(r) in voxel at r location as an average over all fiber orientations assuming N fibers inside the voxel with R k orientation matrix for every fiber k, i.e.
with the complete form of the conductivity tensor Σ again given by (A1).

Appendix C: Wave trajectory integration
To obtain the trajectory of a brain wave with frequency ω emitted at a point r 0 with an initial wave vector k 0 we substitute the inhomogeneous anisotropic conductivity tensor Σ in the wave dispersion relation (imaginary part of eq. (3) of the paper) and obtain wave ray equations as We integrate these equations starting at τ =0, t = 0, r = r 0 and k = k 0 to trace the wave trajectory t(τ ), r(τ ) and k(τ ). For numerical integration it is beneficial to split k into magnitude k and directionk parts (k = kk, and |k| 2 = 1), and rewrite the last two equations as where in the last equation the right hand side is orthogonal tok to guarantee that |k| 2 = 1, and the algebraic expression (C1) was substituted for the differential equation for wave vector magnitude k.

Appendix D: Wave dissipation and excitation
An integral along the wave trajectory t(τ ), r(τ ) and k(τ ) with any appropriate model form of wave excitation γ exc describes a change of wave energy W along its path. This may allow the study of many interesting questions of brain wave dynamics, i.e. to identify potential active area where wave intensity may grow as a result of certain frequency and/or spatial distributions of neuronal spiky activation, or to estimate the spatial extent of coherent activation area as a result of some particular point sources, or to find when and/or how these waves may potentially become important in triggering neuronal firing, that is to study possible mechanisms of synchronization and feedback, etc.

Appendix E: Persistent wave loop patterns
One of the interesting questions is the possibility of persistent pattern formation/selforganization from a randomly emitted distribution of brain waves influenced by inhomogeous anisotropic structure of their dispersion. Considering a simplest case of fixed difference between wave excitation and dissipation, i.e. assuming propagation of wave packets of fixed width (exponentially decaying) and searching for any closed parts of wave trajectories (loops) that are shorter than the packet width can provide a clue about possible answer. Fig. 4 as well as Figs. A1 to A16 (and hyperlinked movies) show formation of persistent loops for a variety of wave packet initial conditions as well as for different models of conductivity tensor inhomogeneity and anisotropy constructed from dMRI or whole brain tractography.
Examples of wave trajectories and emergent persistent loop patterns for the spherical shell cortex model with varying amounts of tensor anisotropy and inhomogeneous shell layer thickness (Figs. A1 to A4).
Examples of wave trajectories and emergent persistent loop patterns for cortical fold geometry with different approaches used for estimation of inhomogeneity and anisotropy are shown in Figs.A5 to A16. Among those examples are several simple cases with variable inhomogeneity and fixed anisotropy (similar to the above spherical shell cortex model) as well as with more complex estimates of anisotropy based on multiple diffusion gradients MRI (dMRI) acquisitions. Assuming linear dependence between diffusion and conductivity tensors [26] several combination of the diffusion tensor eigenvectors were used to describe the conductivity tensor anisotropy. In even more complex approach, the full brain tractography [15] generated fiber distributions [25] were used to infer anisotropy through direct integration of single fiber anisotropy parameters in each voxel.