Impact of stochastic fluctuations in the cell free layer on nitric oxide bioavailability

A plasma stratum (cell free layer or CFL) generated by flowing blood interposed between the red blood cell (RBC) core and the endothelium affects generation, consumption, and transport of nitric oxide (NO) in the microcirculation. CFL width is a principal factor modulating NO diffusion and vessel wall shears stress development, thus significantly affecting NO bioavailability. Since the CFL is bounded by the surface formed by the chaotically moving RBCs and the stationary but spatially non-uniform endothelial surface, its width fluctuates randomly in time and space. We analyze how these stochastic fluctuations affect NO transport in the CFL and NO bioavailability. We show that effects due to random boundaries do not average to zero and lead to an increase of NO bioavailability. Since endothelial production of NO is significantly enhanced by temporal variability of wall shear stress, we posit that stochastic shear stress stimulation of the endothelium yields the baseline continual production of NO by the endothelium. The proposed stochastic formulation captures the natural continuous and microscopic variability, whose amplitude is measurable and is of the scale of cellular dimensions. It provides a realistic model of NO generation and regulation.


INTRODUCTION
Nitric oxide (NO) plays a critical role in the local control of smooth muscle tone and the regulation of blood flow at the microvascular level. Its distribution in the microcirculation is determined by the balance between NO production and consumption in the blood and tissue compartments. The local concentration of NO in blood results from the competition between NO diffusing from the endothelium and NO scavenging by hemoglobin in red blood cells (RBCs) or at times dissolved in plasma. Mathematical modeling was used (Vaughn et al., 1998a,b;Condorelli and George, 2002;Kavdia and Popel, 2003;Lamkin-Kennard et al., 2004;Chen et al., 2006;Ong et al., 2011b;Sriram et al., 2011) to determine the NO distribution in blood vessels simulated by cylindrical and parallelplate compartments, as a function of local transport parameters such as NO production rate, scavenging reaction rate and diffusion coefficients in blood and tissue. All these analyses assume deterministic boundary conditions.
Mechanotransduction generates a significant portion of the NO involved in the regulation of blood flow (Condorelli and George, 2002;Chen et al., 2006). This mechanical effector links the biochemistry of NO production by the endothelium with shear stress induced by blood flow on the vascular wall (WSS). This coupling is tight because NO bioavailability in the vascular wall determines vessel diameter, the anatomical component of vascular flow resistance. WSS is generated at the vascular wall by tangential stresses caused by flowing blood whose RBC concentration (hematocrit or Hct) diminishes from maximal at the blood flow core to zero in the cell free layer (CFL) adjacent to the vessel wall.
Blood in microcirculation can be modeled as a two-layer system consisting of two immiscible fluids. Flow velocity profiles are parabolic within the CFL and plug-like in the RBC core region (Sharan et al., 1997;Martini et al., 2005;Sriram et al., 2011) where the fluid exhibits non-Newtonian behavior. Such models provide a coarse-grained description of more detailed simulations of variations of the WSS induced by a passing column of RBCs (e.g., Dupin et al., 2007). Changes of NO production due to changes in WSS caused by the variation of blood flow and plasma viscosity (Tsai et al., 1998) appear to dominate other factors, such as changes in hemoglobin concentration or Hct (Kavdia and Popel, 2003;Lamkin-Kennard et al., 2004;Chen et al., 2006;Namgung et al., 2011;Sriram et al., 2011). Experimental evidence (e.g., Kanai et al., 1995) suggests that NO production rate by the endothelium varies linearly with WSS (see also Andrews et al., 2010;Ong et al., 2011b).
However, there is a fundamental difference in the rate of NO production between steady and time-varying WSS (Frangos et al., 1996;Ong et al., 2011b), the latter being significantly greater. Temporal variability arises from action of the heart and the respiratory cycle which although attenuated is present in the microcirculation, or vasomotion which is rooted in the microcirculation.
Endothelial response to shear stress is mostly studied in parallel plate flow chambers, where endothelial cells grow to confluent layers. However, it is generally recognized that this experimental setup differs from in vivo conditions, since it relies on perfusion with cell culture media fluid instead of blood. (We are aware of a study by Yalcin et al., 2008 that did use RBC suspension, but it did not compare results with incubation medium). Blood perfusion introduces other, hitherto not considered, sources of temporal variability due to the microscopic spatial fluctuations of the CFL width. Spatial variability is generated by blood flowing over the spatially nonuniform endothelial surface (Barbee et al., 1994;Sato et al., 2000). Spatial and temporal variability is generated at the blood boundary of the CFL due to the changes in motion, position and shape of RBCs in the outer layer of the flowing blood column. As a consequence of these phenomena fluid in the CFL flows between a spatially variable endothelial boundary and spatially and temporally highly irregular RBC surface (Kim et al., 2006(Kim et al., , 2007, both of which randomly affect WSS. These effects are non-existent in the absence of RBCs. Furthermore, as a consequence of their stochastic nature they potentially include all forms of variability. The frequency of stochastic fluctuations of the WSS due to spatio-temporal variability of the CFL may significantly increase the bioavailability of NO. Evidence for this hypothesis is provided by an observation that exposing the tissue damaged by ischemia reperfusion to diagnostic ultrasound improves microvascular functionality, an effect that is significantly reduced by administration of L-NAME (Hightower and Intaglietta, 2008).
These unpredictable fluctuations can be analyzed by treating the CFL width at any given location as a random field that determines the distance between the NO source and its sink. The stochasticity of the CFL width also affects NO concentration due to its repercussion on the local, microscopic variability of the flow field, and therefore on shear rates and shear stresses. While the random spatio-temporal fluctuations of CFL width and its bounding surfaces clearly affect NO bioavailability and production rate, most studies (including those mentioned above) treat the interface between the RBC column and plasma as a smooth deterministic surface. In the present study, we adopt a more realistic approach by treating the surface between flowing RBCs and the CFL as a random field whose statistics are obtained from experimental studies (Ong et al., 2011a,b). Within this conceptual framework, we formulate a model that determines the distribution of NO concentration in the region of the interface between blood and tissue at the blood vessel wall (RBC-rich core, cell free layer, and tissue layer) (Vaughn et al., 1998a).

Model Formulation
Our analysis deals with NO transport in the geometrical configuration associated with the standard Krogh tissue cylinder model. We consider an arteriolar cross-section that consists of the RBC-rich core (E 1 : 0 ≤ r ≤ r 1 ), the CFL (E 2 : r 1 < r ≤ r 2 ), the endothelial-cell region (E 3 : r 2 < r ≤ r 3 ), and the smoothmuscle region (E 4 : r 3 < r ≤ r 4 ). Stochastic fluctuations of the interface formed by flowing RBCs, r 1 (θ, t), are modeled by treating it as a random function of both angular coordinate θ and time t, i.e., r 1 = r 1 (θ, t; ω) with ω ∈ indicating a realization ("coordinate") in the probability space . This renders the CFL width w = r 2 − r 1 random, i.e., w = w(θ, t; ω). Our goal is to capture the effects of stochastic fluctuations of the RBC-CFL interface r 1 (θ, t; ω) on distribution of NO concentration, C NO , in the Krogh tissue cylinder D = {(r, θ) : 0 ≤ r ≤ r 4 , 0 ≤ θ ≤ 2π}.
In each region of the computational domain, E i (i = 1, . . . , 4), the concentration C NO satisfies a reaction-diffusion equation where D i and k i are the diffusion coefficient and degradation (reaction) rate in the i-th region, respectively. These four equations are coupled by the continuity conditions at the interfaces r i (i = 1 . . . , 3), Here the superscripts − and + indicate the left and right limits of the corresponding quantities at the i-th interface, F n = F i · n i is the normal component of Fick's flux F i = −D i ∇C NO at the i-th interface whose outward unit normal is n i , andq i denotes the NO production rates at the interface r = r i . Since only endothelium cells are involved in NO production,q 1 ≡ 0. We assume that no nitric oxide leaves the outer boundary of the smooth-muscle region, r = r 4 , so that n · ∇C = 0, r = r 4 .
The coupling of the reaction-diffusion Equations (1) at the interfaces r = r i (i = 1, 2, 3) propagates uncertainty (randomness) in the topology of the RBC-CFL interface r 1 (θ, t; ω) through the modeling process, leading to randomly varying NO concentration C NO (r, θ, t; ω) throughout the Krogh tissue cylinder. The problem formulation given by Equations (1)-(3) implicitly assumes that blood flow is laminar, fully-developed, and incompressible, vessel walls are impermeable to blood flow, NO concentration at the vessel inlet equals that at the vessel outlet, and all the reaction rates are spatially uniform. In the deterministic setting with a uniform CFL width, these assumptions imply that the radial (u r ) and angular (u θ ) components of flow velocity u = (u r , u θ , V) ⊤ are u r = u θ ≡ 0, its longitudinal component is V = V(r), and NO concentration C NO = C NO (r). This results in advective flux of NO in the blood vessel, u · ∇C NO , that is identically zero. In the stochastic setting with randomly fluctuating CFL width, advective flux u · ∇C NO is zero in the mean (Tartakovsky and Xiu, 2006;Park et al., 2012).

Model Parameterization
While the reaction rates k i in the endothelium (i = 3) and tissue (i = 4) can be considered constant, the reaction rate in the RBC-rich core (k 1 ) is related to hemoglobin levels. The latter depends on hematocrit H(r) and radial component of blood flow velocity V(r, θ, t). Let k s denote a reference rate of NO scavenging by RBCs at a reference level of hematocrit H s . Then, the NO scavenging rate k 1 corresponding to a given hematocrit level H c is given by Ong et al. (2011b) and Chen et al. (2006) The hematocrit ratio H c /H s is determined by mass conservation, In the general stochastic framework we advocate here, blood is a two-phase fluid that exhibits non-Newtonian behavior in the RBC-rich core and Newtonian one in the CFL, with the random surface r 1 (θ, t; ω) separating the two regions. This implies that flow velocity V(r, θ, t; ω) is random as well, being given by a solution of corresponding flow equations in random domains (Park et al., 2012). To focus on NO transport, we simplify the flow calculations by adopting two alternative approximations. The first is based on a lubrication approximation in which random geometry parameterizes an otherwise deterministic velocity profile (Tartakovsky and Xiu, 2006). This approach yields a random velocity profile V(r, θ, t; ω), where V max = Jr 2 2 /(4µ p ) is the (maximum) centerline velocity, J is the externally imposed pressure gradient, and µ p and µ c are the viscosities of the plasma and RBC-rich core. In this formulation, the only source of the non-Newtonian behavior of the RBC-rich core is the dependence of the core viscosity µ c on the (random) CFL width. Following Martini et al. (2005) and many others, we assume a linear relationship µ c = 0.1678H c − 4.348 between a hematocrit level H c and the viscosity of the RBC core µ c . Specifying a (random) radial distribution of hematocrit, H = H(r, t; ω), as a step function enables one to compute the randomly fluctuating NO scavenging rate k 1 (t; ω) by combining Equations (4)-(7). First, the system of Equations (5)-(7) was solved using Matlab function "solve" to compute H c for µ p = 1.2 cP and two values of H s . Then k 1 (t; ω) was obtained from Equation (4). The second alternative for obtaining k 1 (t; ω) treats blood as a single-phase fluid with a parabolic velocity profile This formulation replaces the CFL and the random RBC-CFL interface r 1 (θ, t; ω) with a radial distribution of hematocrit, Substituting Equations (8) and (9) into Equations (4) and (5) yields an alternative expression for the NO scavenging rate k 1 (t; ω). This approach was used by Ong et al. (2011a) in the deterministic context that treated r 1 (θ, t) as constant. Finally, we allow the NO production rates by the endothelium, i.e.,q 2 andq 3 in Equation (2), to vary with the wall shear stress τ w exerted on the endothelium walls by blood flow. Following Ong et al. (2011a,b), Vaughn et al. (1998a) and others, we assume a linear relationq where τ ⋆ w is the reference wall shear stress,q ⋆ is the control NO production rate, and V e is the mean velocity at the outer edge of the RBC core. These production rates fluctuate randomly, i.e., q 2 (t; ω) andq 3 (t; ω), due to their dependence on the random flow velocity V and the CFL width w(θ, t; ω) = r 2 − r 1 (θ, t; ω).
In the numerical results reported below we assume the diffusion coefficients D i in Equation (1) to be the same and equal to D. Its value and the values of the remaining parameters used in our are model are reported in Table 1.
We represent spatio-temporal variations of the RBC-CFL interface, as the product of mutually uncorrelated temporal and angular fluctuations r t (t; ω) and r θ (θ; ω), respectively. A Reynolds decomposition is used to represent each of these fields, r = r + r ′ , as the sum of its ensemble mean r and zero-mean fluctuations r ′ . Setting r θ = 1 yields the mean and variance of the RBC-CFL interface: r 1 = r t and σ 2 r = r 2 t σ 2 θ + σ 2 t (1 + σ 2 θ ). The coefficient of variation of the CFL width, CV w = σ w /w, is given by where w = r 2 − r 1 is the mean CFL width and σ w is its standard deviation.
Since the random field r ′ θ (θ, ω) is periodic, a truncated Fourier-type expansion provides its natural representation. Here the eigenvalues ν n (ω) are complex zero-mean random variables, whose real and imaginary parts are mutually independent for all n. Each has zero mean and variance σ 2 n = C n /4, where are coefficients of the Fourier cosine expansion of a 2π-periodic covariance function C p θ of the random field r ′ θ (θ, ω). It is constructed as follows. First, we note that statistics of r ′ θ (θ; ω) are rotationally invariant on the circle, such that a covariance function C θ is Then C p θ is constructed by extending the covariance function C θ of the random field r ′ θ (θ, ω) to a 2π-periodic periodic domain. We employ a Gaussian covariance function C θ ( θ ) = exp(− 2 θ /l 2 θ ) with the correlation length l θ . The decay of the Fourier cosine coefficients C n determines the number of terms N θ in the expansion in Equation (13) that is required to achieve a given truncation error. As the correlation length l θ decreases, N θ increases.
We represent the random field r ′ t (t; ω) via a truncated Karhunen-Loéve expansion, where Y m (ω) (m ≥ 1) are independent random variables, and λ m and f m (t) are, respectively, the eigenvalues and eigenfunctions of Fredholm equations, For an exponential correlation function ρ t (t, t ′ ) = exp(−|t − t ′ |/l t ) with the correlation length l t > 0, the eigenvalue problems in Equation (17) admit an analytical solution (Lin et al., 2010), where γ m are solutions of (l 2 t γ 2 − 1) sin(γ T) = 2l t γ cos(γ T) and m ≥ 1. The truncation error of the Karhunen-Loéve expansion in Equation (16) depends on the correlation length l t . The smaller the correlation length, the more terms N t are necessary to represent the random field r ′ (t, ω) with a given degree of accuracy.
Within the statistical framework adopted here, the random RBC-CFL interface is characterized by four parameters: variances σ 2 θ and σ 2 t , and correlation lengths l θ and l t . Experimental data, such as those reported by Kim et al. (2007), can be used to estimate these statistics. Table 2 contains the values of these parameters used in our simulations.

Transformed Stochastic Equations
The mapping defined by Equation (19) renders the transformation Jacobian and other related metrics stochastic, i.e., dependent on a set of K = 2N θ + N t independent random variables {Y i (ω)} K i=1 . The first N t variables Y 1 , . . . , Y N t coincide with those introduced in Equation (16) and the remaining 2N θ variables represent their counterparts in Equation (13), such that Y N t +1 = ν −N θ , . . . , Y K = ν N θ . Consequently, the deterministic reactiondiffusion Equations (1) are transformed into stochastic equations of the form (Appendix) where the random coefficients A 11 , A 12 = A 21 , and A 22 are given by Equation (A2) in the Appendix. These stochastic differential equations on the deterministic domain B can be solved with a variety of well-established techniques, including perturbation-based moment equations (Tartakovsky and Winter, 2001), stochastic finite elements (Ghanem and Spanos, 1991), and stochastic collocation on sparse grids (Lin et al., 2010, and the references therein). In the subsequent numerical simulations we employ the latter approach (Appendix).

RESULTS
Solutions of the stochastic system of transport Equations (1)-(3) are given in terms of their statistical moments. Ensemble means, e.g., mean NO concentrationC NO , serve as unbiased predictors of the system behavior; variances, e.g., NO concentration variance σ 2 C , act as a measure of predictive uncertainty.

Data-driven Model Parameterization
The CFL width measurements (Ong et al., 2011a) are used to construct a probabilistic model for the random input parameter w(θ, t; ω) = r 2 − r 1 in which the random RBC-CFL interface r 1 (θ, t; ω) is given by Equation (11). These data, which represent temporal fluctuations of w at a single spatial location (say, θ = 0), give rise to the histogram and auto-correlation reported in Figures 1 and 2, respectively. This histogram (and the obvious fact that the CFL width is both non-negative and smaller than the vessel radius r 2 ) indicates that the random field w(θ, t; ω) is non-Gaussian. We fit the histogram in Figure 1 with a beta distribution f w (W) = B −1 r  complete gamma function, 0 ≤ W ≤ r 2 , and α > 0 and β > 0 are shape parameters. Setting α = 4.358 and β = 32.9 provides the best data fit, resulting in the mean CFL widthw = 2.73 µm. The auto-correlation data in Figure 2 were fitted with an exponential correlation function ρ(t, t ′ ) = exp(−|t − t ′ |/l t ), yielding the correlation length l t = 0.007 s.
Experimental limitations preclude data acquisition at multiple azimuths θ, which requires us to postulate a probabilistic model for r θ (θ; ω). In analogy with its temporal counterpart r t (t; ω), we chose r θ (θ; ω) to have the beta distribution with unit mean and variance σ 2 θ and the exponential correlation function with correlation length l θ . In the formulation provided by Equation (12), the amplitude of spatio-temporal (in the angular coordinate θ and time t) fluctuations of both the RBC-CFL interface r 1 (θ, t; ω) and CFL width w(θ, t; ω) = r 2 − r 1 increases with the variances σ 2 θ and σ 2 t , while the smoothness of these fluctuations increases with the correlation lengths l θ and l t . This behavior, which reflects chaotic motion of RBCs in the blood core, is demonstrated by two representative realizations of the random CFL width shown in Figure 3.

Random Fluctuations of Wall Shear Stress
The CFL width w in Equation (10) is inversely proportional to the wall shear stress (WSS) τ w . Hence the random spatio-temporal fluctuations in w induce corresponding fluctuations in τ w . In the computations of the WSS we use the values of the edge velocity V e = 0.54 mm/s and the corresponding pressure gradient J = 2.15 × 10 4 computed by fitting the smooth-wall model to the experimentally observed peak NO concentration of 11.2 nM.
The statistic commonly available from experimental studies similar to Kim et al. (2006Kim et al. ( , 2007 and Ong et al. (2011a) is the coefficient of variation of the CFL width, CV w = σ w /w. Figure 4 shows how the mean WSSτ w , normalized with the smoothvessel smooth-vessel WSS τ ⋆ w , increases with CV w . (Recall that the fixed/smooth boundaries of the CFL correspond to CV w = 0 andτ w /τ ⋆ w = 1). The rate of growth of the mean WSS depends on the model's statistical parameters, some of which, especially l θ , are not found in the experiments (Kim et al., 2006(Kim et al., , 2007Ong et al., 2011a). Fortunately, Figure 4 reveals that the mean WSS is nearly insensitive to l θ , being dominated by the temporal fluctuations statistics that are more readily measurable.

NO Production Rate
It follows from Equation (10) that the rate of NO production by the endothelium,q 2 , is directly proportional to the WSS. When normalized by the control production rateq ⋆ , it is equal to the ratio τ w /τ ⋆ w . In other words, the statistics of the ratiosq 2 /q ⋆ and τ w /τ ⋆ w coincide. Therefore, Figure 4 also demonstrates how the mean NO production rate by the endothelium,q 2 /q ⋆ , increases with the coefficient of variation of the CFL width, CV w .

Mean Profiles of NO Concentration
Unless specified otherwise, the results reported below correspond to the hematocrit-dependent reaction rate k 1 in Equation (4) given by the constitutive model in Equations (8) and (9). We start by computing a (deterministic) reference NO concentration c ⋆ NO (r) as a solution of Equations (1)-(3) with smooth (constant) interfaces r 1 and r 2 . It serves as an initial condition for transient stochastic simulations.
The mean concentration profiles computed with these simulations,C NO (r), are exhibited in Figure 5. While the NO production rates (q 2 andq 3 ) on both sides of the endothelium (r = r 2 and r 3 ) are the same, the NO scavenging rate in the RBC core (0 ≤ r ≤ r 1 ) is higher than that in the muscle tissue (r > r 3 ). That is why the peak NO concentration is at the endothelium surface facing the tissue (r = r 3 ). Figure 5 also reveals that random fluctuations of the CFL width increase the NO availability relative to that predicted by FIGURE 4 | Mean WSS, normalized with the smooth-vessel WSS τ ⋆ w , (and mean NO production rate, normalized with the control production rateq ⋆ ) as a function of the coefficient of variation (CV w = σ w /w) of the CFL width for several combinations of constitutive statistical parameters. the model that ignores them. This is to be expected, since these fluctuations enhance the NO production by the endothelium (Figure 4). NO production and availability increase with the the degree of roughness of the random RBC-CFL interface r 1 (θ, t; ω): the higher CV w and/or the smaller the correlation lengths l t and l θ , the rougher the interface is.
The simulation results reported in Table 2 demonstrate the relative importance of temporal and angular fluctuations of the CFL width on NO availability. The latter is reported in terms of the ratio of the peak NO concentrations, . Larger values of R indicate stronger impact of the CFL width fluctuations.

Effect of Constitutive Models
The above-made estimates of NO production and availability rely on the NO scavenging rate k 1 (t; ω) given by the constitutive law in Equations (8), (9), which treats blood as a single-phase fluid. The alternative constitutive model for k 1 (t; ω), which explicitly accounts for the CFL presence, is given by Equations (6) and (7). Our simulations demonstrate that the difference between the mean peak NO concentrations predicted with the two models is less than 1% ( Table 2). This provides a confirmation of the robustness of our predictions of expected NO production and availability with respect to model selection for the scavenging rate.

Effect of Dextran Infusion
In the experiments reported by Ong et al. (2011a), infusion of a plasma expander dextran increases the average CFL width from w = 2.73 to 3.22 µm. It also enhances fluctuations of the CFL width, increasing CV w from 0.443 to 0.509 while leaving the correlation length l t practically unchanged (it increases from FIGURE 5 | Radial profile of mean NO concentration for several degrees of spatio-temporal variability of CFL quantified by CV w = 0.442 (l θ = 0.0 and σ θ = 0.0), CV w = 0.515 (l θ = 1.0 and σ θ = 0.07), and CV w = 0.582 (l θ = 1.0 and σ θ = 0.1); in all three cases, l t = 0.008 and σ t = 2.03. Also shown is NO concentration corresponding to constant uniform CFL width. The vertical lines indicate the inner and outer surfaces of the endothelium. 0.007 s before the dextran infusion to 0.008 s after). To match the decrease in the reported peak NO concentration from 11.2 to 9.5 nM, we recalculated the value of the pressure gradient J = 2.15 × 10 4 and 1.72 × 10 4 . Figure 6 demonstrates that changing the mean CFL width (from 2.73 to 3.22 µm) does not significantly change the mean WSS, but has a more pronounced effect on the mean peak NO concentrations. The peak NO concentration ratio R reported in Table 2 further emphasizes this effect.

DISCUSSION
We developed a computational framework to quantify the impact of spatio-temporal fluctuations in the CFL width on the production and transport of NO. This is accomplished by treating the RBC-CFL interface (and the corresponding CFL width) as a space-time correlated random field. This surface is represented via Karhunen-Loéve and Fourier expansions. The differential equations describing blood flow and NO production and transport, defined on random simulation domains, were solved by using a stochastic collocation method.
Our analysis demonstrates that the two-phase nature of microcirculatory flow, partitioned between a central blood column and a peripheral cell free plasma layer, causes stochastic flow variability in the CFL that influences NO bioavailability. Our findings are qualitatively comparable to those obtained with the deterministic analysis conducted by Ong et al. (2011a). However, they differ in an important way since the statical parameters on which they are based can be used to predict NO production, as well as other effects related to the variability of the CFL, with data from other experiments where the same statistical properties can be assessed.
Flow variability has a homeostatic role in the microcirculation where it is a factor in the control of blood flow and inflammation through biochemical mechanotransduction modulation of the production of NO and prostaglandins. Large Reynolds numbers and flow variability in the central circulation can cause the up-regulation of virtually all atherogenic or pro-inflammatory genes ultimately promoting the development of atherosclerotic plaques (Cabrales et al., 2011). Experimentally it is shown that a principal effect of flow (and shear stress) variability is the up and down regulation of the expression of a multitude of genes in endothelial cell cultures (Yee et al., 2008).
The rate of production of NO by mechanotransduction differs between steady and time-dependent flows, the latter resulting in a significantly higher production of NO. This phenomenon has been evidenced in studies using horizontal parallel plates flow chambers in which endothelial cells grow to a confluent layer (Ruel et al., 1995). Continuous changes of flow or "ramp" flow (Frangos et al., 1996) and step changes (Andrews et al., 2010) yielded significantly greater NO production than constant steady flow. Experimental studies also show that the rate of NO production is higher for sinusoidal flow vs. steady flow with the same mean (Noris et al., 1995;Li et al., 2005). There is evidence that flow variation frequency of about 1 Hz, independent of shear, is a determinant of the endothelial responses to pulsatile flow (Balcells et al., 2005).
Although vessel wall shear rates and shear stress (WSS) are similar throughout the circulation, as proposed by the "uniform shear stress hypothesis" (Kassab and Fung, 1995), flow variability is quite different. The main source of flow variability is the periodic action of the heart, which can also cause flow instabilities in locations with large Reynolds numbers. This variability decreases from the systemic blood vessels to the microcirculation where it is attenuated to an amplitude of 1-2% of mean flow (Intaglietta et al., 1971) in the arterioles.
Blood flow in microvessels also undergoes periodic flow changes whose amplitude can reach 100% of mean flow due to the phenomenon of vasomotion. This activity has fundamental frequency of the order of minutes for the larger arterioles, increasing as vessel size decreases (Colantuoni et al., 1984). This phenomenon has received many interpretations since it has been found in some normal conditions but not all, and is elicited by ischemic and low pressure states (Schmidt et al., 1992).
Our analysis and results show the importance of a microscopic random flow variability, whose primary effect is to increase the bioavailability of NO at the vessel/CFL interface. This variability has two components: a vessel wall component due to the unevenness of the endothelium (Barbee et al., 1994;Barbee, 2002) and the blood/CFL interface component. The latter is absent in experimental studies of cell culture flow chambers perfused with culture media, i.e., in the majority of studies.
There is presently no specific evidence of the effect on mediator production rate by stochastically induced mechanotransduction. However, since random excitation contains all the spectral components of sinusoidal, ramp and step flow variability it is likely that random excitation increases the rate of mediator production. These effects are largely unexplored, since endothelial mechanotransduction has been primarily studied with non-blood, steady or varying flow conditions on a scale that is macroscopic relative to the dimensions of endothelial cells (Chien, 2007). Yet, cell-scale variations of shear stress can cause large shear stress gradients (Davies et al., 2003;Leiderman et al., 2008). Experimental studies to determine the physiological role of this variability are not available.
Stochastic variability of WSS is fundamentally distinct from the temporal variability studied in cell cultures. The latter is macroscopic and has primarily a local effect in that it is absent in the microcirculation, where temporal flow changes occur in a time frame of minutes to hours, a condition of quasi steady state. The stochastic fluctuations of the CFL width generate a persistent time-dependent microscopic component of WSS, which is not tested by current experimental methods. Our results and these considerations suggest that endothelial flow chamber experiments aimed at understating the consequences of mechanotransduction on cardiovascular regulation should include the conditions that preserve the stochastic nature of the blood flow/microvessel wall interaction, i.e., use blood as the flow medium. This appears to be the only significant source of timedependent variability of shear stress in the microcirculation. It may be significantly affected by changes in the composition of either blood (due to hemorrhage and anemia) or plasma (possibly due to diabetes and hypertension).

AUTHOR CONTRIBUTIONS
SP carried out numerical implementation of the proposed algorithms, participated in the design of the study and drafted the manuscript; MI formulated the physiological problem, participated in the design of the study, and helped draft the manuscript; DT designed and coordinated the study, and helped draft the manuscript. All authors gave final approval for publication.