# Network-driven anomalous transport is a fundamental component of brain microvascular dysfunction

Florian Goirand, Tanguy Le Borgne, Sylvie Lorthois

► To cite this version:

Florian Goirand, Tanguy Le Borgne, Sylvie Lorthois. Network-driven anomalous transport is a fundamental component of brain microvascular dysfunction. Nature Communications, 2021, 12 (1), pp.7295. 10.1038/s41467-021-27534-8. insu-03557344

HAL Id: insu-03557344

https://insu.hal.science/insu-03557344

Submitted on 4 Feb 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

ARTICLE

https://doi.org/10.1038/s41467-021-27534-8

OPEN

<figure><img src="image_3.png" /></figure>

Check for updates

# Network-driven anomalous transport is a fundamental component of brain microvascular dysfunction

Florian Goirand<sup>1,2</sup>, Tanguy Le Borgne<sup>1✉</sup> & Sylvie Lorthois<sup>2✉</sup>

<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

Blood microcirculation supplies neurons with oxygen and nutrients, and contributes to clearing their neurotoxic waste, through a dense capillary network connected to larger tree-like vessels. This complex microvascular architecture results in highly heterogeneous blood flow and travel time distributions, whose origin and consequences on brain pathophysiology are poorly understood. Here, we analyze highly-resolved intracortical blood flow and transport simulations to establish the physical laws governing the macroscopic transport properties in the brain micro-circulation. We show that network-driven anomalous transport leads to the emergence of critical regions, whether hypoxic or with high concentrations of amyloid-β, a waste product centrally involved in Alzheimer's Disease. We develop a Continuous-Time Random Walk theory capturing these dynamics and predicting that such critical regions appear much earlier than anticipated by current empirical models under mild hypoperfusion. These findings provide a framework for understanding and modelling the impact of microvascular dysfunction in brain diseases, including Alzheimer's Disease.

<footnote n="1">¹University of Rennes, CNRS, Géosciences Rennes, UMR 6118 Rennes, France. ²Institut de Mécanique des Fluides de Toulouse, UMR 5502, CNRS, University of Toulouse, Toulouse, France. ✉email: tanguy.le-borgne@univ-rennes1.fr; sylvie.lorthois@imft.fr</footnote>
The brain microvascular network provides an efficient, highly integrated, and dynamic infrastructure for the distribution of blood¹⁻³: it supplies oxygen, nutrients, and, if needed, drugs to all cells in the brain, and ensures the removal of their metabolic waste⁴,⁵. As the brain lacks any substantial energy reserve, it also acts as a short-term regulation system, which responds quickly and locally to the metabolic needs of neurons²,⁶. In ageing and disease, however, the progressive appearance of abnormal vessel architectures, including reduced capillary diameters or stalling, and the decrease in regulation efficiency together reduce blood flow and the availability of oxygen²,⁵,⁷⁻¹¹. This also alters the clearance of metabolic waste, including neurotoxic forms of amyloid-β centrally involved in the pathogenesis of Alzheimer's Disease (AD)⁴,⁹,¹⁰. Thus, understanding the links between the microvascular architecture, reduced blood flow, and impaired oxygen delivery and metabolic waste clearance is a key challenge to decipher the role of microvascular dysfunction in brain disease.

The microvascular architecture is structured by tree-like arterioles and venules that connect to a dense capillary network³,¹²,¹³. While this organization ensures a large surface of exchange between blood and the brain tissue, it also induces strong spatial heterogeneities of vessel flows and capillary transit times, leading to heterogeneous oxygenations¹⁴,¹⁵. Even in normal conditions, some vessels with low blood flow rates approach the hypoxic threshold¹⁶,¹⁷. These critical vessels may be particularly vulnerable to further pathological stress¹⁷, consistent with the appearance of small hypoxic regions in the cortex of ageing and AD mice⁹,¹⁸. Since reduced capillary flow also compromises metabolic waste clearance, critical vessels with abnormally high intravascular concentrations of amyloid-β may also be expected. Yet, it is unknown how such critical vessels may appear under normal conditions nor how they may progress in response to pathological stress, such as hypoperfusion¹⁰,¹⁹,²⁰.

In fact, the dynamics governing oxygen or amyloid-β distributions in such networks are fundamentally non-local: the solute concentration in a given vessel depends on the blood travel time from penetrating arterioles, where blood enters into the brain cortex, to this vessel, which integrates all blood velocities along the corresponding pathways. Furthermore, the distribution of blood flow is also non-local, i.e., driven by the whole vascular architecture. Because the impact of such non-local dynamics on relevant network scale processes is difficult to resolve explicitly, the blood travel time distribution through the microvascular network has been represented by phenomenological models²¹, including early indicator dilution analysis models²¹⁻²⁴ and the recent Capillary Transit Time Heterogeneity (CTH) model¹⁴. These models rely on empirical blood travel time distributions, following mathematical functions chosen to match the experimentally observed distributions²⁵. A key property emerging from CTH models is that increased transit time heterogeneities induce a decreased efficiency of oxygen supply¹⁴. Consistently, reduced transit time heterogeneity has been experimentally confirmed in the cortical layers with the highest levels of metabolic activity²⁶ or in response to neuronal activation²⁷ while increased transit time heterogeneity has been inferred, based on such models, from clinical imaging data in AD patients²⁸. Besides these few examples, CTH models have been used to interpret a large amount of experimental data¹⁵,¹⁸,²⁹⁻³¹ thus helping to identify increased vascular heterogeneity as a key general mechanism of neuronal injury.

A fundamental limitation of current phenomenological models is that they do not quantitatively relate the transport dynamics to the underlying network architecture and flow distributions. This makes it difficult to understand and predict how changes in vessel architecture may influence blood travel time heterogeneity and

thus alter oxygen supply and metabolic waste clearance. Furthermore, in vivo measurements are limited to time scales smaller than the blood recirculation time (~5 s)³², which limits the range available for the calibration of empirical models. Hence, their predictive power for the statistics of longer time scales, which likely control the appearance of vessels with critical oxygen or amyloid-β concentrations, is strongly dependent on the mathematical functions chosen to parametrize travel time distributions²²,²³. Yet, the physical mechanisms shaping these distributions and how they depend on the network structure²¹ are poorly understood.

Theoretical analyses of transport in model random networks have shown that these systems can exhibit anomalous transport dynamics, i.e., characterized by slow power-law decays of the large travel time probabilities³³. The latter has been successfully described by Continuous Time Random Walk (CTRW) theories, providing analytical expressions of travel time distributions as a function of the microscale structures and flow distributions³³⁻³⁶. Although microvascular networks fundamentally differ from such random networks, their complex structure potentially contains the fundamental ingredients for anomalous transport dynamics to develop. To explore this hypothesis, we use highly resolved simulations of blood flow in such networks, validated from in vivo measurements, that provide access to the full statistics of blood flow and transport dynamics in realistic microvascular networks. We use these insights to uncover the scaling laws of blood flow distributions arising from the microvascular architecture and develop an effective transport model at the scale of the network that captures these properties. This provides analytical solutions for the blood travel time distributions inferred from the physics of transport in these networks. Our model predicts that the interplay between the spatial distribution of arterioles and venules and the mesh-like architecture of the capillary bed¹² leads to the emergence of anomalous transport dynamics. This implies that the occurrence probability of large blood travel times is significantly larger than predicted by current models. We couple this model to the kinetics of oxygen consumption and amyloid-β production in brain cells to show that these anomalous transport properties control the early development of critical vessels with low oxygen or large amyloid-β concentrations under hypoperfusion. These findings hence provide a framework to measure, understand and model the onset and development of brain diseases, such as AD³⁷.

## Results

Blood travel times through microvascular networks, i.e., from any network inlet to any network outlet, are also referred to as transversal times²⁴ or, simply, as transit times¹⁴, whereas the distribution of travel times is sometimes referred to as the impulse response function, as the transfer function of the system or as the microvascular response function³⁸. To avoid any ambiguity, we use below the terms travel times for transport across the whole network and transit times for transport across a single vessel (as defined by Eq. (2)).

## Vessel flow rates and transit times follow broad distributions.

