Numerical solution of airway closure via airway smooth muscle interaction
The results presented in this chapter have appeared in (65).
In this chapter, we rst present a detailed analysis of the comparison between the DM approximation (as developed in the Chapter 2) and the full PDE based (MoC) solution of the ASM cross-bridge system as in (58) and described in section (1.2). This is done under length uctuations, simulating the dynamically changing cross-bridge distributions that occur in ASM during breathing, as performed in Mijailovich et al (58). This serves as a foundation and measure of the accuracy of the DM in capturing ASM dynamics. We then proceed to further test the validity of the DM approximation by making appropriate comparisons to the coupled ASM-airway wall system for both full PDE-based and DM approximated systems. In making this comparison, we are aided by several aspects of the problem: rst, that airways are typically subject to an oscillatory pressure environment, and second, that airway bistability (1, 2, 16) can provide a natural metric by discriminat-ing between \open » and \close » (or e ectively closed) airway states. In this way, we are able to establish regions of validity for the DM approximation to the coupled airway-ASM system. We show that most of the parameter space exhibits good quantitative agreement between the DM approximation model and its full PDE-based counterpart. Furthermore, although signi cant quantitative deviations do begin to appear when the airway is sub-jected to high-amplitude, high-frequency oscillations near airway closure, it is important to note that the qualitative behaviour is maintained, suggesting that a relatively simple correction factor might render the DM approximation acceptable throughout parameter space. Thus, the DM approximation to cross-bridge dynamics within a coupled ASM-airway system provides a viable alternative to the solution of the explicit cross-bridge PDEs and an avenue for future inclusion of (systematically reduced) ASM dynamics in lung-scale interacting models.
Numerical solution of Mijailovich’s crossbridge model
Figure 3.1 shows the comparison between the full Mijailovich cross-bridge PDE model and the DM approximation for the imposed length oscillations. Panels A and B show the bulk observables of force, sti ness and phosphorylation for each model whereas Panels C and D give a more detailed breakdown by cross-bridge state. While di erences are certainly in evidence, broadly speaking the ASM-only model agreement suggests that relatively modest di erences in ASM during forced oscillations, particularly in the form and struc-ture of the myosin bond populations as well as the macroscopic muscle behaviour, are ampli ed by non-linear interactions with the airway wall into more complex di erences, especially in terms of frequency dependence (and at higher oscillation frequencies/ampli-tudes). There is a rapid increase in phosphorylation, force and sti ness during the rst 0-5s, which coincides with early Ca2+ transitions. This increase also corresponds to an increase in the AM and AM p cross-bridges. The increased in force and sti ness slowed whilst phosphorylation decreased and the system obtained equilibrium before the impo-sition of length uctuations. The DM also attained a lower maximum sti ness than the MoC as a result of a lower second moment of the sum of the AM and AM p populations. These oscillatory perturbations disrupt the isometric binding equilibrium and causes the system to adapt to a new perturbed equilibrium. The new steady state is characterised by fewer attached bridges, lower mean force and lower sti ness particularly at 4% length con-tractions. These e ects are caused by an increase of cross-bridges when they are exposed to higher rates of detachment when convected out of the binding region 0 < x < h.
Comparison between the full model (the PDE cross-bridge-based system) and its DM-based approximation is complicated by the conditions imposed on the system. For ex-ample, the accuracy of the DM approximation might be very di erent in an iso-pressure contraction of a large, cartilaginous central airway, as opposed to contraction of a small, peripheral airway under periodic breathing. If one allows an arbitrary external pres-sure waveform, the space over which a comparison must be made is, in fact, of in nite dimension.
In order to simplify this comparison to workable form, we make two key reductions:
We assume a three-parameter external pressure waveform consisting of triangular waves; thus the frequency, amplitude, and minimum pressure (three parameters) are su cient to de ne this space. This is designed to mimic both common experimental protocols, and also in vivo considerations.
Rather than attempting to quantify the entire timeseries solutions to the system, we de ne a key measure based on the relevant physiological concerns. These are that, rst, the airway luminal radius is the primary quantity of interest; and second, that the existence of airway bistability can be exploited. Namely, for a given airway size, degree of activation, and xed external pressure waveform, an airway can respond in only three qualitatively distinct ways: A) remain in the open state; B) close, and remain in the closed state; or C) switch between open and closed states inde nitely. The basic idea is illustrated in Figure 2, showing the underlying bistability (green) as well as sample trajectories for each model (red-MoC, blue-DM).
Using this latter idea we de ne as the degree of activation (e.g. the value of the force activation parameter ) at which the A ! B transition occurs (e.g. airway closure). This critical transition value then is itself a function of the external pressure waveform, e.g. it is a function of the external pressure frequency, amplitude, and minimum external pressure. We then de ne a range of these parameters of physiological interest, and com-pute over this parameter space for both the full model and the DM approximation; by then comparing these we are able to understand the situations under which the DM approximation is most accurate, and those in which deviations are larger.
In order to determine , we consider a threshold radius (denoted r) within the unstable region of the open to close transition (see Figure (3.2)) and then used bisection method to determine the value of required to constrict the airway lumen past r, within the set integration time.
Inde nite switching between open and closed states occurs when the airway is subject to large amplitude, slow frequency oscillations and minimal . At these values the radius of the airway lumen is reduced past r but the force provided by the ASM is not su cient to keep the airway at the closed state and it thus returns to the opened state. Because this occurs in only a very small portion of the physiological parameter space, we concentrate here on the (A ! B) transition.
Results and discussion
Static pressure conditions
The simplest case for comparison between the distribution-moment approximation model and its full cross-bridge counterpart is simply to narrow the airway against a static external pressure. As described in the previous section, we compute as a function of this imposed pressure for both models and compare these critical closure values; the results are given in Figure 3.3 with given for DM ( ) and MoC ( ). In this simplest case, both models show good agreement both with expectations and with one another. Speci cally, one would qualitatively expect monotonic increase in as the external pressure increases, as shown by both models. Furthermore, DM shows excellent quantitative agreement with the full model with a 2% to 5% percentage di erence as well as having absolute di erences in the ranges of 0 to 1.2. This suggests that the DM approximation is suitable at least for this simplest, static case; however, as the entire protocol is quasi-isometric, one might have expected that accurate estimation of the rst bond-moment is all that is necessary here. Further comparison of the bound state distributions is given in panels (B) and (C) for both low and high static external pressure, respectively. In (85), while the distributions themselves showed signi cant departures, the moments agreed rather better; here a similar observation can be made regarding the relatively large deviation in the distributions versus the rather better performance of the coupled system (in terms of ).
External pressure oscillations
A more realistic and challenging scenario is to impose external pressure oscillations, mimicking either the breathing environment or common experimental protocols. Hence we impose triangular pressure oscillations deter-mined by three parameters: amplitude, frequency, and minimal external pressure. For xed pressure oscillations, then, the critical value can be found in that environment, for each model, and again compared. The central idea is illustrated in Figure 3.4, which shows individual radius trajectories for each model for xed external pressure oscilla-tions; the transition from trajectories which remain in the open state, to those which transition to closed, occurs at . To characterize this relationship for not just one xed pressure waveform, but more generally, we compute for a range of pressure waveforms within the 3-parameter space famplitudeg ffrequencyg fminimal external pressureg throughout the physiological range. The percentage di erence between the values (DM vs MoC) is given in the main panel of Figure 3.5 chosen at an amplitude of 25 so as to illustrate maximum percentage di erences. Broadly speaking, there is good qualitative agreement between the DM and MoC models for much of this parameter space, with differences of < 10% existing for broad regions; in particular for larger minimal pressures, lower amplitudes, and lower frequencies. On the other hand, di erences of 40% appear for high-frequency, large-amplitude oscillations with low minimal pressures. In order to better illustrate the behavior of these di erences, we also show selected \slices » of the parameter space in the sub-panels A-F (Figure 3.5); the location of each \slice » is also annotated in the main panel.
Several observations are apparent. First, parametric dependence on minimal external pressure (slices A & D) and oscillation amplitude (slices B & E) are relatively simple and in good qualitative agreement, even within the region of relatively large quantita-tive departure (25% – 40% relative error). On the other hand, the frequency response (slices C & F) is more complex, being both non-linear and non-monotonic with relative percentage di erences ranging from 5% to 40% in C and 5% to 25% in F. It is worth a mention that although the percentage di erences are relatively large, the absolute error between predicted values at the largest deviation is 5 and a value of 0.5 at the lowest percentage di erences. Even so, because of the persistent qualitative agreement through parameter space, these di erences can be reduced to < 15% everywhere with a relatively simple second-order polynomial correction factor
Numerical e ciency
In order to compare the e ciency of the DM to the MoC, we considered the number of ODEs required to be solved by each method. For the DM method, this is xed: 7 ODEs. For the MoC approach, one must choose the spatial discretization; speci cally, the N points over the spatial domain. We consider N = 2000 as our \gold standard » solution. The key idea is that one might also reduce the computational cost by reducing N , though doing so introduces additional error. Thus our approach to quantifying the numerical e ciency is to determine the value of N at which DM achieves the same error (relative to the gold standard). In particular, we determined the maximum absolute error in the force produced between the MoC at N = 2000 (gold standard) and the DM approximation. The number of ODEs for the MoC was then systematically reduced and at each point, the maximum absolute error between this and the MoC at gold standard was found. The number of ODEs for which the error in the MoC is approximately equal to the error between the DM and gold standard MoC, was then considered equivalent to the DM approximation. This process is represented in Fig 3.6. For these parameter values, the MoC requires 690 = 3N ODEs compared with 7 ODEs required by the DM approximation to obtain the same relative error. These parameter values were chosen as they represent a reasonable mid-range value in the main panel of Fig 3.5. Other points in parameter space were similarly explored for numerical e ciency and displayed similar improvements in DM vs MoC. Thus using this approach, the DM method is approximately 100 times more e cient numerically than an equivalent MoC scheme. Additionally it should be noted that when numerically integrating these models, both the DM and MoC methods display loss of order when solved using Runge-Kutta methods, and so must be carefully convergence tested.
In this chapter we presented results associated with the application, analysis and com-parison of the full PDE based cross-bridge model as well as its DM-based approximation, to the coupled ASM-airway wall system. By making suitable comparisons between the distribution-moment model, and the full coupled cross-bridge-airway model, we show that:
the DM approximation is accurate when static external pressures are imposed; 2) for much of the oscillatory waveform parameter space, the DM approximation again provides a satisfactory approximation, with the possible exception of high-frequency, large ampli-tude oscillations with low minimal pressures; 3) even within this region of relatively large quantitative departure, persistent qualitative agreement allows a relative simple correc-tion factor to reduce relative errors to less than 15% throughout the physiological range. Taken together, this suggests that the DM approximation provides a viable method for vast reduction in computational complexity, while retaining key ASM dynamics through a systematic reduction.
Rate of ASM force recovery post deep inspiration.
A deep inspiration (DI) is a breath taken from functional residual capacity (volume of air in lungs at the end of a full expiration) to total lung capacity (TLC). It has been the subject of extensive research and many of these studies have focused on the di erences in the responses between asthmatic and non asthmatic patients. Speci cally, DIs are known to be e ective in reversing bronchoconstriction in healthy individuals, but may fail to prevent bronchoconstriction in asthmatic patients. Asthmatic patients also show limited bronchodilation of constricted airways as compared with normals (22, 66). Fish et al (23) in their study of response of patients with asthma and allergic rhinitis (irri-tation and in ammation of the mucous membrane inside the nose), rst suggested that important clues of the mechanism of airway hyperresponsiveness (increased sensitivity of the airways to an inhaled constrictor) may be found in the di erent response to DI experienced by asthmatics and non-asthmatics. They studied the response of patients to inhaled methacholine (MCh) before and after showed that asthmatic subjects exhibits greater airway closure than non-asthmatics. The works of Skloot et al (71) con rmed these results, and showed that without periodic DI, normal subjects resembled asthmat-ics, in that they appear to be hyperresponsive. Additionally, King et al (47) reported that preventing DI during a complete dose of constrictor accentuates airway narrowing in normals. The results associated with these studies seems to suggest that preventing ASM from stretching for a signi cant period of time (20-30 minutes), increases its ability to shorten and decreases the e ectiveness of subsequent DIs (78). Mechanically during a DI, a stretch is applied to the ASM and parenchyma, both of which produce less tension than pre-stretched values post DI. In asthmatic subjects, ASM strain may be decreased potentially due to sti er airway walls or decreased stress transmitted to the ASM than in normals with a DI. When strain acts on activated ASM there is premature detachment of cross-bridges. As a result, Fredberg et al (25) suggested that during a DI, ASM strain can detach the slowly cycling latch bridges and increase the cycling of cross-bridges, which would lead to a decrease in the muscle force produced by the ASM. Therefore the greater the percentage of latch bridges, the sti er the airway smooth muscle will be and thus less strain will be associated with the DI (78). Additionally, periodic length uctuations mimicking tidal breathing are known to perturb the binding of myosin to actin, thereby limiting the tension that ASM can generate.
In vitro experiments performed by Mailhot-Larouche et al (57) investigated whether the duration between deep inspirations a ected the rate of ASM force recovery post deep inspiration. Their results reported that the shorter the time interval between successive DIs, the faster the rate of ASM force recovery post DI. The duration between successive deep inspirations on the rate of maximum ASM force recovery following DI in contracted guinea pig bronchi as performed by Mailhot-Larouche et al (57) is depicted in Figure 4.1. These DIs were simulated at 2, 5, 10 and 30 minutes duration between successive DIs. Their work also demonstrated that varying the time between successive DIs had no e ect on the measured ASM force or bronchodilatory e ect (ASM force before DI-ASM force after DI) of the DI. Following a period of equilibration, the ASM preparations were exposed to length oscillation under di ering contractile states; contraction triggered by methacholine (MCh) and relaxation triggered by isoproterenol (ISO). Each protocol consisted of a sequence of length oscillations at 0:2Hz for approximately 60 minutes. This simulated tidal breathing interspersed with DI. The ASM preparations were then subjected to a series of four sequences of tidal breathing of varying durations, each ending with a simulated DI. Tidal breathing lasted two minutes after the last DI.
1.1 Thesis structure
1.2 Airway smooth muscle physiology
1.3 Numerical methods for airway smooth muscle
1.4 Airway wall models
1.5 ASM and airway wall coupled system
2 The Distribution Moment Approximation
2.1 Zahalak’s approximation to the Huxley model
2.2 Extension of DM method to coupled 4-state ASM-airway wall model
3 Numerical solution of airway closure via airway smooth muscle interaction
3.1 Numerical solution of Mijailovich’s crossbridge model
3.2 Inter-model comparison
3.3 Results and discussion
4 Rate of ASM force recovery post deep inspiration.
4.1 Numerical simulation of deep inspiration
4.2 Results and discussion
5 Clustered dynamics of branched airways
5.2 Numerical simulation of branched airways
GET THE COMPLETE PROJECT
Airway smooth muscle and airway wall coupled dynamics: A numerical study using the distribution moment approach