We propose a model of lymphatic drainage from the brain via the cerebral BM based on Darcy's law for flow in porous media. By exploiting the discrepancy of scales between the width of the BM and its length we approximate the flow by a one-dimensional lubrication model (see for example Ockendon and Ockendon, 1995). Furthermore, we allow the permeability of the BM to depend on the pressure gradient of the ISF. Such a pressure gradient dependent permeability could, for instance, model proteins in the BM protein network that have structural properties that allow them to bend easily only in one direction, thereby acting as a valve, and thus provide a higher resistance to flow in the direction of the blood flow. A similar idea has been investigated in Sharp et al. (2015), where proteins were explicitly modeled as cilia. However, here we model the active component of the BM generically in this model and thus allow for alternative explanations of the phenomenon.

The BM consists of a complex matrix of proteins (extracellular matrix, ECM) that can be interpreted as a porous medium through which soluble metabolites drain via the ISF. We model flow through this medium, in the standard fashion, by Darcy's law in a 2D cylindrical coordinate system. Darcy's law states that the flux q is related to the fluid pressure gradient ∇p by

where k is the intrinsic permeability of the ECM and μ is the viscosity of ISF. Since ISF flows up the artery, counter to the direction of arterial pulsations, a valve-like mechanism in the ECM is required if arterial pulsations are indeed the motive force for this flow. This valve-like property of the ECM can be modeled by a Darcy flow Equation (1), in which the intrinsic permeability is a function of the axial pressure gradient ∂p/∂z along the vessel. On assuming an artery aligned along the z-axis with net blood flow and arterial pulsations propagating the positive z-direction we model valve-like properties of the ECM in the BM by a pressure gradient dependent permeability

In order to obtain a net flux up the vessel, counter to the blood flow, we require this function to be increasing in p_z, representing the notion that it is harder to push a flow in the positive z-direction (p_z < 0 and K(p_z) small) than in the negative z-direction (p_z ≥ 0 and K(p_z) large).

The governing equations in cylindrical coordinates are thus

where e_z and e_r describe the unit vectors in the z and r direction. Incompressibility of the fluid is represented by Equation (4).

To determine the boundary conditions of the system consider the geometry in Figure 1. Under the assumption that the vessel deformations are radially symmetric the positions of the inner and outer sheaths bounding the BM can be described by the two expressions r = R_i(z, t) and r = R_o(z, t), respectively. Assuming that the flows through the surfaces of the BM are negligible the following kinematic boundary condition applies

and analogously for R_o(z, t).

Because the BM thickness is much smaller than the vessel length it is possible, on using the lubrication approximation, to derive a 1D model (in z) for the BM thickness h_bm(z, t); this takes the form

where γ is the fluid volume fraction in the BM (see also Supplementary Material 1). Equation (6) describes the evolution of the BM thickness and can be solved using a finite difference approximation. For a detailed derivation see Supplementary Material 1. Another method to arrive at the same equation is derived by considering an infinitesimally small volume of BM dV = 2π · R_i(z, t) · h_bm(z, t) · dz. Then the difference between the flux into and out of the element dV is equal to the rate of change of the fluid volume ∂(γdV)/∂t, which can be written as

which is equal to the model Equation (6) in the limit dz → 0. The pressure p(z, t) and the displacement of the BM wall R_i(z, t) are, to a very good approximation, determined solely by arterial pressure (i.e., almost independent of the BM flow), see Supplementary Material 2 for more details. R_i(z, t) and p_z are obtained from a numerical artery model described in the next section. To model the protein valve that is being proposed here, K(p_z) is defined as a step function

such that K₁ > K₀ > 0 so that the BM is more permeable to ISF flow in the negative z-direction than in the positive one.