Our analysis is based on highly resolved simulations of blood flow in anatomic microvascular networks, validated by comparison with in vivo measurements (see Methods and Supplementary Note 1). We first present the results obtained in a microvessel network (~15,000 vessels) digitized from a 1 mm³ of the mouse cortex (see Supplementary Figure S1) and then compare the results to another mouse intracortical microvessel network and to bio-mimetic networks (Supplementary Notes 3 and 8). The dense capillary bed in the first network is homogeneous and space-
<figure><img src="image_7.png" /><figcaption>Fig. 1 The blood flow organization in cortical microvessel networks induces broadly distributed flow rates and vessel transit times. a Simulated blood flow map in a microvascular network feeding 1mm³ of the mouse cortex (see Supplementary Figure S1a, b and Supplementary Movie S1). Vessel flow rates are represented with blue shades for $Q < Q_c$ and red shades for $Q > Q_c$ (see also Supplementary Figure S1c, d, f and Supplementary Movie S2). b PDF of simulated flow rates (black dots) compared with the approximation of Eq. (1) (continuous blue line). The characteristic flow rate $Q_c$ separating the uniform and power-law regimes are indicated as a dotted line. Inset: PDF of vessel transit times. The diffusion coefficients for oxygen and amyloid-$\beta$ ($D_{O_2} = 2 \times 10^{-9}$ m² s⁻¹ and $D_{A\beta} = 6 \times 10^{-11}$ m² s⁻¹) yield diffusion times ($t_{D}^{O_2}$ and $t_{D}^{A\beta}$) indicated by the green and red dashed lines, respectively. c Example of trajectories visiting <30 vessels (orange) and >70 vessels (blue), originating from the arteriole shown by the arrow (see Supplementary Figure S2a-d and Supplementary Movies S3).</figcaption></figure>

filling¹⁷,³⁹, with a narrow distribution of vessel diameters (4.8 ± 0.9 µm). The network is fed and drained by ~15 arteriolar and ~30 venular trees (Fig. 1a, Supplementary Figure S1 and Supplementary Movie S1). Simulations integrate the non-linear blood rheology and red blood cell repartition at diverging vessel bifurcations (see Methods). The average flow rate is $\langle Q \rangle_0 = 4 \times 10^{-2}$ nL/s. The probability density function (PDF) of blood flow rates in vessels shows a large spread across about seven decades around the average flow rate (Fig. 1b). This PDF exhibits two different regimes: it is uniform in the low flow rate range and decays as a power law with exponent -2 above a characteristic value $Q_c = 10^{-2}$ nL/s. The number of vessels in these regimes is approximately $P(Q < Q_c) = 60\%$ and $P(Q > Q_c) = 40\%$, with large flow rates developing preferentially in the neighborhood of arterioles and venules (Fig. 1a and Supplementary Figure S1f, Supplementary Movie S2).

The flow rate PDF is well approximated by a Cauchy distribution (Fig. 1b), consistent with the theoretical asymptotic behavior obtained in large Random Regular Graphs⁴⁰.

$$p_Q(Q) = \frac{2}{\pi Q_c} \frac{1}{1 + (Q/Q_c)^2}, \quad (1)$$

A direct consequence of the broad distribution of flow rates Q is the broad distribution of advective transit times in vessels (inset of Fig. 1b), which are defined as:

$$t = \frac{l\pi d^2}{4Q}, \quad (2)$$

where $l$ is the vessel arc length and $d$ is the vessel diameter. As the flow rate PDF, the vessel transit time PDF is characterized by two regimes separated by the characteristic time $t_c = \frac{(l)_{\text{cap}} \pi (d)_{\text{cap}}^2}{4Q_c} = 10^{-1}$ s, where $\langle l \rangle_{\text{cap}} = 50$ µm and $\langle d \rangle_{\text{cap}} = 5$ µm are taken as characteristic capillary length and diameter, respectively. Since the vessel length $l$ and diameter $d$ vary weakly compared to the flow rate Q, the transit time variability is mainly driven by the flow rate fluctuations. We thus estimate the PDF of vessel transit times from the PDF of flow rates as $p_t(t) = p_Q(Q)dQ/dt$, also yielding a Cauchy distribution

$$p_t(t) = \frac{2}{\pi t_c} \frac{1}{1 + (t/t_c)^2}, \quad (3)$$

in good agreement with the simulations (inset of Fig. 1b). The

scaling $p_t(t) \sim t^{-2}$ for long times, $t > t_c$ is induced by the uniform flow rate PDF $p_Q(Q) \sim \text{cst}$ at low flow rates, $Q < Q_c$, leading to $p_t(t) \sim dQ/dt \sim t^{-2}$ (Eq. (2)). Equation (3) implies that the probability of a vessel belonging to this power-law regime is $p(t>t_c) = \int_{t_c}^{\infty} dt p_t(t) = 1 - 2/\pi \tan^{-1}(1) = 0.5$. Hence, half of the vessels belong to this regime. Since our anatomical networks have about $N = 15,000$ vessels, the lowest probability that we can measure is $1/N = 7e^{-5}$. This corresponds to a maximum transit time $t_m = t_c \tan(\pi/2(1 - 1/N)) \approx 10^4 t_c$. Therefore, the size of the network allows observing this power-law behavior across up to four orders of magnitudes in transit times from $t_c$ to $10^4 t_c$ (see inset of Fig. 1b). This power-law scaling yields a non-negligible probability of extremely long vessel transit times (e.g., 1% of vessels have a transit time $t > 50 t_c = 5$ s), which may lead to the emergence of anomalous transport properties at the network scale³³,³⁴.

This vessel transit time PDF $p_t(t)$ characterizes the vessel transport statistics in the absence of diffusion. Diffusion introduces a maximum cutoff time $t_D = \langle l \rangle^{2/D}$ in the transit time PDF, corresponding to the diffusive transport time over a vessel length. For oxygen and amyloid-$\beta$, the two species considered here, the range of times $t_c < t < t_D$ over which the power-law $p_t(t) \sim t^{-2}$ holds thus covers respectively one and three orders of magnitude (Fig. 1b, inset). Therefore, this power-law regime affects significantly more amyloid-$\beta$ clearance than oxygen supply.

**Network trajectory lengths and travel times show anomalous transport statistics.** In addition to the broad distribution of vessel transit times, solute transport at the scale of the microvascular network is also controlled by the distribution of trajectory lengths from arterioles to venules. Our particle tracking simulations (see Methods, Fig. 1c and 2a and Supplementary Movie S4) show that trajectories lengths $L$, expressed in the number of visited vessels, vary from <10 to ~80 (Fig. 2b). The trajectory length PDF is characterized by a power-law scaling $p_L(L) \sim L^{-2}$ between two characteristic lengths $L_0 = 12$ and $L_c = 50$. Above $L_c$, the PDF decays sharply as $p_L(L) \sim \exp(-L/L^*)$ with $L^* = 5$ (Fig. 2b):

$$\begin{cases} p_L(L) \sim L^{-2} & L_0 < L \le L_c \\ p_L(L) \sim \exp(-L/L^*) & L > L_c \end{cases} \quad (4)$$

The scale $L^*$ characterizes the exponential decay of the PDF of trajectory lengths at large lengths and is estimated by fitting Eq.
<figure><img src="image_8.png" /><figcaption>**Fig. 2** Trajectory lengths and travel times are broadly distributed at network scale. a 3D view of typical particle trajectories (see also Supplementary Movie S4) with growing numbers of visited vessels (Red: L = 20; Yellow: L = 30; Orange: L = 40; Green: L = 50; Cyan: L = 60; Blue: L = 70; Violet: L = 80). b Probability Density Function (PDF) of trajectory lengths. Inset: Local vessel transit times $\bar{t}^L$ averaged over all trajectories of same length L as a function of the number n of vessels visited since the inlet arteriole, normalized by L. Results for different trajectory lengths are shown with the same color conventions as in a). c PDF of travel times through the network for purely advective (black dots), oxygen (green dots), and amyloid-β (red dots) transport. Inset: Average travel time $\bar{T}^L$ as a function of trajectory length L (see also Supplementary Note 4 and Supplementary Figure S4a, b). The linear and power-law tendencies (Eq. (5)) are shown respectively as dashed and dotted blue lines. Note that the linear tendency does not appear as a straight line in logarithmic coordinates because of the constant $T_0$ in Eq. (5) (See Supplementary Note 4).</figcaption></figure>

(4) to the simulated PDF (Fig. 2b). This cutoff in trajectory lengths is likely controlled by a complex interplay between the topology of the network and the density of venules/arterioles, which prevents arbitrary large trajectories to develop.

Averaging over multiple trajectories of equal length, we computed the evolution of the local average transit time along the trajectory, as a function of the number n of vessels visited since the inlet arteriole (Fig. 2b, inset). These average local transit times $\bar{t}^L$ are of the order of $10^{-2}$ s close to the inlet arterioles, then increase, up to two orders of magnitudes for the longest trajectories, then decrease again in the vicinity of the outlet venules. The network travel time averaged over all trajectories of equal length $\bar{T}^L$ thus increases with the trajectory length L following two different trends (Inset of Fig. 2c, Supplementary Note 4 and Supplementary Figure S4a, b),

