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.

1. 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 parallel-plate 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 non-uniform 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, 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).

2. Mathematical Model of NO Transport

2.1. 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\le r\le r_1$), the CFL ($E^2:r_1<r\le r_2$), the endothelial-cell region ($E^3:r_2<r\le r_3$), and the smooth-muscle region ($E^4:r_3<r\le r_4$). Stochastic fluctuations of the interface formed by flowing RBCs, r₁(θ, t), are modeled by treating it as a random function of both angular coordinate θ and time t, i.e., r₁ = r₁(θ, t; ω) with ω ∈ Ω indicating a realization (“coordinate”) in the probability space Ω. This renders the CFL width w = r₂ − r₁ random, i.e., w = w(θ, t; ω). Our goal is to capture the effects of stochastic fluctuations of the RBC-CFL interface r₁(θ, t; ω) on distribution of NO concentration, C_NO, in the Krogh tissue cylinder ${D}$ = {(r, θ):0 ≤ r ≤ r₄, 0 ≤ θ ≤ 2π}.

In each region of the computational domain, Eⁱ (i = 1, …, 4), the concentration C_NO satisfies a reaction-diffusion equation

$\begin{array}{lll}\hfill \frac{\partial C_{NO}}{\partial t}=D_i{\nabla }^2{C}_{NO}-k_i{C}_{NO},(r,\theta )\in E^i,& & (1)\end{array}$

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),

$\begin{array}{lll}\hfill C_{NO}^{-}=C_{NO}^{+},F_n^{-}-F_n^{+}={\dot{q}}_i,r=r_i.& & (2)\end{array}$

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, and ${\dot{q}}_i$ denotes the NO production rates at the interface r = r_i. Since only endothelium cells are involved in NO production, ${\dot{q}}_1\equiv 0$. We assume that no nitric oxide leaves the outer boundary of the smooth-muscle region, r = r₄, so that

$\begin{array}{lll}\hfill \mathbf{\text{n}}·\nabla C=0,r=r_4.& & (3)\end{array}$

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₁(θ, 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 $\mathbf{\text{u}}={(u_r,u_{\theta },V)}^{\top }$ 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).

2.2. 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₁) 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₁ corresponding to a given hematocrit level H_c is given by Ong et al. (2011b) and Chen et al. (2006)

$\begin{array}{lll}\hfill k_1=\frac{H_c}{H_s}{k}_s.& & (4)\end{array}$

The hematocrit ratio H_c/H_s is determined by mass conservation,

$\begin{array}{lll}\hfill {\displaystyle \underset{0}{\overset{2\pi }{\int }}}{\displaystyle \underset{0}{\overset{r_2}{\int }}}H(r,t)V(r,\theta ,t)rdrd\theta =H_s{\displaystyle \underset{0}{\overset{2\pi }{\int }}}{\displaystyle \underset{0}{\overset{r_2}{\int }}}V(r,\theta ,t)rdrd\theta .& & (5)\end{array}$

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₁(θ, 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; ω),