In order to obtain the input functions R_i(z, t) and p_z for the analytical model developed in the previous section, a section of human MCA was modeled as a 1D axisymmetric tube in a cylindrical coordinate system with cross-sectional lumen area A(z, t) and flux Q(z, t) using the Python package VaMpy (Diem and Bressloff, 2016). The geometry parameters are listed in Table 1. Blood flow through arteries and the propagation of arterial pulsations are commonly described by a one-dimensional averaged model (Olufsen et al., 2000) (for more details see Diem, 2016). It is represented by the conservation law

where

The equation is solved for cross-sectional area A(z, t) of the vessel lumen and flux through the vessel Q(z, t), where R(z, t) = (A(z, t)/π)^1/2 is the lumen radius, δ_b is the boundary layer thickness and Re is the Reynolds number. Elasticity of the artery wall is modeled via the quantity f(r₀) = 4Eh/(3r₀), where E is the Young's modulus and h is the wall thickness and r₀(z) is the radius at rest. The relationship Eh/r₀ is based on compliance estimates

with k₁, k₂, k₃, see Table 2 (Olufsen et al., 2000). Blood pressure p(z, t) is related to A(z, t) and f(r₀) via the state equation

where p₀ is the diastolic pressure and A₀(z) is the lumen radius at rest. Blood pressure provides the basis for determining ISF pressure inside the BM via the calculations in Supplementary Material 2.

To derive Equation (9) it was necessary to make assumptions about the velocity profile of blood flow through an artery. Blood flow is considered pulsatile laminar and vessels are considered slightly tapered, therefore the velocity profile is assumed to be mostly flat with a thin boundary layer δb=(νT/(2π))0.5, such that δ_b ≪ R(z, t) (Olufsen et al., 2000).

This system of equations can be solved numerically using Richtmyer's two-step version of the Lax-Wendroff method, as done by others (LeVeque, 1992; Olufsen et al., 2000; Kolachalama et al., 2007). We use the Python module from Diem and Bressloff (2016) to implement a model of the MCA and its first bifurcation. All simulation parameters are listed in Table 2. The MCA geometry data was taken from Cousins and Gremaud (2014), estimating lengths by 50 · R_u. The radii of the daughter vessels were estimated using a scaling factor of 0.91 (Cousins and Gremaud, 2014). The length of the right daughter vessel was shortened slightly to avoid reflective waves canceling one another out (Kolachalama et al., 2007). A three-element Windkessel (3WK) model is used at the outlet whose parameters were obtained from Olufsen et al. (2002).

At the inlet flux is prescribed using patient-specific flow velocity measurements collected from the MCA of a healthy adult male using Doppler Sonography ultrasound. The velocity data was converted to volumetric flux values within the range reported in Olufsen et al. (2002). To obtain a smooth inlet function the final ten peaks of the time series were averaged in the Fourier space. The resulting inlet boundary condition is shown in Figure 2.

In the previous sections the models governing blood flow through the MCA and ISF flow through the artery wall have been introduced. Here, the results of the previous sections are coupled and used to calculate ISF flow through the BM with and without a valve mechanism to show its necessity to achieve net reverse drainage. At the same time it is shown that blood pressure driven perivascular drainage is far too slow to provide any meaningful flow, thereby disproving the popular hypothesis that arterial pulsations provide the major driving force for perivascular drainage. MCA simulations were run for eight cardiac cycles to allow the system to settle. An average human heartbeat lasts 0.85 s and the total simulation time was 6.8 s (see Table 2). In the results presented here only the final heartbeat is shown. Figure 3 shows blood pressure (a), wall displacement (b), and ISF pressure (c) in the MCA. Because the radius of the MCA is very small its wall stiffness is high. Pressure gradients in time are therefore steep. This data provides the basis for estimating perivascular drainage through the cerebral vasculature.

Figure 3B shows the displacement of the artery wall due to the pressure pulse. To calculate the width of the BM using the model Equation (6), the input function R_i(z, t) is the sum of the radius at rest and the wall displacement, as depicted in Figure 1. Figure 4A shows the width of the BM h_bm(z, t) as calculated from Equation (6) over the length of the artery at a number of time points, while Figure 4B shows ∂h_bm(z, t)/∂z for the same time points. The BM width only varies minimally around its initial width and thus flow rates of ISF are very small.