$$ \left\{ \begin{array}{l} \bar{T}^L(L) \approx \tau_1(L - L_0) + T_0 \quad L_0 < L \le L_c \\ \bar{T}^L(L) \approx T_c \left(\frac{L}{L_c}\right)^4 \quad L > L_c \end{array} \right. \qquad (5) $$

where $T_0 = 0.2$ s, $T_c = 1.8$ s, and $\tau_1 = 0.04$ s is on the order of the mean transit time over particle trajectories.

The broad distribution of vessel transit times, together with the distribution of trajectory lengths, leads to a broad range of travel times at the network scale $p_T(T)$ (Fig. 2c). Without diffusion, simulated travel times vary over four orders of magnitudes. Accounting for the diffusive cutoff transit time (inset of Fig. 1b), the travel time distribution still covers over two orders of magnitude for oxygen and three orders of magnitude for amyloid-β (Fig. 2c).

**Flow and transport properties emerge from the physics of dipole flows in networks.** The statistical properties driving these transport dynamics can be understood as arising from different topological properties of the flow field, as schematized in Fig. 3. The measured scaling $p_Q(Q) \sim Q^{-2}$ is a characteristic of dipole flows$^{41}$ (see Supplementary Note 3, Supplementary Figure S3d). In the present system, high flow rates are localized around arterioles and venules that act as multiple sources and sinks (Fig. 1a and Supplementary Figure S1, Supplementary Movie S1), driving the flow in the network. Hence, in the large flow range, $Q > Q_c$, the flow field behaves statistically as a superposition of

<figure><img src="image_9.png" /><figcaption>**Fig. 3** Schematic illustration of the blood flow organization in the cortical microcirculation inferred from the observed flow and transport properties. Arterioles and venules are represented in red and black respectively. Blood flow paths through the capillary bed are represented as light purple lines. The characteristic trajectory length $L_c$ is the maximum length of trajectories that can travel directly from one arteriole to a neighboring venule (green line).</figcaption></figure>

dipoles. At low flow rates, $Q < Q_c$, the flow rate statistics differ from the continuous dipole model and become uniform, a signature of the network structure (see Supplementary Note 3, Supplementary Figure S3d). In random networks, the flow rate PDF is indeed theoretically expected to follow exponential distributions, driven by the random additions and divisions of flow at vertices$^{42}$. In the low flow rate range, this asymptotically leads to a uniform distribution of flows. Here, the low flow vessels are far from arterioles and venules. Therefore, their statistics are dominated by the random fluctuations induced by the network structure. Noteworthy, similar flow rate PDFs are observed for dipole flows on space-filling networks with architectures of increasing complexity, see Supplementary Figure S3b-d. Thus, the observed flow rate statistics are generic and arise from the interplay between a structured dipole-like topology for large flow rates and random network topology for small flow rates.
<figure><img src="image_10.png" /><figcaption>**Fig. 4** Long travel times, captured by our mechanistic model, drive the emergence of critical areas with impaired oxygen delivery and amyloid-β clearance under hypoperfusion. a CTRW model predictions (full black line) compared to purely advective, amyloid-β, and oxygen transport simulations (black, red, and green dots, respectively). The prediction of the mean-field transport model (Eq. (6)) and of the reference CTH model (Supplementary Note 7) are shown as dashed black and blue lines, respectively. b Fraction $f_p$ of network travel times larger than $\tau_c^{O_2}$ for oxygen delivery as a function of average flow rate $\langle Q \rangle / \langle Q \rangle_0$. Simulations (green dots) are compared to the predictions of the mean-field model (dashed black line) and of the reference CTH model (dashed blue line). The inset shows the fraction $f_v$ of vessels in the network that are only reached by flow paths with travel times to these vessels larger than $\tau_c^{O_2}$, as a function of the flow rate. c Fraction $f_p$ of network travel times larger than $\tau_c^{A\beta} = 8$ s (squares), $\tau_c^{A\beta} = 16$ s (dots) and $\tau_c^{A\beta} = 40$ s (triangles), as a function of average flow rate for amyloid-β clearance. The predictions of the CTRW model (Eq. (7)) for each value of $\tau_c^{A\beta}$ are shown as continuous black lines. The predictions of the reference CTH model (see Supplementary Note 7) for each value of $\tau_c^{A\beta}$ are shown as blue dashed lines.</figcaption></figure>

The distribution of trajectory lengths (Fig. 2b), is also consistent with this dipole flow interpretation (see Supplementary Note 3, Supplementary Figure S3e). For dipole flows on a finite-size network, the power-law regime $p_L(L) \sim L^{-2}$ develops until a cutoff length, corresponding to the network size, which sets the maximum trajectory length. In the present system, the characteristic trajectory length $L_c$ that sets the transition to the exponential cutoff (Eq. (4)) is the maximum length of trajectories that can travel directly from one arteriole to a neighboring venule<sup>43</sup> (Fig. 1c and Supplementary Figure S2a, b, Supplementary Movie S3 (orange trajectories)). These trajectories typically reach the bottom of the simulated domain and include on average 20 arteriolar steps, 20 venular steps, and 10 capillary steps in between (see Fig. 3 and Supplementary Note 2, Supplementary Figure S2f). Longer trajectories connect more distant arterioles and venules via the deep capillary bed, where most vessels are in the low flow regime ($Q < Q_c$) (Figs. 1c, 3 and Supplementary Figure S2c, d, Supplementary Movie S3 (blue trajectories)).

The relationship between average travel times and trajectory lengths $\overline{T}^L(L)$ (Eq. (5)) shows a transition from a linear to a power-law scaling (Inset of Fig. 2c and Supplementary Figure S4a, b). In the first linear regime, the number of visited capillaries remains approximately constant when the trajectory length increases (see Supplementary Note 2 and Supplementary Figure S2f). This implies that the additional visited vessels belong to arterioles and venules. The linear scaling of $\overline{T}^L(L)$ indicates that the average transit time $t$ in such vessels remains approximately constant when trajectories explore deeper sections of the network. As arterioles penetrate at depth, their flow rate decreases but so does their diameter and length, which may explain this finding. This may constitute an evolutionary advantage contributing to spatially uniformizing the supply and clearance of solutes across the network while minimizing total dissipation and blood volume<sup>44–46</sup>. In the second regime, the power-law scaling $\overline{T}^L(L) \sim L^4$ is characteristic of dipole flow in 3D systems (see Supplementary Note 3): in the deep network, blood flows from the pre-capillary arterioles to the post-capillary venules through a complex 3D network of capillaries (Fig. 1c, 3, Supplementary Figure S2c, d and Supplementary Movie S3 (blue trajectories)) where the number of visited capillaries increases steeply with the trajectory length (Supplementary Figure S2f).

**Mean-field transport dynamics are governed by dipolar trajectory length distributions.** The distribution of trajectory lengths (Eq. (4)), coupled with the relationship between average travel time and trajectory length (Eq. (5)), provides a mean-field transport model for the travel time PDF across the network (Supplementary Note 4). This mean-field description, which neglects random fluctuations due to the network structure but captures the dipole-driven trajectory length distribution, is characterized by a transition from power law to a stretched exponential behavior:

$$ 
\begin{align} 
p_{\overline{T}}(\overline{T}) &\sim ((\overline{T} - T_0)/\tau_1 + L_0)^{-2} \nonumber \\ 
&\quad \text{for } \overline{T} < T_c \nonumber \\ 
p_{\overline{T}}(\overline{T}) &\sim \left(\frac{\overline{T}}{T_c}\right)^{-3/4} \exp\left(-\frac{L_c}{L^*}\left(\frac{\overline{T}}{T_c}\right)^{1/4}\right) \tag{6} \\ 
&\quad \text{for } \overline{T} > T_c \nonumber 
\end{align} 
$$

with $T_c = 1.8$ s (Eq. (5)). This model captures the travel time distribution over the first two orders of magnitude (Fig. 4a and Supplementary Note 4, Supplementary Figure S4c). We compare this prediction to that of a reference CTH model, that assumes a Gamma distribution of travel times, whose parameters are calibrated from in vivo data<sup>14</sup> (see Supplementary Note 7). Since experimentally measured travel time distributions are limited to times smaller than ~5 s due to blood recirculation<sup>27,32</sup>, this CTH model serves here as a reference to assess the effect of neglecting the experimentally inaccessible longest travel times. While it captures relatively well the shape of the travel time distribution in the low range, the reference CTH model significantly underestimates the probability of late times (Fig. 4a). Accounting for the trajectory length distribution via the mean-field model allows capturing a significant part of this long-time dynamics. The mean-field model, however, does not capture the power-law behavior of the longest travel times, driven by vessels with transit times $t > t_c$ in the deep capillary network (Fig. 1c). Since the oxygen cutoff diffusion time $t_D^{O_2}$ is close to $t_c$ (inset of Fig. 1b), this long-time regime does not affect much oxygen transport, which is well represented by the mean-field model (Fig. 4a).

**Random flow fluctuations in the capillary network control large blood travel times.** To obtain a full description of the transport dynamics, in particular for low diffusivity solutes such
<figure><img src="image_11.png" /><figcaption>Fig. 5 Simulated distribution of blood oxygen concentration relative to oxygen arteriolar concentration in the microvascular network represented in Fig. 1. See also Supplementary Movie S5. a Side view. b Top view. c Evolution of the mean vessel oxygen concentration along trajectories of size L, concentrations are normalized by the arteriolar concentration $c_{O_2}^A$. Red: L = 20; Yellow: L = 30; Orange: L = 40; Green: L = 50; Cyan: L = 60; Blue: L = 70; Violet: L = 80.</figcaption></figure>

as amyloid-β, we seek an effective transport model that captures the heterogeneity of vessel transit times around the average behavior described above. Fluid elements transported along a given trajectory in the network move from one vessel to the next with a broadly varying time step (see Supplementary Notes 5 and 6 and Supplementary Figure S5a)). This closely corresponds to the conceptual framework of CTRW$^{34,36}$, which has been proved relevant as an effective representation of transport in random networks$^{33,35}$. However, as a consequence of the dipole structure of the flow at a large scale, two characteristics of the transport dynamics differ from a conventional CTRW representation. First, the number of steps in a given trajectory is broadly distributed according to trajectory lengths (Eq. (4)). Second, the mean local transit time at the nth step within a trajectory of length L depends on the trajectory length (Eq. (5)), inset of Fig. 2b and Supplementary Note 4). Fluctuations of local transit times around these mean local transit times define a noise term that is independent of trajectory lengths and follows a power-law scaling controlled by the flow distribution (see Supplementary Note 5). We, therefore, developed a CTRW framework capturing these two properties (see Supplementary Note 6).