$\begin{array}{cc}\begin{array}{l}\frac{V}{V_{\text{max}}}=\{\begin{array}{ll}1-\frac{{\mu }_p}{{\mu }_c}\frac{r^2}{r_2^2}-\left(1-\frac{{\mu }_p}{{\mu }_c}\right)\frac{r_1^2}{r_2^2}\hfill & 0\le r\le r_1\hfill \\ \hfill & \hfill \\ 1-\frac{r^2}{r_2^2}\hfill & r_1\le r\le r_2\hfill \end{array}\hfill \end{array}& (6)\end{array}$

where $V_{max}=J{r}_2^2∕(4{\mu }_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

$\begin{array}{cc}\frac{H}{H_c}=\{\begin{array}{ll}1\hfill & 0\le r\le r_1\hfill \\ 0\hfill & r_1<r\le r_2\hfill \end{array}& (7)\end{array}$

enables one to compute the randomly fluctuating NO scavenging rate k₁(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₁(t; ω) was obtained from Equation (4).

The second alternative for obtaining k₁(t; ω) treats blood as a single-phase fluid with a parabolic velocity profile

$\begin{array}{lll}\hfill \frac{V}{V_{max}}=1-\frac{r^2}{r_2^2},0\le r\le r_2.& & (8)\end{array}$

This formulation replaces the CFL and the random RBC-CFL interface r₁(θ, t; ω) with a radial distribution of hematocrit,

$\begin{array}{cc}\frac{H}{H_c}=\{\begin{array}{ll}1\hfill & 0\le r\le r_1\hfill \\ \hfill & \hfill \\ {\left(\frac{r_2-r}{r_2-r_1}\right)}^2\hfill & r_1\le r\le r_2\hfill \end{array}.& (9)\end{array}$

Substituting Equations (8) and (9) into Equations (4) and (5) yields an alternative expression for the NO scavenging rate k₁(t; ω). This approach was used by Ong et al. (2011a) in the deterministic context that treated r₁(θ, t) as constant.

Finally, we allow the NO production rates by the endothelium, i.e., ${\dot{q}}_2$ and ${\dot{q}}_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 relation

$\begin{array}{lll}\hfill {\dot{q}}_2={\dot{q}}_3=\frac{{\tau }_w}{{\tau }_w^{\star }}{\dot{q}}^{\star },{\tau }_w={\mu }_p\frac{V_e}{w},& & (10)\end{array}$

where ${\tau }_w^{\star }$ is the reference wall shear stress, ${\dot{q}}^{\star }$ 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., ${\dot{q}}_2(t;\omega )$ and ${\dot{q}}_3(t;\omega )$, due to their dependence on the random flow velocity V and the CFL width w(θ, t; ω) = r₂ − r₁(θ, 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.

TABLE 1

Table 1. Model parameters and their values.

We represent spatio-temporal variations of the RBC-CFL interface,

$\begin{array}{lll}\hfill r_1(\theta ,t;\omega )=\left[{\overline{r}}_t+r_t^{\prime }(t;\omega )\right]\left[{\overline{r}}_{\theta }+r_{\theta }^{\prime }(\theta ;\omega )\right],& & (11)\end{array}$

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=\overline{r}+r^{\prime }$, as the sum of its ensemble mean $\overline{r}$ and zero-mean fluctuations r′. Setting ${\overline{r}}_{\theta }=1$ yields the mean and variance of the RBC-CFL interface: ${\overline{r}}_1={\overline{r}}_t$ and ${\sigma }_r^2={\overline{r}}_t^2{\sigma }_{\theta }^2+{\sigma }_t^2(1+{\sigma }_{\theta }^2)$. The coefficient of variation of the CFL width, $C{V}_w={\sigma }_w∕\overline{w}$, is given by

$\begin{array}{lll}\hfill C{V}_w^2={\left(\frac{{\overline{r}}_t}{\overline{w}}\right)}^2{\sigma }_{\theta }^2+C{V}_t^2(1+{\sigma }_{\theta }^2)& & (12)\end{array}$

where $\overline{w}=r_2-{\overline{r}}_1$ is the mean CFL width and σ_w is its standard deviation.

Since the random field $r_{\theta }^{\prime }(\theta ,\omega )$ is periodic, a truncated Fourier-type expansion

$\begin{array}{lll}\hfill r_{\theta }^{\prime }(\theta ;\omega )\approx {\sigma }_{\theta }{\displaystyle \sum _{n=-N_{\theta }}^{N_{\theta }}}{\nu }_n(\omega )e^{-in\theta }& & (13)\end{array}$

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 ${\sigma }_n^2=C_n∕4$, where

$\begin{array}{lll}\hfill C_n=\frac{1}{\pi }{\displaystyle \underset{0}{\overset{2\pi }{\int }}}{C}_{\theta }^{p}cos(n\theta )d\theta ,-N\le n\le N& & (14)\end{array}$

are coefficients of the Fourier cosine expansion of a 2π-periodic covariance function $C_{\theta }^p$ of the random field $r_{\theta }^{\prime }(\theta ,\omega )$. It is constructed as follows. First, we note that statistics of $r_{\theta }^{\prime }(\theta ;\omega )$ are rotationally invariant on the circle, such that a covariance function C_θ is

$\begin{array}{ll}\langle r_{\theta }^{\prime }({\theta }_1;\omega )r_{\theta }^{\prime }({\theta }_2;\omega )\rangle =C_{\theta }({\Delta }_{\theta }),{\Delta }_{\theta }=\left|{\theta }_1-{\theta }_2\right|.& (15)\end{array}$

Then $C_{\theta }^p$ is constructed by extending the covariance function C_θ of the random field $r_{\theta }^{\prime }(\theta ,\omega )$ to a 2π-periodic periodic domain. We employ a Gaussian covariance function $C_{\theta }({\Delta }_{\theta })=exp(-{\Delta }_{\theta }^2∕l_{\theta }^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^{\prime }(t;\omega )$ via a truncated Karhunen-Loéve expansion,

$\begin{array}{lll}\hfill r_t^{\prime }(t,\omega )={\sigma }_t{\displaystyle \sum _{m=1}^{N_t}}\sqrt{{\lambda }_m}{f}_m(t)Y_m(\omega ),& & (16)\end{array}$

where Y_m(ω) (m ≥ 1) are independent random variables, and λ_m and f_m(t) are, respectively, the eigenvalues and eigenfunctions of Fredholm equations,

$\begin{array}{lll}\hfill {\displaystyle \underset{0}{\overset{T}{\int }}}{\rho }_t(t,t^{\prime })f_m(t^{\prime })d{t}^{\prime }={\lambda }_m{f}_m(t),m\ge 1.& & (17)\end{array}$

For an exponential correlation function ${\rho }_t(t,t^{\prime })=exp(-|t-t^{\prime }|∕l_t)$ with the correlation length l_t > 0, the eigenvalue problems in Equation (17) admit an analytical solution (Lin et al., 2010),

$\begin{array}{lll}\hfill {\lambda }_m=\frac{2l_t}{l_t^2{\gamma }_m^2+1},f_m=\frac{l_t{\gamma }_{m}cos({\gamma }_{m}t)+sin({\gamma }_{m}t)}{\sqrt{(l_t^2{\gamma }_m^2+1)T∕2+l_t}}& & (18)\end{array}$

where γ_m are solutions of $(l_t^2{\gamma }^2-1)sin(\gamma T)=2l_t\gamma cos(\gamma 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 ${\sigma }_{\theta }^2$ and ${\sigma }_t^2$, 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.

TABLE 2

Table 2. Statistical parameters and summary of simulation results.

2.3. Numerical Solution

2.3.1. Mapping onto Deterministic Domain

We introduce a new coordinate system (ξ₁, ξ₂), which maps the original stochastic domain ${D}$ onto a rectangle B = {(ξ₁, ξ₂):−1 ≤ ξ₁ ≤ 1, 0 ≤ ξ₂ ≤ 2π}. A mapping ${D}={\cup }_{i=1}^4{E}^i\to B$ is accomplished analytically by the coordinate transformation

$\begin{array}{lll}\hfill r=r_{i-1}+\frac{{\xi }_1+1}{2}(r_i-r_{i-1}),\theta ={\xi }_2;(r,\theta )\in E^i& & (19)\end{array}$

where i = 1, …, 4 and r₀ = 0. The random RBC-CFL interface r₁(θ, t; ω) is represented by the expansions described above.

2.3.2. Transformed Stochastic Equations

The mapping defined by Equation (19) renders the transformation Jacobian

$\begin{array}{lll}\hfill J({\xi }_1,{\xi }_2,\omega )\equiv \frac{\partial (r,\theta )}{\partial ({\xi }_1,{\xi }_2)}=J\left[{\xi }_1,{\xi }_2,Y_1(\omega ),\cdots ,Y_K(\omega )\right].& & (20)\end{array}$

and other related metrics stochastic, i.e., dependent on a set of K = 2N_θ + N_t independent random variables ${\left\{Y_i(\omega )\right\}}_{i=1}^K$. The first N_t variables Y₁, …, Y_Nt coincide with those introduced in Equation (16) and the remaining 2N_θ variables represent their counterparts in Equation (13), such that Y_{Nt + 1} = ν_−Nθ, …, Y_K = ν_Nθ. Consequently, the deterministic reaction-diffusion Equations (1) are transformed into stochastic equations of the form (Appendix)

$\begin{array}{ll}\frac{\partial C_{NO}}{\partial t}={\displaystyle \sum _{i,j=1}^2}\frac{\partial }{\partial {\xi }_i}\left(D{A}_{ij}\frac{\partial C_{NO}}{\partial {\xi }_j}\right)-k{C}_{NO}& (21)\end{array}$

where the random coefficients A₁₁, A₁₂ = A₂₁, and A₂₂ 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).

3. 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 concentration ${\bar{C}}_{NO}$, serve as unbiased predictors of the system behavior; variances, e.g., NO concentration variance ${\sigma }_C^2$, act as a measure of predictive uncertainty.

3.1. 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₂ − r₁ in which the random RBC-CFL interface r₁(θ, 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₂) 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}_2^{1-\alpha -\beta }{W}^{\alpha -1}{(r_2-W)}^{\beta -1}$, where B(α, β) = Γ(α + β)∕[Γ(α)Γ(β)] is the beta function, Γ(·) is the complete gamma function, 0 ≤ W ≤ r₂, and α > 0 and β > 0 are shape parameters. Setting α = 4.358 and β = 32.9 provides the best data fit, resulting in the mean CFL width $\bar{w}=2.73$ μm. The auto-correlation data in Figure 2 were fitted with an exponential correlation function $\rho (t,t^{\prime })=exp(-|t-t^{\prime }|∕l_t)$, yielding the correlation length l_t = 0.007 s.

FIGURE 1

Figure 1. Frequency distribution of the temporally fluctuating CFL width w reported by Ong et al. (2011a) and the fitted β-distribution.

FIGURE 2

Figure 2. Auto-correlation of the temporally fluctuating CFL width w reported by Ong et al. (2011a) (solid line) and the fitted exponential correlation function ρ(Δt) = exp(−Δt∕l_t) with the correlation length l_t = 0.007 s (dashed line).

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 ${\sigma }_{\theta }^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₁(θ, t; ω) and CFL width w(θ, t; ω) = r₂ − r₁ increases with the variances ${\sigma }_{\theta }^2$ and ${\sigma }_t^2$, 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.

FIGURE 3

Figure 3. Realizations of temporal fluctuations of the CFL width w at angular coordinate θ = 0.0 corresponding to two sets of statistical parameters.

3.2. 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⁴ 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. (2006, 2007) and Ong et al. (2011a) is the coefficient of variation of the CFL width, ${CV}_w={\sigma }_w∕\bar{w}$. Figure 4 shows how the mean WSS ${\bar{\tau }}_w$, normalized with the smooth-vessel smooth-vessel WSS ${\tau }_w^{\star }$, increases with CV_w. (Recall that the fixed/smooth boundaries of the CFL correspond to CV_w = 0 and ${\bar{\tau }}_w∕{\tau }_w^{\star }=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, 2007; Ong 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.

FIGURE 4

Figure 4. Mean WSS, normalized with the smooth-vessel WSS ${\tau }_w^{\star }$, (and mean NO production rate, normalized with the control production rate ${\dot{q}}^{\star }$) as a function of the coefficient of variation (${\text{CV}}_w={\sigma }_w/\overline{w}$) of the CFL width for several combinations of constitutive statistical parameters.

3.3. NO Production Rate

It follows from Equation (10) that the rate of NO production by the endothelium, ${\dot{q}}_2$, is directly proportional to the WSS. When normalized by the control production rate ${\dot{q}}^{\star }$, it is equal to the ratio ${\tau }_w∕{\tau }_w^{\star }$. In other words, the statistics of the ratios ${\dot{q}}_2∕{\dot{q}}^{\star }$ and ${\tau }_w∕{\tau }_w^{\star }$ coincide. Therefore, Figure 4 also demonstrates how the mean NO production rate by the endothelium, ${\stackrel{°}{\stackrel{̄}{q}}}_2∕{\dot{q}}^{\star }$, increases with the coefficient of variation of the CFL width, CV_w.

3.4. Mean Profiles of NO Concentration

Unless specified otherwise, the results reported below correspond to the hematocrit-dependent reaction rate k₁ in Equation (4) given by the constitutive model in Equations (8) and (9). We start by computing a (deterministic) reference NO concentration $c_{NO}^{\star }(r)$ as a solution of Equations (1)–(3) with smooth (constant) interfaces r₁ and r₂. It serves as an initial condition for transient stochastic simulations.

The mean concentration profiles computed with these simulations, ${\bar{C}}_{NO}(r)$, are exhibited in Figure 5. While the NO production rates (${\dot{q}}_2$ and ${\dot{q}}_3$) on both sides of the endothelium (r = r₂ and r₃) are the same, the NO scavenging rate in the RBC core (0 ≤ r ≤ r₁) is higher than that in the muscle tissue (r > r₃). That is why the peak NO concentration is at the endothelium surface facing the tissue (r = r₃).

FIGURE 5

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.

Figure 5 also reveals that random fluctuations of the CFL width increase the NO availability relative to that predicted by 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₁(θ, 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, ${R}=({\overline{C}}_{\text{max}}-C_{\text{max}}^{\star })/C_{\text{max}}^{\star }$, where $C_{max}^{\star }=C_{NO}^{\star }(r_3)$ and ${\bar{C}}_{max}={\bar{C}}_{NO}(r_3)$. Larger values of ${R}$ indicate stronger impact of the CFL width fluctuations.

3.5. Effect of Constitutive Models

The above-made estimates of NO production and availability rely on the NO scavenging rate k₁(t; ω) given by the constitutive law in Equations (8), (9), which treats blood as a single-phase fluid. The alternative constitutive model for k₁(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.

3.6. 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 $\bar{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 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⁴ and 1.72 × 10⁴.

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.

FIGURE 6

Figure 6. Mean peak NO concentration (A) and mean WSS ratio (B) corresponding to temporal and spatial variations of CFL whose mean value is $\bar{w}=2.73$ μm or 3.22 μm.

4. 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 time-dependent 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.

Funding

This research was supported in part by NIH under award numbers R01-HL064395 and R24-HL64395, Air Force Office of Scientific Research under award number DE-FG02-07ER25815, and by National Science Foundation under award number DMS-1522799.

Conflict of Interest Statement

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.

References

Andrews, A. M., Jaron, D., Buerk, D. G., Kirby, P. L., and Barbee, K. A. (2010). Direct, real-time measurement of shear stress-induced nitric oxide produced from endothelial cells in vitro. Nitric Oxide 23, 335–342. doi: 10.1016/j.niox.2010.08.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Balcells, M., Fernández Suárez, M., Vázquez, M., and Edelman, E. R. (2005). Cells in fluidic environments are sensitive to flow frequency. J. Cell Physiol. 204, 329–335. doi: 10.1002/jcp.20281

PubMed Abstract | CrossRef Full Text | Google Scholar

Barbee, K. A. (2002). Role of subcellular shear-stress distributions in endothelial cell mechanotransduction. Ann. Biomed. Eng. 30, 472–482. doi: 10.1114/1.1467678

PubMed Abstract | CrossRef Full Text | Google Scholar

Barbee, K. A., Davies, P. F.,and Lal R. (1994). Shear stress-induced reorganization of the surface topography of living endothelial cells imaged by atomic force microscopy. Circ. Res. 74, 163–171. doi: 10.1161/01.RES.74.1.163

PubMed Abstract | CrossRef Full Text | Google Scholar

Cabrales, P., Zanini, G. M., Meays, D., Frangos, J. A., and Carvalho, L. J. (2011). Nitric oxide protection against murine cerebral malaria is associated with improved cerebral microcirculatory physiology. J. Infect. Dis. 203, 1454–1463. doi: 10.1093/infdis/jir058

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, X., Jaron, D., Barbee, K. A., and Buerk, D. G. (2006). The influence of radial RBC distribution, blood velocity profiles, and glycocalyx on coupled NO/O2 transport. J. Appl. Physiol. 100, 482–492. doi: 10.1152/japplphysiol.00633.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

Chien, S. (2007). Mechanotransduction and endothelial cell homeostasis: the wisdom of the cell. Am. J. Physiol. Heart Circ. Physiol. 292, H1209–H1224. doi: 10.1152/ajpheart.01047.2006

PubMed Abstract | CrossRef Full Text | Google Scholar

Colantuoni, A., Bertuglia, S.,and Intaglietta M. (1984). Quantitation of rhythmic diameter changes in arterial microcirculation. Am. J. Physiol. 246, H508–H517.

PubMed Abstract | Google Scholar

Condorelli, P., and George, S. C. (2002). Free nitric oxide diffusion in the bronchial microcirculation. Am. J. Physiol. Heart. Circ. Physiol. 283, 2660–2670. doi: 10.1152/ajpheart.00003.2002

PubMed Abstract | CrossRef Full Text | Google Scholar

Davies, P. F., Zilberberg, J., and Helmke, B. P. (2003). Spatial microstimuli in endothelial mechanosignaling. Circ. Res. 92, 359–370. doi: 10.1161/01.RES.0000060201.41923.88

PubMed Abstract | CrossRef Full Text | Google Scholar

Dupin, M. M., Halliday, I., Care, C. M., Alboul, L., and Munn, L. L. (2007). Modeling the flow of dense suspensions of deformable particles in three dimensions. Phys. Rev. E 75:066707. doi: 10.1103/PhysRevE.75.066707

PubMed Abstract | CrossRef Full Text | Google Scholar

Frangos, J. A., Huang, T. Y., and Clark, C. B. (1996). Steady shear and step changes in shear stimulate endothelium via independent mechanisms – superposition of transient and sustained nitric oxide production. Biochem. Biophys. Res. Commun. 224, 660–665. doi: 10.1006/bbrc.1996.1081

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghanem, R., and Spanos, P. (1991). Stochastic Finite Elements: A Spectral Approach (New York, NY: Springer-Verlag).

Google Scholar

Hightower, C. M. Intaglietta M. (2008). Diagnostic frequency continuous ultrasonography directly mitigates venular ischemia reperfusion damage. J. Am. Coll. Surg. 206, 540–547. doi: 10.1016/j.jamcollsurg.2007.09.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Intaglietta, M., Richardson, D. R., and Tompkins, W. R. (1971). Blood pressure, flow, and elastic properties in microvessels of cat omentum. Am. J. Physiol. 221, H922–H928.

PubMed Abstract | Google Scholar

Kanai, A. J., Strauss, H. C., Truskey, G. A., Crews, A. L., Grunfeld, S.,and Malinski T. (1995). Shear stress induces ATP-independent transient nitric oxide release from vascular endothelial cells, measured directly with a porphyrinic microsensor. Circ. Res. 77, 284–293. doi: 10.1161/01.RES.77.2.284

PubMed Abstract | CrossRef Full Text | Google Scholar

Kassab, G. S., and Fung, Y. C. (1995). The pattern of coronary arteriolar bifurcations and the uniform shear hypothesis. Ann. Biomed. Eng. 23, 13–20. doi: 10.1007/BF02368296

PubMed Abstract | CrossRef Full Text | Google Scholar

Kavdia, M., and Popel, A. S. (2003). Wall shear stress differentially affects NO level in arterioles for volume expanders and Hb-based O2 carriers. Microvasc. Res. 66, 49–58. doi: 10.1016/S0026-2862(03)00008-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, S., Kong, R. L., Popel, A. S., Intaglietta, M., and Johnson, P. C. (2006). A computer-based method for determination of the cell-free layer width in microcirculation. Microcirculation 13, 199–207. doi: 10.1080/10739680600556878

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, S., Kong, R. L., Popel, A. S., Intaglietta, M., and Johnson, P. C. (2007). Temporal and spatial variations of cell-free layer width in arterioles. Am. J. Physiol. Heart Circ. Physiol. 293, H1526–H1535. doi: 10.1152/ajpheart.01090.2006

PubMed Abstract | CrossRef Full Text | Google Scholar

Kuo, L., Chilian, W. M., and Davis, M. J. (1990). Coronary arteriolar myogenic response is independent of endothelium. Circ. Res. 66 860–866. doi: 10.1161/01.RES.66.3.860

PubMed Abstract | CrossRef Full Text | Google Scholar

Lamkin-Kennard, K. A., Jaron, D., and Buerk, D. G. (2004). Impact of the Fahraeus effect on NO and O2 biotransport: a computer model. Microcirculation 11, 337–349. doi: 10.1080/10739680490437496

PubMed Abstract | CrossRef Full Text | Google Scholar

Leiderman, K. M., Miller, L. A., and Fogelson, A. L. (2008). The effects of spatial inhomogeneities on flow through the endothelial surface layer. J. Theor. Biol. 252, 313–325. doi: 10.1016/j.jtbi.2008.01.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, Y., Zheng, J., Bird, I. M., and Magness, R. R. (2005). Effects of pulsatile shear stress on signaling mechanisms controlling nitric oxide production, endothelial nitric oxide synthase phosphorylation, and expression in ovine fetoplacental artery endothelial cells. Endothelium 12, 21–39. doi: 10.1080/10623320590933743

PubMed Abstract | CrossRef Full Text | Google Scholar

Lin, G., Tartakovsky, A. M., and Tartakovsky, D. M. (2010). Uncertainty quantification via random domain decomposition and probabilistic collocation on sparse grids. J. Comput. Phys. 229, 6995–7012. doi: 10.1016/j.jcp.2010.05.036

CrossRef Full Text | Google Scholar

Martini, J., Carpentier, B., Negrete, A. C., Frangos, J. A.,and Intaglietta M. (2005). Paradoxical hypotension following increased hematocrit and blood viscosity. Am. J. Physiol. Heart Circ. Physiol. 289, 2136–2143. doi: 10.1152/ajpheart.00490.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

Namgung, B., Ong, P. K., Johnson, P. C., and Kim, S. (2011). Effect of cell-free layer variation on arteriolar wall shear stress. Ann. Biomed. Eng. 39, 359–366. doi: 10.1007/s10439-010-0130-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Noris, M., Morigi, M., Donadelli, R., Aiello, S., Foppolo, M., Todeschini, M., et al. (1995). Nitric oxide synthesis by cultured endothelial cells is modulated by flow conditions. Circ. Res. 76, 536–543. doi: 10.1161/01.RES.76.4.536

PubMed Abstract | CrossRef Full Text | Google Scholar

Ong, P. K., Jain, S.,and Kim S. (2011a). Temporal variations of the cell-free layer width may enhance NO bioavailability in small arterioles: effects of erythrocyte aggregation. Microvasc. Res. 81, 303–312. doi: 10.1016/j.mvr.2011.02.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Ong, P. K., Jain, S., and Kim, S. (2011b). Modulation of NO bioavailability by temporal variation of the cell-free layer width in small arterioles. Ann. Biomed. Eng. 39, 1012–1023. doi: 10.1007/s10439-010-0216-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Park, S. W., Intaglietta, M., and Tartakovsky, D. M. (2012). Impact of endothelium roughness on blood flow. J. Theor. Biol. 300, 152–160. doi: 10.1016/j.jtbi.2012.01.017

PubMed Abstract | CrossRef Full Text | Google Scholar

Ruel, J., Lemay, J., Dumas, G., Doillon, C.,and Charara J. (1995). Development of a parallel plate flow chamber for studying cell behavior under pulsatile flow. ASAIO J. 41, 876–883. doi: 10.1097/00002480-199541040-00011

PubMed Abstract | CrossRef Full Text | Google Scholar

Sato, M., Nagayama, K., Kataoka, N., Sasaki, M.,and Hane K. (2000). Local mechanical properties measured by atomic force microscopy for cultured bovine endothelial cells exposed to shear stress. J. Biomech. 33, 127–135. doi: 10.1016/S0021-9290(99)00178-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Schmidt, J. A., Intaglietta, M., and Borgström, P. (1992). Periodic hemodynamics in skeletal muscle during local arterial pressure reduction. J. Appl. Physiol. 73, 1077–1083.

PubMed Abstract | Google Scholar

Sharan, M., Singh, B.,and Kumar P. (1997). A two-layer model for studying the effect of plasma layer on the delivery of oxygen to tissue using a finite element method. Appl. Math. Modell. 21, 419–426. doi: 10.1016/S0307-904X(97)00035-8

CrossRef Full Text | Google Scholar

Sriram, K., Vázquez, B. Y. S, Yalcin, O., Johnson, P. C., Intaglietta, M., and Tartakovsky, D. M. (2011). The effect of small changes in hematocrit on nitric oxide transport in arterioles. Antioxid. Redox Signal. 14, 175–185. doi: 10.1089/ars.2010.3266

PubMed Abstract | CrossRef Full Text | Google Scholar

Tartakovsky, D. M., and Winter, C. L. (2001). Dynamics of free surfaces in random porous media. SIAM J. Appl. Math. 61, 1857–1876. doi: 10.1137/S0036139999358180

CrossRef Full Text | Google Scholar

Tartakovsky, D. M. Xiu D. (2006). Stochastic analysis of transport in tubes with rough walls. J. Comput. Phys. 217, 248–259. doi: 10.1016/j.jcp.2006.02.029

CrossRef Full Text | Google Scholar

Tsai, A. G., Friesenecker, B., McCarthy, M., Sakai, H.,and Intaglietta M. (1998). Plasma viscosity regulates capillary perfusion during extreme hemodilution in hamster skin fold model. Am. J. Physiol. Heart. Circ. Physiol. 275, H2170–H2180.

PubMed Abstract | Google Scholar

Vaughn, M. W., Kuo, L., and Liao, J. C. (1998b). Effective diffusion distance of nitric oxide in the microcirculation. Am. J. Physiol. Heart Circ. Physiol. 274, 1705–1714.

PubMed Abstract | Google Scholar

Vaughn, M. W., Kuo, L., and Liao, J. C. (1998a). Estimation of nitric oxide production and reaction rates in tissue by use of a mathematical model. Am. J. Physiol. Heart Circ. Physiol. 274, 2163–2176.

PubMed Abstract | Google Scholar

Yalcin, O., Ulker, P., Yavuzer, U., Meiselman, H. J., and Baskurt, O. K. (2008). Nitric oxide generation by endothelial cells exposed to shear stress in glass tubes perfused with red blood cell suspensions: role of aggregation. Am. J. Physiol. Heart Circ. Physiol. 294, H2098–H2105. doi: 10.1152/ajpheart.00015.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

Yee, A., Bosworth, K. A., Conway, D. E., Eskin, S. G., and McIntire, L. V. (2008). Gene expression of endothelial cells under pulsatile non-reversing vs. steady shear stress: comparison of nitric oxide production. Ann. Biomed. Eng. 36, 571–579. doi: 10.1007/s10439-008-9452-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, J., Johnson, P. C., and Popel, A. S. (2009). Effects of erythrocyte deformability and aggregation on the cell free layer and apparent viscosity of microscopic blood flows. Microvas. Res. 77, 265–272. doi: 10.1016/j.mvr.2009.01.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Appendix A

Transformed Transport Equations

A general transformation of coordinates ξ_i(x₁, x₂) transforms the Fickian flux F(x) = −D∇_xC_NO into a vector F(ξ) whose contravariant components F¹ and F² are given by

$\begin{array}{ll}{F}^i=D\frac{1}{J}{\displaystyle \sum _{k=1}^2}{A}_{ik}\frac{\partial C_{NO}}{\partial {\xi }_k},i=1,2& (\mathrm{A1})\end{array}$

where J is the transformation Jacobian in Equation 20 and

$\begin{array}{llll}\hfill A_{11}=& {\left(\frac{\partial x_2}{\partial {\xi }_2}\right)}^2+{\left(\frac{\partial x_1}{\partial {\xi }_2}\right)}^2,& & (\mathrm{A2a})\end{array}$

$\begin{array}{llll}\hfill A_{12}=& A_{21}=-\frac{\partial x_2}{\partial {\xi }_1}\frac{\partial x_2}{\partial {\xi }_2}-\frac{\partial x_1}{\partial {\xi }_1}\frac{\partial x_1}{\partial {\xi }_2},& & (\mathrm{A2b})\end{array}$

$\begin{array}{llll}\hfill A_{22}=& {\left(\frac{\partial x_2}{\partial {\xi }_1}\right)}^2+{\left(\frac{\partial x_1}{\partial {\xi }_1}\right)}^2.& & (\mathrm{A2c})\end{array}$

Furthermore,

$\begin{array}{lll}\hfill {\nabla }_{\mathbf{\text{x}}}·\mathbf{\text{F}}=\frac{1}{J}\left(\frac{\partial F^1}{\partial {\xi }_1}+\frac{\partial F^2}{\partial {\xi }_2}\right)\equiv \frac{1}{J}{\nabla }_{\xi }·\mathbf{\text{F}}.& & (\mathrm{A3})\end{array}$

Substituting Equations (A1)–(A3) into Equation (1) leads to Equation (21).

Appendix B

Numerical Implementation

Let ${\left\{{\text{Ψ}}_m(\mathbf{\text{Y}})\right\}}_{m=0}^M$ denote a set of multidimensional orthogonal polynomials of the random vector $\mathbf{\text{Y}}(\omega )\equiv {(Y_1,\dots ,Y_K)}^T$. The polynomials are chosen to have the ensemble means ${\overline{\text{Ψ}}}_0=1$ and ${\overline{\text{Ψ}}}_k=0$ (k ≥ 1) and to satisfy the orthogonality condition $\langle {\text{Ψ}}_i{\text{Ψ}}_j\rangle =\langle {\text{Ψ}}_i^2\rangle {\delta }_{ij}$, where the 〈·〉 operator is defined by

$\begin{array}{lll}\hfill \langle {\text{Ψ}}_i{\text{Ψ}}_j\rangle \equiv \int {\text{Ψ}}_i(\mathit{\text{Y}}){\text{Ψ}}_j(\mathit{\text{Y}})W(\mathit{\text{Y}})d{Y}_1\dots d{Y}_K,& & (\mathrm{A4})\end{array}$

δ_ij is the Kronecker delta, and W(Y) is a weight function corresponding to a given polynomial type. The size of the polynomial set, M, is determined by the “stochastic dimension” K and the order P of polynomials Ψ_k, according to

$\begin{array}{ll}M=\frac{(K+P)!}{K!P!}-1.& (\mathrm{A5})\end{array}$

Polynomial chaos expansions represent a system state, e.g., u(ξ, ω), a random field whose ensemble statistics are to be determined, as a series

$\begin{array}{lll}\hfill u(\xi ,\omega )={\displaystyle \sum _{k=0}^M}{û}_k(\xi ){\text{Ψ}}_k\left[\mathit{\text{Y}}(\omega )\right].& & (\mathrm{A6})\end{array}$

Instead of stochastic Galerkin method, we adopt a stochastic collocation method to solve the random boundary problem. This method is described below.

Appendix C

Stochastic Collocation Method

Consider a stochastic partial differential equation

$\begin{array}{cc}{L}(x,\xi (\omega );u)=f(x,\xi (\omega )).& (\text{A}7)\end{array}$

Its solution is approximated by using a Lagrange formula

$\begin{array}{lll}\hfill u(\mathit{\text{x}},\xi (\omega ))={\displaystyle \sum _{k=0}^{Nq}}{u}_k(\mathit{\text{x}},{\xi }^k){{L}}_k(\xi ),& & (\mathrm{A8})\end{array}$

where $u_k(x,{\xi }^k)$ is the solution at the set of collocation points ${\left\{{\xi }^k\right\}}_{k=0}^{N_q}$, L_k is the Lagrange polynomial of order N_q + 1, and L_i(ξ_j) = δ_ij. Next, Equation (A8) is substituted into the inner product formula given by Equation (A4) with the basis functions chosen to be Dirac delta functions δ(ξ − ξ^k). Then N_q + 1 uncoupled deterministic problems are solved at the collocation points,

$\begin{array}{cc}L(x,{\xi }^k;u_k)=f(x,{\xi }^k),\text{for}k=0,\dots ,N_q.& (\text{A}9)\end{array}$

Once the solutions at a set of collocation points are obtained, their statistics are computed by using the corresponding quadrature rule,

$\begin{array}{lll}\hfill \overline{u}(\mathit{\text{x}})={\displaystyle {\int }_{\Gamma }}u(\mathbf{\text{x}},\xi )p(\xi )d\xi \approx {\displaystyle \sum _{k=1}^{Nq}}u(\mathbf{\text{x}},{\xi }^k)w_k& & (\mathrm{A10})\end{array}$

and

$\begin{array}{lll}\hfill {\sigma }_u^2(\mathit{\text{x}})={\displaystyle {\int }_{\Gamma }}{\left[u^{\prime }(\mathit{\text{x}},\xi )\right]}^{2}p(\xi )d\xi \approx {\displaystyle \sum _{k=1}^{Nq}}u{(\mathit{\text{x}},{\xi }^k)}^2{w}_k-{\overline{u}}^2& & (\mathrm{A11})\end{array}$

where ${\left\{w_k\right\}}_{k=1}^{Nq}$ is a set of weights corresponding to the set of quadrature points ${\left\{{\xi }^k\right\}}_{k=0}^{Nq}$.