The results show that a valve mechanism is required to drive perivascular drainage in the reverse direction of the blood flow with an average flux of −1.12e−3 μm3 s⁻¹ for a ratio of K₀/K₁ = 0.01, while it is 6.33e−2 μm3 s⁻¹ for K₀/K₁ = 1.0. While the results confirm the necessity of a strong valve mechanism under pulse driven flow they also indicate that arterial pulsations have too long a wavelength to drive significant perivascular drainage. Drainage velocity is fastest without a valve mechanism (1.81e−2 μm s⁻¹, and then it is in the wrong direction), but even in that case it is almost three orders of magnitude slower than the expected value of 8.33 μm s⁻¹ (Carare et al., 2008). The valve mechanism only guarantees net reverse drainage for a permeability ratio K₀/K₁ < 2.74e−2, where K₁ is chosen such that k = 1e−10 cm² is consistent with values for other tissues in the body and μ = 1.5e−3 Pa s (Heppell et al., 2013). Figure 5 illustrates flux through the BM for a number of ratios K₀/K₁ and for different positions of the BM r = r₀ + ηh. Reverse flow is achieved for K₀/K₁ < 2.74e−2. Drainage flux values are fastest when the BM is placed close to the arterial lumen (η = 0) due to the elasticity of the artery wall.

The minute flux values are a result of the very small variations in BM width (see Figure 4) over the length of the MCA, in fact they are less than 1 Å. This value is far too small to be realistic and demonstrates that arterial pulsations cannot be strong enough to drive perivascular drainage through the BM. To illustrate it in another way, the total ISF volume in a human brain is 280 ml (Weller et al., 2015). At a drainage rate of −1.12e−3 μm3 s⁻¹ per cerebral artery the total turnover time for ISF would be 4.83e11 days, much more than a human life span.

Brains were cut into 20 μm coronal sections using a freezing microtome (Leica) and the injection site located by examining on a Zeiss Axioskop 2 fitted with a rhodamine filter. Sections 400 and 800 μm both posterior and anterior to the injection site were selected for immunohistochemistry. Sections were blocked in 15% goat serum (Sigma 9023) for 1 h at room temperature (RT). They were then incubated in primary antibodies diluted in 0.01 M phosphate buffered saline 0.1% triton x100 (PBSt). Primary antibodies comprised rabbit anti-collagen IV (AbCam, ab6586) 1/400 and mouse anti-smooth muscle actin (SMA) FITC conjugated (Sigma F3777) 1/200. Sections were incubated overnight at 4 °C in a moist chamber. Following washing 3 × 10 min 0.01 M PBS sections were incubated in goat anti-rabbit Alexa Fluor 633 (ThermoFisher Scientific A-21070) 1/200 PBSt.

Images of the whole hippocampus were captured on a Leica SP8 confocal microscope fitted with LASX software in tile scanning mode with sequential imaging. Once captured, images were exported as TIFF files and uploaded into Adobe Photoshop CS6. An area of hippocampus was drawn extending transversely from CA3 pyramidal cells to the apex of the blades of the granule cell layer and from the ventricle to stratum oriens in the superior/inferior axis. Aβ containing vessels were defined according to size and presence of SMA: capillaries had a diameter of ≤ 10 μm, veins had no SMA and a diameter of > 10 μm. The number of vessels was expressed per mm2 of hippocampus. Statistical analysis was performed by an independent t-test with a confidence interval of 95%.

Within 5 min of injections into the hippocampus, fluorescent Aβ was observed diffusely in the hippocampus and co-localizing with collagen IV with collagen IV in the walls of capillaries and arteries, with very few veins involved in both Ate-treated and control mice (Figure 6). Application of an independent t-test to the data revealed no significant differences between the saline and Ate treated mice either anterior or posterior to the injection site (Figure 7).

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest. The reviewer AA and handling Editor declared their shared affiliation.