This framework allows deriving an analytical solution for the travel time distribution $p_T(T)$ in the network,

$$p_T(T) = p_T^1(T) + p_T^2(T), \quad (7)$$

where $p_T^1(T)$ and $p_T^2(T)$ are defined in Laplace space as,

$$\tilde{p}_T^1(s) = e^{-sT_0} \sum_{L=L_0+1}^{L_c} p_L(L) (\mathcal{P}(s\tau_1))^{L-L_0} \quad (8)$$

and

$$\tilde{p}_T^2(s) = e^{-sT_0} \mathcal{P}^{L_c-L_0}(s\tau_1) \sum_{L=L_c+1}^{\infty} p_L(L) \mathcal{P}^{L-L_c} \left( s \frac{T_c(L/L_c)^4 - T_c}{L - L_c} \right) \quad (9)$$

where s is the Laplace variable, $\mathcal{P}(s) = 2sK_2(2\sqrt{s})$ with $K_2$ the Bessel function of the second kind and $p_L(L)$ is given by Eq. (4).

This CTRW model is fully determined from the trajectory length PDF (Eq. (4)) and the relationship between average time and trajectory length (Eq. (5)). It provides an accurate prediction over a broad range of travel times with no fitting parameter (Fig. 4a). In particular, it captures the late-time power-law decay:

$$p_T(T) \sim T^{-3}. \quad (10)$$

The emergence of this power law is consistent with the CTRW

theory that predicts stable power-law distributions for noises characterized by power-law exponents equal or smaller than 3 (see Supplementary Figure S5)$^{34,36}$. Hence, the late-time transport behavior is characterized here by stable anomalous transport, driven by the broad distribution of flow rates in the network. As far as we know, no experimental data have been obtained to confirm this result in brain microcirculation. However, injections of flow-limited tracers in the coronary network of isolated rabbit hearts, which avoids blood recirculation, yielded a power-law decay of the late-time regime with exponent $-3.21 \pm 0.27^{47}$. This suggests that, despite the variability of microvascular architecture between organs$^{48,49}$, the combination of tree-like structures with a dense capillary network in their global organization is sufficient to drive this late-time power-law decay.

**Simulation of three-dimensional oxygen distribution in microvascular networks.** Assuming that oxygen consumption in vessels follows first-order kinetics$^{25}$ (see Methods and Supplementary Note 9), we use our particle tracking simulations to compute three-dimensional fields of oxygen concentration in vessels (Fig. 5a, b and Supplementary Movie S5). The characteristic reaction time of the first-order kinetics $\tau_{T_0}^A = 1.5$ s is estimated by matching the simulated ratio between oxygen concentrations at venular outlets and those at arteriolar inlets to typical measured values (see Methods). The evolution of the oxygen concentration along trajectories, averaged over trajectories of equal length (Fig. 5c), is in qualitative agreement with previous experimental observations$^{17}$. Oxygen decays close to the network inlets, as trajectories penetrate down in the cortex, then reaches a minimum in the capillary bed (except for the shortest trajectories $L = 20$), and increases again as blood flows up to the venular outlets (Fig. 5c). The minimum oxygen value along trajectories decreases with the trajectory length L, as longer trajectories penetrate deeper in the network. This directly results from the increase of blood travel time with increasing trajectory length (Eq. (5) and inset of Fig. 2c). The minimal oxygen concentration in the network, reached in the deep capillaries of the longest trajectories, is close to 1/12 of the inlet concentration, i.e., 10 mm Hg of oxygen partial pressure, assuming 120 mm Hg at the inlets. This value is typically used to identify hypoxic brain regions in animal experiments$^{50,51}$. Hence, our results are consistent with the recent experimental observation that some vessels approach the hypoxic threshold in the cortex of normal mice$^{17}$. The increase of oxygen towards the outlet venules is also consistent with in vivo observations showing an increase of oxygen
concentrations with increasing venous diameters¹⁷. Although its origins are debated¹⁷, the two minimal ingredients included in our simulations, i.e., first-order decay of oxygen with travel time and mixing at vessel intersections, are sufficient to capture this behavior. Hence, the travel time and trajectory length statistics explored here are key to explaining the oxygen dynamics.

**Anomalous transport drives the early appearance of hypoxic regions under conditions of hypoperfusion.** Hypoperfusion, i.e., the decrease of the average blood flow, is a major pathological stress associated with many diseases¹⁹, including early stages of AD²⁰,³⁷. Schematically, it may be due to reduced perfusion pressure or increased cerebrovascular resistance, e.g., induced by capillary occlusions⁸ or reduced vessel diameters⁹ in AD. In the latter case, the total resistance of the capillary bed has been estimated to increase more than threefold, leading to a ~50% blood flow reduction. This large hypoperfusion level significantly increases tissue hypoxia⁹ and is believed to favor amyloid-β accumulation in the brain¹⁰, thus participating in disease progression. However, whether lower levels of hypoperfusion (5–30%), such as induced by capillary occlusions ⁸, may be involved in disease onset is debated¹⁰. In this context, we use our modeling framework to quantify the appearance of hypoxic regions, focusing on the occurrence probability of trajectories with travel times above the critical travel time $f_p(T > \tau_c^{O_2})$, which corresponds to the hypoxic threshold $c_{O_2} = 1/12c_{O_2}^{A, 50, 51}$. Using the first-order kinetic model (Supplementary Note 9, equation S28), we estimate this critical travel time to be $\tau_c^{O_2} = \log(12)\tau_r^{O_2} = 3.7$ s. The probability of exceeding $\tau_c^{O_2}$ along a given flow pathway is derived from our transport model as $f_p(T > \tau_c) = \int_{\tau_c}^{\infty} dT p_T^{(Q)}(T)$, where $p_T^{(Q)}(T)$ is the predicted network travel time distribution for any flow rate $\langle Q \rangle$. The latter is derived from the travel time distribution in normal conditions as $p_T^{(Q)}(T) = p_T^{(Q)0}(T \langle Q \rangle_0/\langle Q \rangle)$, where, $\langle Q \rangle_0$ is the blood flow in normal conditions, deduced from the above simulations, for which all parameters, including the perfusion pressure, correspond to physiological data (see Methods and Supplementary Note 1). In other words, changing the average flow rate $\langle Q \rangle$ is equivalent to rescaling the critical time as $\tau_c \langle Q \rangle / \langle Q \rangle_0$. This also holds true when the flow rate change is induced by capillary occlusions as the shape of the travel time PDFs is similar, not only for different sets of mouse microvascular networks, but also when up to 10% of capillaries are occluded (see Supplementary Note 8 and Supplementary Figure S6).

Under normal perfusion, the fraction of travel times larger than $\tau_c^{O_2}$ is equal to $f_p(T > \tau_c^{O_2}) = 1.3\%$ (Fig. 4b). This fraction increases non-linearly as the average flow decreases, e.g., by a factor 4 when $\langle Q \rangle$ decreases by a factor 2, an evolution accurately captured by the mean-field model (Eq. (6)). As discussed above, oxygen transport is mostly controlled by the dipole-driven transport regime and its travel time distribution is relatively well represented by the mean-field transport model (Fig. 4a). We thus approximate $p_T^{(Q)0}$ by the mean-field travel time distribution (Eq. (6)), which provides an accurate prediction of the simulated fraction of critical travel times for all flow rates (Fig. 4b). Since the latter integrates the broad distribution of trajectory lengths, it predicts a significantly larger occurrence probability $f_p(T > \tau_c^{O_2})$ than the reference CTH model (Fig. 4b), which is two orders of magnitude below the simulated probabilities under baseline perfusion. Noteworthy, the fractions of trajectories $f_p(T > T')$ measured at the outlet are very close to the fractions of critical vessels $f_v(T > T')$ within the network, i.e.,

vessels that are only visited by fluid elements taking a time larger than $T'$ to reach them from the inlet arterioles (Supplementary Note 2, Supplementary Figure S2e). Hence, the travel time statistics at the outlet venules offer a surrogate for the transport statistics within the network. As a result, the mean-field model also provides a good prediction for the probability of occurrence of hypoxic vessels within the network (inset of Fig. 4b). When hypoperfusion is induced by capillary occlusions, these hypoxic vessels may likely add up to the occluded vessels, thus enhancing their impact at the early stages of AD⁸.

**The weak diffusivity of amyloid-β amplifies the impact of anomalous transport.** We use a similar modeling method as for oxygen to relate the blood travel time statistics to amyloid-β concentrations (see Methods and Supplementary Note 10). Owing to its low diffusivity, amyloid-β is highly sensitive to random flow fluctuations in the capillary network as discussed above (Fig. 2a). Therefore we approximate $p_T^{(Q)0}$ by the CTRW model (Eqs. (7–9)). The metabolism and neurotoxicity of amyloid-β involve multiple soluble and insoluble isoforms and are still poorly understood⁵²⁻⁵⁴. To account for this uncertainty, we consider a range of critical times $\tau_c^\beta = 8, 16, \text{ or } 40$ s. This leads respectively to a three-, five- or tenfold arterio-venous increase of the total intravascular amyloid concentration (see Methods and Supplementary Note 10), which is much larger than the measured physiological increase of ~20%⁵⁵, thus yielding different degrees of compromised clearance. As expected, the probabilities of occurrence of trajectories with travel times above these critical times $f_p(T > \tau_c^\beta)$ vary significantly with $\tau_c$ (Fig. 4c). Although the CTH model predicts an exponential evolution of the fraction of critical travel times $f_p(T > \tau_c^\beta)$ with $\langle Q \rangle$ (linear in semi-log in Fig. 4c), the CTRW model captures the clearly non-exponential trend observed from the simulations. This is a signature of the power-law scalings of travel time distributions, driven by anomalous transport. Hence, the reference CTH model underestimates the fractions of travel times with inefficient amyloid-β clearance by several orders of magnitude, while the CTRW model accurately predicts the probability of these critical travel times for all $\langle Q \rangle$ and $\tau_c$.

## Discussion

By reducing the complexity of the transport problem in anatomically realistic networks while keeping the essential physics of transport emerging from the network architecture, our analysis reveals the physical mechanisms by which the microvascular architecture shapes the blood travel time distribution. This is a major fundamental open question in microvascular physiology and a bottleneck for accurate quantification of hemodynamic parameters from brain imaging data in a broad range of applications, from clinical studies aimed at improving the diagnosis and/or staging of brain disease to fundamental studies on cerebral blood flow and metabolism, neurovascular coupling, cerebral autoregulation and/or blood-brain barrier function in health and disease²¹.

We have demonstrated that the blood travel time distributions are driven by two fundamental mechanisms constitutive of anomalous transport dynamics³⁴: broadly distributed blood trajectory length and broadly distributed capillary transit times. As schematized in Fig. 3, the former is determined by dipolar flow patterns resulting from the localized connections with upstream and downstream surface vessels, while the latter is driven by random-like fluctuations within the capillary network. For high diffusivity species, such as oxygen, travel time distributions are cutoff by diffusive transport, a third fundamental mechanism that
dampens the random-like fluctuations. These fundamental insights yield the first physics-based analytical solutions for transport at the scale of cortical microvascular networks, accurately predicting the statistical distributions of travel times of different solutes, including oxygen (Eq. (6) and amyloid-β (Eqs. (7–9)).

This offers an alternative to current effective models of transport at the scale of microvascular networks. Such models either consider simplified networks as combinations of parallel elements (e.g., ⁵⁶), overlooking the link between the transport dynamics and the underlying microvascular architecture, or use empirical functions (e.g., ¹⁴). The latter are calibrated from in vivo measurements, which, due to blood recirculation, are limited to travel times of 5s. Extrapolation to larger time scales is constrained by the underlying parametrization chosen for the distribution of travel times, most often gamma distributions¹⁴,²⁵,⁵⁷. The resulting exponential decay underestimates the probability of large travel times that follow a slow power-law decay induced by the anomalous transport dynamics uncovered here. Such models predict that a significant level of hypoperfusion, where blood flow is reduced by ~20% compared with normal perfusion, should be reached for hypoxic vessels to appear. Even in normal conditions, we found that a small proportion (~1.2%) of hypoxic vessels develop in the microvascular network. This finding is consistent with experimental measurements¹⁶,¹⁷. Furthermore, accounting for anomalous transport leads to a regular non-linear increase of the proportion of hypoxic vessels with the decreased flow. The impact of anomalous transport is stronger for amyloid-β, which has a smaller diffusion coefficient. Hence, the probabilities of occurrence of critical vessels with inefficient metabolic waste clearance is orders of magnitude larger than predicted by current empirical models under normal conditions and their increase under hypoperfusion occurs much earlier than anticipated by these models.

Overall, these findings underpin the physical mechanisms by which moderate levels of hypoperfusion yield a non-linear expansion of hypoxic regions and lead to increased accumulation of amyloid-β in the brain tissue. Crucially, hypoxia and amyloid-β accumulation are two key ingredients of the amyloid cascade, the positive feedback loop linking hypoperfusion and amyloid-related pathways of AD (see, e.g., Fig. 3 in ref. ¹⁰). Combined with the recent discovery that, at early stages of AD, every single capillary occlusion has a similar, and cumulative, impact on blood flow, without any threshold effect⁸, this suggests that capillary occlusions, even in small proportions, may trigger this positive feedback loop, leading to pericyte activation and capillary constrictions. Current numerical simulations⁸ neglecting these pathological stresses as well as any compensatory effects, e.g., variations of the cerebral perfusion pressure, tend to underestimate the magnitude of hypoperfusion compared to experiments in animal models of AD. Hence, we speculate that pathological stresses induced by the amyloid cascade have an outsized impact on disease progression. Our results suggest that such cascading events may be initiated by spatial heterogeneities of blood flow. In particular, the critical regions with abnormally high amyloid-β concentration, which can be interpreted as hotspots for amyloid species accumulation, are deterministically located in regions fed by long trajectories. This could possibly explain why no association between capillaries with abnormal/stalled flow and amyloid deposition has been found experimentally⁸. Furthermore, as larger travel times are expected for trajectories feeding the sub-cortical regions, this may explain their specific vulnerability⁵⁸,⁵⁹, as highlighted by the appearance of white matter hyper-intensities in clinical imaging, in both cerebrovascular disease and AD²⁸,⁶⁰. This suggests that, while the capillary bed is a continuously adapting

network, which remodels to reduce the extent of hypoxic regions⁶¹, strong constraints may be imposed by the architecture of the vascular network and the associated locations of inlets and outlets. Besides shedding light on the impact of capillary occlusions at early stages of AD, and more generally on the role of hypoperfusion in AD onset and progression, these findings contribute to explaining the considerable overlap between vascular and neurodegenerative factors in the pathogenesis of brain disease²,⁷.

These results will also contribute to improving the quantification of physiological parameters from brain perfusion or functional imaging data, whether acquired by optical imaging, computed tomography, or magnetic resonance. Such quantification generally relies on the choice of mathematical functions to represent the distribution of intravascular travel times below the scale of imaging resolution. Our findings hence establish the physical grounds for defining these travel time functions and relate them to the microvascular architecture. This may help to account for vascular alterations, which are currently overlooked when interpreting human clinical imaging data in patient populations⁶², contributing to bridging the gap with knowledge acquired from animal experiments.

The analysis of transport dynamics in highly-resolved simulations of blood flow in anatomically realistic microvascular networks is complementary to in vivo experiments in that it provides access to the full flow and transport statistics, opening a window to deciphering the underlying non-local physics. With the fast progress of numerical simulations and imaging capacities, the presented framework may be further improved. For example, larger simulations domains, possibly up to the whole mice brain (see, e.g., ⁶³) may be considered. The diphasic nature of blood or more realistic reactive intravascular transport dynamics (see, e.g., ⁶⁴), simplified here as first-order kinetics and neglecting the oxygen-binding cooperativity to hemoglobin, could also be considered in the future, as well as transport and metabolic processes within the brain tissue (see, e.g., ²⁵,⁵⁴,⁶⁵).

Our theoretical framework hence opens perspectives for the development of predictive, physics-based, transport models at the scale of brain microvascular networks that account for the complexity of microvascular architectures. The resulting scaling laws are generic to a large variety of networks, from simplified ones to accurate anatomical representations, with or without capillary occlusions. This suggests that the uncovered anomalous transport mechanisms are general, even if the parameter values and pre-factors of scaling laws may be slightly dependent on the specific assumptions in our blood flow computational scheme. Thus, the variability of the network architecture, including differences between brain areas and between species³⁹,⁶⁶, differences due to long term vessel remodeling in hypoxia, aging or disease⁶¹,⁶²,⁶⁷, as well as passive or active diameter variations resulting from changes in pressure, blood flow, brain auto-regulation and/or neurovascular coupling²,⁶,²¹,⁶⁸, is unlikely to fundamentally alter the nature of the statistical laws that we have derived. Blood is ultimately transported in the brain capillary network, the structure of which is highly similar between species³⁹, so that velocity fluctuations are expected to follow similar distributions as described here. Interestingly, the predicted late-time power-law decay of travel time probabilities, with exponent -3, has even been observed with a similar exponent in the coronary microcirculation⁴⁷, despite the much stronger variability of microvascular architecture between organs⁴⁸,⁴⁹. Hence, anomalous transport induced by velocity fluctuations in the capillary network provides an alternative fundamental mechanism for this scaling, previously interpreted as arising from an underlying fractal organization of blood flow²¹,⁶⁹.
## Methods

**Brain microvascular networks.** We use a first large postmortem data set (~15,000 vessel segments in a ~1 mm³ region) from the mouse vibrissa primary sensory (vS1) cortex obtained by³,⁷⁰ (see Supplementary Note 1 and Supplementary Movie S1) and previously used for simulation studies by⁸,⁴³,⁶⁴,⁷¹,⁷². Details on the procedures used for correcting the vessel's diameters to match the in vivo distributions and for classifying vessels into arterioles, capillaries, and venules are given in the Supplementary Material from⁸. For comparison, we use a second dataset from the same cortical region, obtained from another mouse³,⁷⁰, as well as dense, space-filling networks of increasing complexity, including randomly generated networks, which reproduce well the spatial and functional statistics of the mouse capillary bed³⁹ (Supplementary Notes 3 and 8).

**Simulation of blood flow and Lagrangian transport.** Blood flow is modeled using a non-linear network approach accounting for the complex rheology of blood in the microcirculation⁸,⁷³, with prescribed boundary conditions (see Supplementary Note 1 for details on methods and validation procedures). This yields the pressure P at each vertex and the flow rate Q and haematocrit H in each vessel. We then simulate the transport of $5 \times 10^7$ passive Lagrangian particles, injected from arteriolar inlets with a probability proportional to their flow rate. We advect these particles at local average vessel velocity ($v = \langle Q \rangle / \langle Q \rangle_0 / \pi r^2$), where $\langle Q \rangle / \langle Q \rangle_0$ represents the global variation of cerebral blood flow compared with physiological baseline conditions. We distribute them in downstream vessels using flux-weighted fractions, i.e., with a probability proportional to the local flow rate (Supplementary Note 2). For each particle trajectory, we computed the corresponding trajectory length L, i.e., the total number of visited vessels, the Lagrangian transit time series, i.e., the succession of local advective transit times t from vessel to vessel in the trajectory (Fig. S5a, Supplementary Note 5), as well as the travel time τ, i.e., the sum of these local transit times. The effect of diffusion on intra-vascular transport is taken into account by introducing a maximum transit time equal to the diffusion time over the vessel length: the local transit times t is set equal to the local diffusion times ($t_D = \pi^2/D$) if $t_D < t_a$, where $t_a$ is the local vessel advection time. Note that to compute the advection time, we considered the average velocity in each vessel and did not resolve the Poiseuille flow profile in the direction transverse to the flow in each vessel. Resolving these velocity profiles within capillaries does not change much the travel time distributions and only slightly modifies the early travel times⁶⁴. Hence, at the network scale, the intra-vessel flow variability is negligible compared to the inter-vessel variability.

All numerical procedures for blood flow simulations have been implemented in a custom-built C++ code⁷⁴ and the post-processing algorithms for the Lagrangian analysis have been developed in Python (V2.7).

**First-order kinetics models of oxygen supply and amyloid-β clearance.** In order to evaluate the characteristic times associated with oxygen supply and amyloid-β clearance, and following¹⁴,⁵⁴,⁷⁵,⁷⁶, we assume that they can both be described by the following first-order kinetics models (see Supplementary Notes 9 and 10):

$$ \frac{\partial c_{O_2}}{\partial t} = -k_{O_2} c_{O_2} \quad (11) $$

and

$$ \frac{d(c_{A_β})}{dt} = \frac{d(c_{A_β} - c_{A_β}^T)}{dt} = -k_{A_β}(c_{A_β} - c_{A_β}^T) \quad (12) $$

where $k_{O_2}^{-1} = \tau_{O_2}^T$ and $k_{A_β}^{-1} = \tau_{A_β}^T$ are, respectively, the characteristic times for oxygen consumption in the brain tissue and for amyloid-β clearance.

Coupling these kinetics with the travel time distribution $p(T)$, we obtain the ratio between oxygen concentration at the outlet (venules), $c_{O_2}^V$, and at the inlet (arterioles), $c_{O_2}^A$,

$$ c_{O_2}^V/c_{O_2}^A = \int_0^\infty dT \exp(-k_{O_2} T) p_T(T) \quad (13) $$

In the same way, we obtain the ratio between the amyloid-β concentration at the outlet (venules), $c_{A_β}^V$, and its tissue concentration, $c_{A_β}^T$,

$$ c_{A_β}^V/c_{A_β}^T = 1 - \int_0^\infty dT \left\{ \begin{aligned} & \left(1 - \frac{c_{A_β}^T}{c_{A_β}^A}\right) \\ & \exp\left(-k_{A_β} T\right) p_T(T) \end{aligned} \right\} \quad (14) $$

Using the measured travel time distribution for oxygen (Green dots in Fig. 2c), we match Eq. (13) to the typical resting oxygen extraction fraction of ≈30%¹⁴, which gives $c_{O_2}^V/c_{O_2}^A \approx 0.7^{14}$, yielding $\tau_{O_2}^T = k_{O_2}^{-1} = 1.5$ s. Interestingly, this time is of the same order as the decay time measured in a single cell in vitro measurements of RBC oxygen desaturation dynamics⁷⁷, equal to 800 ms. Similarly, using typical amyloid-β concentration values (see, e.g., ⁵⁴,⁷⁸) and an arterio-venous increase of 20%⁵⁵ yields $\tau_{A_β}^T = k_{A_β}^{-1} = 97$ s. This estimated time is of the same order as can be deduced by abluminal-to-luminal permeability measurements in an in vitro blood-brain barrier model (hCMEC/D3 endothelial monolayers)⁷⁹, which yielded

a permeability $P \approx 18 \times 10^{-5}$ cm/min. Assuming an endothelial thickness of ~1 µm, this leads to a characteristic time of 33 s. By contrast, previous estimates based on a compartmental model yielded $\tau_{A_β}^T = 3000$ s⁷⁶, i.e., two orders of magnitudes above the in vitro results. Hence, despite the simplified reaction kinetics considered here, capturing the full range of travel times appears to be a key element for modeling reactive processes in the cortex.

**Reporting summary.** Further information on research design is available in the Nature Research Reporting Summary linked to this article.

## Data availability

The flow simulation data analyzed in this manuscript and the corresponding Lagrangian trajectories have been made available on Harvard Dataverse (https://doi.org/10.7910/DVN/QJDUUA).

## Code availability

The C++ flow simulation code and the custom Python (V2.7) post-processing algorithms are available within 2 months from the corresponding authors upon request subject to a nonexclusive, revocable, non-transferable, and limited right to use for research and evaluation purposes only, excluding any commercial use.

Received: 17 February 2021; Accepted: 18 November 2021; Published online: 15 December 2021

## References

1. Duvernoy, H. M., Delon, S. & Vannson, J. Cortical blood vessels of the human brain. *Brain Res. Bull.* **7**, 519–579 (1981).

2. Iadecola, C. Neurovascular regulation in the normal brain and in Alzheimer's disease. *Nat. Rev. Neurosci.* **5**, 347–360 (2004).

3. Blinder, P. et al. The cortical angiome: an interconnected vascular network with noncolumnar patterns of blood flow. *Nat. Neurosci.* **16**, 889–897 (2013).

4. Tarasoff-Conway, J. M. et al. Clearance systems in the brain–implications for Alzheimer disease. *Nat. Rev. Neurol.* **11**, 457 (2015).

5. Kisler, K., Nelson, A. R., Montagne, A. & Zlokovic, B. V. Cerebral blood flow regulation and neurovascular dysfunction in Alzheimer disease. *Nat. Rev. Neurosci.* **18**, 419–434 (2017).

6. Attwell, D. et al. Glial and neuronal control of brain blood flow. *Nature* **468**, 232–243 (2010).

7. Zlokovic, B. V. Neurovascular pathways to neurodegeneration in Alzheimer's disease and other disorders. *Nat. Rev. Neurosci.* **12**, 723–738 (2011).

8. Cruz-Hernández, J. C. et al. Neutrophil adhesion in brain capillaries reduces cortical blood flow and impairs memory function in Alzheimer's disease mouse models. *Nat. Neurosci.* **22**, 413–420 (2019).

9. Nortley, R. et al. Amyloid beta oligomers constrict human capillaries in Alzheimer's disease via signaling to pericytes. *Science* **365**, eaav9518 (2019).

10. Korte, N., Nortley, R., & Attwell, D. Cerebral blood flow decrease as an early pathological mechanism in Alzheimer's disease. *Acta Neuropathol.* **140**, 793–810 (2020).

11. Cortes-Canteli, M. & Iadecola, C. Alzheimer's disease and vascular aging. *J. Am. Coll. Cardiol.* **75**, 942–951 (2020).

12. Lorthois, S. & Cassot, F. Fractal analysis of vascular networks: insights from morphogenesis. *J. Theor. Biol.* **262**, 614–633 (2010).

13. Hirsch, S., Reichold, J., Schneider, M., Szekely, G. & Weber, B. Topology and hemodynamics of the cortical cerebrovascular system. *J. Cereb. Blood Flow. Metab.* **32**, 952–967 (2012).

14. Jespersen, S. N. & Østergaard, L. The roles of cerebral blood flow, capillary transit time heterogeneity, and oxygen tension in brain oxygenation and metabolism. *J. Cereb. Blood Flow. Metab.* **32**, 264–277 (2012).

15. Erdener, Ş. E. & Dalkara, T. Small vessels are a big problem in neurodegeneration and neuroprotection. *Front. Neurol.* **10**, 889 (2019).

16. Parpaleix, A., Houssen, Y. G. & Charpak, S. Imaging local neuronal activity by monitoring PO₂ transients in capillaries. *Nat. Med.* **19**, 241–246 (2013).

17. Sakadźić, S. et al. Large arteriolar component of oxygen delivery implies a safe margin of oxygen supply to cerebral tissue. *Nat. Commun.* **5**, 5734 (2014).

18. Moeini, M. et al. Compromised microvascular oxygen delivery increases brain tissue vulnerability with age. *Sci. Rep.* **8**, 8219 (2018).

19. Dong, S. Cerebral hypoperfusion and other shared brain pathologies in ischemic stroke and Alzheimer's disease. *Transl. Stroke Res.* **9**, 13 (2011).

20. Bracko, O., Cruz-Hernandez, J. C., Park, L., Nishimura, N. & Schaffer, C. B. Causes and consequences of baseline cerebral blood flow reductions in Alzheimer's disease. *J. Cereb. Blood Flow Metab.* **41**, 1501–1516 (2021).

21. Zierler, K. Indicator dilution methods for measuring blood flow, volume, and other properties of biological systems: a brief history and memoir. *Ann. Biomed. Eng.* **28**, 836–848 (2000).
22. Lassen, N. A., Henriksen, O. & Sejrsen, P. Indicator methods for measurement of organ and tissue blood flow in *Comprehensive Physiology* (ed. Terjung, R.) (John Wiley & Sons, Inc., 2011).

23. Lassen, N. A. & Perl, W. Tracer kinetic methods in medical physiology. pp 189 (Raven Press, 1979).

24. Lipowsky, H. H., Mckay, C. B. & Seki, J. Transit time distributions of blood flow in the microcirculation in *Microvascular Mechanics* (eds Lee J.-S. & Skalak, T. C.), 13–27, (Springer New York, 1989).

25. Angleys, H., Østergaard, L. & Jespersen, S. N. The effects of capillary transit time heterogeneity (CTH) on brain oxygenation. *J. Cereb. Blood Flow. Metab.* **35**, 806–817 (2015).

26. Li, B. et al. More homogeneous capillary flow and oxygenation in deeper cortical layers correlate with increased oxygen extraction. *eLife* **8**, e42299 (2019).

27. Gutierrez-Jimenez, E. et al. Effect of electrical forepaw stimulation on capillary transit-time heterogeneity (CTH). *J. Cereb. Blood Flow. Metab.* **36**, 2072–2086 (2016).

28. Eskildsen, S. F. et al. Increased cortical capillary transit time heterogeneity in Alzheimer's disease: a DSC MRI perfusion study. *Neurobiol. Aging* **50**, 107–118 (2017).

29. Østergaard, L. et al. The capillary dysfunction hypothesis of Alzheimer's disease. *Neurobiol. Aging* **34**, 1018–1031 (2013).

30. Mundiyanapurath, S. et al. Capillary transit time heterogeneity is associated with modified Rankin scale score at discharge in patients with bilateral high grade internal carotid artery stenosis. *PloS One* **11**, e0158148 (2016).

31. Gutierrez-Jimenez, E. et al. Disturbances in the control of capillary flow in an aged APPswe/PS1δE9 model of Alzheimer's disease. *Neurobiol. Aging* **62**, 82–94 (2018).

32. Merkle, C. W. & Srinivasan, V. J. Laminar microvascular transit time distribution in the mouse somatosensory cortex revealed by dynamic contrast optical coherence tomography. *NeuroImage* **125**, 350–362 (2016).

33. Kang, P. K., Dentz, M., Le Borgne, T. & Juanes, R. Spatial Markov model of anomalous transport through random lattice networks. *Phys. Rev. Lett.* **107**, 180602 (2011).

34. Metzler, R. & Klafter, J. The random walk's guide to anomalous diffusion: a fractional dynamics approach. *Phys. Rep.* **339**, 1–77 (2000).

35. De Anna, P. et al. Flow intermittency, dispersion, and correlated continuous time random walks in porous media. *Phys. Rev. Lett.* **110**, 184502 (2013).

36. Dentz, M., Kang, P. K., Comolli, A., Le Borgne, T. & Lester, D. R. Continuous time random walks for the evolution of lagrangian velocities. *Phys. Rev. Fluid* **1**, 074004 (2016).

37. Iturria-Medina, Y. et al. Early role of vascular dysregulation on late-onset Alzheimer's disease based on multifactorial data-driven analysis. *Nat. Commun.* **7**, 11934 (2016).

38. Van Osch, M. J. et al. Model of the human vasculature for studying the influence of contrast injection speed on cerebral perfusion MRI. *Magn. Reson. Med.* **50**, 614–622 (2003).

39. Smith, A. F. et al. Brain capillary networks across species: a few simple organizational requirements are sufficient to reproduce both structure and function. *Front. Physiol.* **10**, 233 (2019).

40. Goirand, F., Georgeot, B., Giraud, O. & Lorthois, S. Network community structure and resilience to localized damage: application to brain microcirculation. *Brain Multiphysics* **2**, 100028 (2021).

41. Kurowski, P., Ippolito, I., Hulin, J., Koplik, J. & Hinch, E. Anomalous dispersion in a dipole flow geometry. *Phys. Fluids* **6**, 108–117 (1994).

42. Alim, K., Parsa, S., Weitz, D. A. & Brenner, M. P. Local pore size correlations determine flow distributions in porous media. *Phys. Rev. Lett.* **119**, 144501 (2017).

43. Schmid, F., Tsai, P. S., Kleinfeld, D., Jenny, P. & Weber, B. Depth dependent flow and pressure characteristics in cortical microvascular networks. *PLoS Comput. Biol.* **13**, e1005392 (2017).

44. Murray, C. The physiological principle of minimum work I. The vascular system and the cost of blood volume. *Proc. Natl. Acad. Sci.* **12**, 207–214 (1926).

45. Durand, M. Architecture of optimal transport networks. *Phys. Rev. E* **73**, 016116 (2006).

46. Meigel, F. J. & Alim, K. Flow rate of transport network controls uniform metabolite supply to tissue. *J. R. Soc. Interface* **15**, 20180075 (2018).

47. Bassingthwaighte, J. B. & Beard, D. A. Fractal ¹⁵O-labeled water washout from the heart. *Circ. Res.* **77**, 1212–1221 (1995).

48. Augustin, H. G. & Koh, G. Y. Organotypic vasculature: from descriptive heterogeneity to functional pathophysiology. *Science* **357**, eaal2379 (2017).

49. Waters, S. L. et al. Theoretical models for coronary vascular biomechanics: progress and challenges. *Prog. Biophys. Mol. Biol.* **104**, 49–76 (2011).

50. Nishimura, N. et al. Targeted insult to subsurface cortical blood vessels using ultrashort laser pulses: three models of stroke. *Nat. Methods* **3**, 99–108 (2006).

51. Shih, A. Y. et al. The smallest stroke: occlusion of one penetrating vessel leads to infarction and a cognitive deficit. *Nat. Neurosci.* **16**, 55–63 (2012).

52. Dahlgren, K. N. et al. Oligomeric and fibrillar species of amyloid beta peptides differentially affect neuronal viability. *J. Biol. Chem.* **277**, 32046–32053 (2002).

53. Blennow, K., Hampel, H., Weiner, M. & Zetterberg, H. Cerebrospinal fluid and plasma biomarkers in Alzheimer's disease. *Nat. Rev. Neurol.* **6**, 131–144 (2010).

54. Potter, R. et al. Increased in vivo amyloid 42 production, exchange, and loss in presenilin mutation carriers. *Sci. Transl. Med.* **5**, 189ra77–189ra77 (2013).

55. Xiang, Y. et al. Physiological amyloid beta clearance in the periphery and its therapeutic potential for Alzheimer's disease. *Acta Neuropathol.* **130**, 487–499 (2015).

56. King, R. B., Raymond, G. M. & Bassingthwaighte, J. B. Modeling blood flow heterogeneity. *Ann. Biomed. Eng.* **24**, 352–372 (1996).

57. Park, C. S. & Payne, S. J. A generalized mathematical framework for estimating the residue function for arbitrary vascular networks. *Interface Focus* **3**, 20120078 (2013).

58. Li, B. et al. Two-photon microscopic imaging of capillary red blood cell flux in mouse brain reveals vulnerability of cerebral white matter to hypoperfusion. *J. Cereb. Blood Flow. Metab.* **40**, 501–512 (2020).

59. Hartung, G. et al. Voxelized simulation of cerebral oxygen perfusion elucidates hypoxia in aged mouse cortex. *PLoS Comput. Biol.* **17**, e1008584 (2021).

60. Gorelick, P. B. et al. Vascular contributions to cognitive impairment and dementia: a statement for healthcare professionals from the american heart association/american stroke association. *Stroke* **42**, 2672–2713 (2011).

61. Secomb, T. W., Alberding, J. P., Hsu, R., Dewhirst, M. W. & Pries, A. R. Angiogenesis: an adaptive dynamic biological patterning problem. *PLoS Comput. Biol.* **9**, e1002983 (2013).

62. D'Esposito, M., Deouell, L. Y. & Gazzaley, A. Alterations in the BOLD fMRI signal with ageing and disease: a challenge for neuroimaging. *Nat. Rev. Neurosci.* **4**, 863–872 (2003).

63. Kirst, C. et al. Mapping the fine-scale organization and plasticity of the brain vasculature. *Cell* **180**, 780–795 (2020).

64. Berg, M., Davit, Y., Quintard, M. & Lorthois, S. Modelling solute transport in the brain microcirculation: is it really well mixed inside the blood vessels? *J. Fluid Mech.* **884**, 1–43 (2020).

65. Holter, K. E. et al. Interstitial solute transport in 3D reconstructed neuropil occurs by diffusion rather than bulk flow. *Proc. Natl. Acad. Sci.* **114**, 9894–9899 (2017).

66. Hartmann, D. A., Hyacinth, H. I., Liao, F.-F. & Shih, A. Y. Does pathology of small venules contribute to cerebral microinfarcts and dementia? *J. Neurochem.* **144**, 517–526 (2018).

67. Steinman, J., Sun, H.-S. & Feng, Z.-P. Microvascular alterations in Alzheimer's disease. *Front. Cell. Neurosci.* **14**, 618986 (2021).

68. Lorthois, S., Cassot, F. & Lauwers, F. Simulation study of brain blood flow regulation by intra-cortical arterioles in an anatomically accurate large human vascular network. Part II: Flow variations induced by global or localized modifications of arteriolar diameters. *NeuroImage* **54**, 2840–2853 (2011).

69. Beard, D. A. & Bassingthwaighte, J. B. The fractal nature of myocardial blood flow emerges from a whole-organ model of arterial network. *J. Vasc. Res.* **37**, 282–296 (2000).

70. Tsai, P. S. et al. Correlations of neuronal and microvascular densities in murine cortex revealed by direct counting and colocalization of nuclei and vessels. *J. Neurosci.* **29**, 14553–14570 (2009).

71. Gould, I. G., Tsai, P., Kleinfeld, D. & Linninger, A. The capillary bed offers the largest hemodynamic resistance to the cortical blood supply. *J. Cereb. Blood Flow Metab.* **37**, 52–68 (2016).

72. Sweeney, P. W., Walker-Samuel, S. & Shipley, R. J. Insights into cerebral haemodynamics and oxygen utilising in vivo mural cell imaging and mathematical modelling. *Sci. Rep.* **8**, 1373 (2018).

73. Pries, A. R., Secomb, T. W., Gaehtgens, P. & Gross, J. Blood flow in microvascular networks. Experiments and simulation. *Circ. Res.* **67**, 826–834 (1990).

74. Peyrounette, M., Davit, Y., Quintard, M. & Lorthois, S. Multiscale modelling of blood flow in cerebral microcirculation: details at capillary scale control accuracy at the level of the cortex. *PLoS One* **13**, e0189474 (2018).

75. Gagnon, L. et al. Modeling of cerebral oxygen transport based on in vivo microscopic imaging of microvascular network structure, blood flow, and oxygenation. *Front. Comput. Neurosci.* **10**, 82 (2016).

76. Shibata, M. et al. Clearance of Alzheimer's amyloid beta-40 peptide from brain by LDL receptor related protein 1 at the blood-brain barrier. *J. Clin. Invest.* **106**, 1489–1499 (2005).

77. Di Caprio, G., Stokes, C., Higgins, J. M. & Schonbrun, E. Single cell measurement of red blood cell oxygen affinity. *Proc. Natl Acad. Sci.* **112**, 9984–9989 (2015).

78. Roberts, K. F. et al. Amyloid efflux from the central nervous system into the plasma: Brain efflux of amyloid. *Ann. Neurol.* **76**, 837–844 (2014).

79. Swaminathan, S. K. et al. Insulin differentially affects the distribution kinetics of amyloid beta 40 and 42 in plasma and brain. *J. Cereb. Blood Flow. Metab.* **38**, 904–918 (2018).
## Acknowledgements

Research reported in this publication was supported by the European Research Council under ERC grant agreements 615102 (BrainMicroFlow) and 648377 (ReactiveFronts) and by the NIH (awards R21CA214299 and 1RF1NS110054). We gratefully acknowledge P. Blinder, P. Tsai, and D. Kleinfeld for sharing anatomical networks, C.B. Schaffer and N. Nishimura for inspiring discussions, N. Nishimura for critically reading a previous version of this manuscript, and M. Berg who developed the network flow solver. The funders had no role in the study design, data collection, and analysis, decision to publish, or preparation of the manuscript.

## Author contributions

All authors contributed equally to the conceptual aspects of this work and to manuscript writing. F.G. developed the numerical methods and analyzed the data. T.L.B. and S.L. designed the research, jointly supervised the work, and contributed to data analysis and interpretation.

## Competing interests

The authors declare no competing interests.

## Additional information

**Supplementary information** The online version contains supplementary material available at https://doi.org/10.1038/s41467-021-27534-8.

**Correspondence** and requests for materials should be addressed to Tanguy Le Borgne or Sylvie Lorthois.

**Peer review information** *Nature Communications* thanks Stephen Payne and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available.

**Reprints and permission information** is available at http://www.nature.com/reprints

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_12.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2021