Mathematical Model of Nucleotide Regulation on Airway Epithelia

In the airways, adenine nucleotides support a complex signaling network mediating host defenses. Released by the epithelium into the airway surface liquid (ASL) layer, they regulate mucus clearance through P2 (ATP) receptors, and following surface metabolism through P1 (adenosine; Ado) receptors. The complexity of ASL nucleotide regulation provides an ideal subject for biochemical network modeling. A mathematical model was developed to integrate nucleotide release, the ectoenzymes supporting the dephosphorylation of ATP into Ado, Ado deamination into inosine (Ino), and nucleoside uptake. The model also includes ecto-adenylate kinase activity and feed-forward inhibition of Ado production by ATP and ADP. The parameters were optimized by fitting the model to experimental data for the steady-state and transient concentration profiles generated by adding ATP to polarized primary cultures of human bronchial epithelial (HBE) cells. The model captures major aspects of ATP and Ado regulation, including their >4-fold increase in concentration induced by mechanical stress mimicking normal breathing. The model also confirmed the independence of steady-state nucleotide concentrations on the ASL volume, an important regulator of airway clearance. An interactive approach between simulations and assays revealed that feed-forward inhibition is mediated by selective inhibition of ecto-5′-nucleotidase. Importantly, the model identifies ecto-adenylate kinase as a key regulator of ASL ATP and proposes novel strategies for the treatment of airway diseases characterized by impaired nucleotide-mediated clearance. These new insights into the biochemical processes supporting ASL nucleotide regulation illustrate the potential of this mathematical model for fundamental and clinical research.

Mucociliary clearance (MCC) 4 constitutes the first line of defense against airway infection (1). Inhaled pathogens are trapped by a mucus layer, positioned above ciliated epithelia by a periciliary (PCL) layer. Together, mucus and PCL form the ASL layer, which is continuously transported cephalad by coordinated ciliary beating. Through their interactions with P2 receptors, adenine nucleotides regulate all major epithelial functions supporting MCC, including Cl Ϫ /liquid secretion via the Ca 2ϩ -activated Cl Ϫ channel and the cystic fibrosis transmembrane regulator (CFTR) (2), mucin secretion (3,4), and ciliary beat activity (5,6). Nucleotides are released by the epithelium into the ASL layer under basal conditions (7,8), and their release rate increases in response to mechanical stress imparted by tidal breathing (8 -12). Nucleotide-mediated epithelial responses are modified by surface enzymes (ectonucleotidases) that convert a fraction of the ATP into Ado (13)(14)(15). This nucleoside activates airway epithelial signaling pathways through P1 receptors, typically the A 2B receptor, which stimulates ciliary beat activity (16) and Cl Ϫ /liquid secretion by CFTR (17). In healthy airways, the purinergic regulation of MCC is mediated by both P1 and P2 receptor-mediated responses (10,18). In chronic lung diseases, the accelerated metabolism of ATP (14) can lead to severe complications related to airway obstruction and inflammation (2,19). The cystic fibrosis (CF) patients are particularly vulnerable to infection because of the sole dependence of MCC on P2 receptor-mediated signaling, reflecting the absence of the final effector for P1 signaling, i.e. functional CFTR. This information led to the current interest in aerosolized nucleotide analogs for the improvement of MCC (19).
Despite the importance of purinergic signaling for airway defenses, it has not been possible, until recently, to develop a mathematical model of ASL nucleotide regulation because of the absence of quantitative information regarding ASL nucleotide and Ado levels, the rates of nucleotide release and Ado uptake, and the kinetic properties of the ectonucleotidases. In the past decade, extensive work has been conducted to characterize ATP release from the apical surface of HBE cells (7, 8, 10 -12, 20, 21). Members of the three ecto-ATPase families (alkaline phosphatases, ecto-nucleoside triphosphate diphosphohydrolases (E-NTPDases), and nucleotide pyrophosphatase/phosphodiesterases (E-NPPs)) (22,23) have been identified on HBE cells (13-15, 24 -26), along with ecto-AK (2ADP 7 ATP ϩ AMP) (27,28). Also, Ado levels are regulated on HBEs by a balance between the production rate (AMP 3 Ado) by ecto-5Ј-NT and NSAP (13) and the elimination rate by Ado deaminase 1 (ADA1; Ado 3 Ino) and cellular uptake by concentrative nucleoside transporter 3 (CNT3) (29).
Based on these findings, we developed a mathematical model of ASL nucleotide regulation. The kinetic coefficients were optimized by fitting the output of the model to experimental values of ASL nucleotide and Ado concentrations, nucleotide release rates (12,20), ATP metabolism (13,26,28), and Ado uptake (29). The capacity of the model to reproduce the impact of physical parameters (ASL volume and mechanical stress) and functional interactions between network components was verified by biochemical assays on primary cultures of HBE cells. This study provides the first mathematical approach to the investigation of an important defense mechanism in human airways.

EXPERIMENTAL PROCEDURES
Cell Culture-Model validation was conducted on HBE cultures maintained at air-liquid interface, as described previously (30). The epithelial cells were obtained from healthy donors at the time of transplantation by the Tissue Culture Core of the Cystic Fibrosis Center (University of North Carolina) under the auspices of protocols approved by the Institutional Board of Review. They were isolated by protease digestion, seeded on 12-mm permeable supports (0.45 m pores, T-Clear; Costar) coated with human type IV collagen, and kept in Ham's F-12based medium supplemented with 5 g⅐ml Ϫ1 transferrin, 30 nM triiodothyronine, 1 M hydrocortisone, and 3.75 g⅐ml Ϫ1 endothelial cell growth substance. The confluent cultures were maintained at air-liquid interface with medium on the basolateral side. After 4 weeks, they were composed of ciliated and secretory cells, with transepithelial electrical resistance Ն300 ohms/cm 2 .
Cyclic Compressive Stress (CCS)-The HBE cultures were subjected to 1 h of CCS mimicking normal tidal breathing (20 cm of H 2 O; 20 cycles/min) using a system enclosed in an incubator (37°C; 5% CO 2 ) (10). Two pieces of tubing, each containing a high speed solenoid valve, were inserted into a conical silicon plug fitted to the top of the T-Clear, one connected to a pressurized chamber and the other connected to the incubator atmosphere. A microprocessor regulated the amplitude and frequency of the chamber pressure applied to the culture by controlling the open/closed state of the solenoid valves (10).
Biochemical Measurements of Steady-state Nucleotide Concentrations-The steady-state ASL ATP concentration was measured on HBE cultures by luminometry using soluble luciferase or a chimeric Staphylococcus aureus protein A-luciferase attached to the apical surface via keratan sulfate-specific antibodies (20). Steady-state concentrations of ATP, ADP, AMP, and Ado were determined on HBE cultures by ethenyl derivatization and fluorescence high pressure liquid chromatography (HPLC) (12).
Transient Nucleotide Concentration Profiles and Ectonucleotidase Activities-HBE cultures were preincubated 30 min (37°C; 5% CO 2 ) in Krebs buffer (pH 7.5) (0.35 ml; bilateral). Then 0.1-1000 M ATP was added to the apical surface to determine the transient concentration profile and initial rate of ATP hydrolysis by HPLC analysis of buffer samples collected over 60 min (14). The contribution of each ectoenzyme was calculated from paired assays conducted on the same culture wells before and after pretreatment with selective noncompetitive inhibitor(s). The kinetic parameters of NSAP for ATP and ADP were determined as described previously for AMP (13), but in the presence of 0.5 mM Ap 5 A to inhibit ecto-AK (ATP ϩ AMP 7 2ADP) (28). NSAP activity was calculated from the difference in linear rate of substrate hydrolysis measured before and after exposure to 10 mM levamisole, at each substrate concentration (0.1-1000 M). Michaelis constants (K M ) and maximal velocities (V max ) were calculated as shown previously (13).
Kinetic Model of Airway ATP/Ado Regulation-This section outlines the theory and methodology used to build the mathematical model of airway nucleotide regulation on the apical surface of human airway epithelia. A biochemical network was derived using quantitative data gathered from the literature on nucleotide release (12,20) and metabolism (13-15, 24, 26) by HBE cultures (Fig. 1). Because UTP levels are about 15% that of ATP (25) in the ASL layer, the influence of uridine nucleotides on ATP regulation was considered negligible. Therefore, the reactions governing ATP dephosphorylation to Ado and Ado deamination are as shown in Reactions 1-5, Each of the above reactions was assumed to follow Michaelis-Menten kinetics of the form shown in Equation 1, where [S] is the substrate concentration, and V max (j) and K M (j) are the maximum velocity and Michaelis constant. The superscript j is used to index the multiple enzymes involved in the pathway. In contrast, the reversible ecto-AK Reaction 6, ATP ϩ AMP 7 2ADP REACTION 6 which involves two binding sites, has been shown to follow a Bi Bi mechanism (32) as shown in Equation 2, This expression is a generalization of Michaelis-Menten kinetics that relies on the same underlying assumptions. The model also incorporates the feed-forward inhibition of AMP hydrolysis by ATP and ADP. Experimental evidence indicates that ATP and ADP competitively bind to enzymes that catalyze AMP dephosphorylation (13,33). Equation 3 was used to model this effect, Although it is possible that the inhibition constants K ATP in and K ADP in may differ for ecto-5Ј-NT and NSAP, these values were assumed identical for simplicity. The remaining elements that were incorporated into the model are the rates of ATP release and Ado/Ino uptake. In accordance with experimental observations (12), the basal rate of ATP release (J ATP ) was considered constant. Because Ado and Ino are taken back into the cell through CNT3, and therefore compete for access to the transporter, this effect is modeled using the following rates for Ado and Ino uptake as shown in Equations 4 and 5, which allows for the possibility of different V max and K M values for the two molecules.
One of the functions of nucleotide release is to stimulate liquid secretion, which in turn may affect ASL nucleotide concentrations. Under static conditions, the basal rate of release does not maintain ATP at sufficiently high levels to activate P2 receptors (12). However, Ado levels generated from ATP metabolism under static conditions stimulate liquid secretion into the ASL layer via the A 2B receptor-CFTR pathway (17). Therefore, the model equations were derived for the general condition of a time-dependent ASL volume. To illustrate how the ASL volume enters the model, consider a simple case in which ATP is released into the ASL at a constant rate of J ATP molecules/s, and then hydrolyzed by enzyme E. The surface area is denoted by A, and [E] T represents the surface density of E, expressed in mol⅐cm Ϫ2 . If V indicates ASL volume, the equation governing ATP concentration is shown in Equation 6, where (k cat A)/(V)[E] T corresponds to maximum velocity (V max ) and the last term models changes in ATP concentration because of the changing ASL volume. The considerations discussed above led to the formulation of five equations to describe the regulation of adenine nucleotides in the ASL layer (Equations 7-11). All model simulations were performed in MatLab (Mathworks, Natick, MA) using the ordinary differential equation solver of 15 s.

Linear Model for the Steady-state Nucleotide Concentrations-
The biochemical data showed that, under physiological conditions, steady-state nucleotide concentrations are orders of magnitude lower than the K M values of the enzymatic reactions (12). Hence, to a good approximation, the Michaelis-Menten reaction rates can be simplified (linearized) as shown in Equation 12, where the parameter P (k) denotes the sum of the ratios of V max (j) to K M (j) for all enzymes catalyzing the same reaction. The expressions for Ado and Ino uptake can be simplified in a similar way. The reversible ecto-AK reaction (Reaction 6) is a second order reaction requiring two substrates and therefore does not contribute significantly to the regulation of ASL nucleotide concentrations under steady-state conditions. Also, steady-state ATP and ADP levels are 1000-fold lower than their K i values for the feed-forward inhibition reaction (13) and therefore insufficient to inhibit AMP. Therefore, these equations were excluded from the linear model. Using these simplifications, Equations 7-11 can be approximated under steady-state conditions by the following set of linear Equations 13-17, J ATP Ϫ ͑P ͑1͒ ϩ P ͑5͒ ͓͒ATP͔ ss ϭ 0 (Eq. 13) P ͑1͒ ͓ATP͔ ss Ϫ P ͑2͒ ͓ADP͔ ss ϭ 0 (Eq. 14) P ͑5͒ ͓ATP͔ ss ϩ P ͑2͒ ͓ADP͔ ss Ϫ P ͑3͒ ͓AMP͔ ss ϭ 0 where U (1) and U (2) result from linearizing the expressions for Ado (Equation 4) and Ino (Equation 5) uptake, respectively.
Statistical Analysis-All biochemical assays were conducted on HBE cultures from Ն5 subjects (ϮS.E.), and the mean experimental value was considered nonsignificantly different from the corresponding simulated value when Ͻ15% different.

Generation of the Experimental Data Set for Model Parameter Optimization-
The strength of this model resides in the fact that all biochemical data were acquired from a single experimental system, polarized primary cultures of HBE cells maintained under air-liquid interface conditions and assayed in a physiological buffered solution (pH 7.5). Such consistency allowed us to combine data gathered over a decade into one experimental data set (Table 1, column A).
The first kinetic parameters of ATP and ADP hydrolysis measured on HBE cultures represent the sum of all apical ecto-ATPase and ecto-ADPase activities (26). More recently, we determined the kinetic properties of ecto- Table 1, column A) (28) and E-NPPs (V max (10) and K M (10) ; Table 1, column A) (24) on these cells. Based on this information, the kinetic measurements of ATP and ADP hydrolysis were repeated in the presence of the ecto-AK inhibitor, Ap 5 A (28), to eliminate the influence of the nonlinear reaction (ATP ϩ AMP 7 2ADP). Under these conditions, the ecto-ATPase and ecto-ADPase activities described by V max (1)(2)(3)(4)(5) (Table  1, column A) combined the activities of two enzyme families: E-NTPDases and alkaline phosphatases. We also recently Ϫ V max (7) [AMP] showed that HBE cells express NTPDase 1 and NTPDase 3 (15).
It is important to note that the K M corresponds to the substrate affinity of an enzyme, which is independent of cell type or expression level, whereas V max constitutes a quantitative measure of the expression level of an enzyme (33). Because specific inhibitors of NTPDase 1 are currently unavailable, NTPDase 1 and NTPDase 3 were assigned to the parameter sets with K M values approximating (Ϯ10%) their K M (Table 1,  We recently demonstrated that NSAP is the alkaline phosphatase isoform functionally expressed by human airway epithelia (13). Clinical assays, generally performed in basic solutions (pH 9 -10), depict NSAP as a low affinity enzyme (K M Ͼ400 M) (23). However, NSAP has been shown to exhibit both low (lowNSAP; K M Ͼ400 M) and high (highNSAP; K M ϭ  (13). Because the kinetic parameters of airway NSAP for ATP and ADP were unavailable, we conducted these experiments, as described previously (13) but in the presence of the ecto-AK inhibitor, Ap 5 A (28). Two kinetic parameter sets were identified for ATP and ADP hydrolysis (Table 2), with K M values within the ranges reported for lowNSAP and highNSAP activities (36 -38). These values also approximated the K M values of the high affinity and low affinity ecto-ATPase (K M (1) and K M (3) ) and ecto-ADPase (K M (4) and K M (5) ) activities (Table 1, column A). Because NSAP, NTPDase 1, and NTPDase 3 all support the dephosphorylation of ATP and ADP, enzyme activities with similar K M values were combined and designated as the high affinity (NTPDase 1/highNSAP) and the low affinity (NTPDase 3/lowNSAP) enzymes (Table 1, column A).
The list of experimental data were completed by adding the kinetic parameters of ADA1 (29), the cellular uptake of Ado by CNT3 (29), the basal rate of ATP release (12), and the feedforward inhibition of Ado production by ATP and ADP (13). Altogether, these assays and studies provided the experimental data set that was used to develop the mathematical model of nucleotide regulation on human airway epithelia (Table 1, column A).
Model Validation with Experimental Data-The model was tested by comparing simulated and previously reported experimental measurements of transient (26) and steady-state (12) nucleotide and nucleoside concentrations on the apical surface of HBE cells. We started our investigations using the multicomponent model ( Fig. 1) described by Equations 7-11. Experimental values were available for all parameters from earlier studies (Table 1), except the K M of ecto-AK for AMP and the K M and V max values of CNT3 for Ino. When used as free parameters to fit the data, considerable discrepancies were found between the simulated and experimental transient profiles of ATP hydrolysis ( Fig. 2A) and steady-state concentrations (Fig.  2B). The performance of the model was considerably improved by using the experimental data set as initial guesses in the nonlinear least square fitting routine in MatLab (Mathworks, Natick, MA). Using this approach, the model closely reproduced the transient and steady-state profiles (Fig. 2, C and D). Among the 30 estimated parameters, only 6 were Ͼ2-fold different from the experimental values (Table 1, column B) as follows: high affinity ecto-ATPase (V max (1) ), E-NPPs (V max (10) and K M (10) ), ecto-AK (V max (f) ) and the feed-forward inhibition constants (K ATP (in) and K ADP (in) ). However, a 12-fold discrepancy persisted between measured and predicted ATP release rates. Similar experimental values were reported by different research groups using various techniques (10,12,20), ruling out experimental error as the source of discrepancy. On the other hand, if the rate of ATP release was held fixed at the experimental value and the other parameters free to be determined by data fitting,

TABLE 2 Kinetic parameters of ATP and ADP dephosphorylation by NSAP on primary cultures of HBE cells
Michaelis constants (K M ; M) and maximal velocities (V max : nmol ⅐ min Ϫ1 ⅐ ml Ϫ1 ) were calculated from Woolf-Augustinson Hoftsee plots.

Reactions
Enzyme identity Parameters Experimental data the model reproduced the transient behavior of the system (Fig.  2E) but not the steady-state profile (Fig. 2F). Steady-state Analysis with a Linear Model-In an attempt to resolve the differences between experimental and predicted steady-state ATP concentrations, we developed a linear model to describe the system near steady state, as described by Equations 13-17. The model was analyzed using two strategies. First, the experimental steady-state concentrations were used as inputs, and the model was solved for the parameter values. Because the steady-state Ino level cannot be measured by fluorescent HPLC, it was considered an unknown value in this procedure. In this case, solving the linear equations generated a negative Ino concentration, which is biologically impossible. Second, the experimental data set (Table 1, column A) was used as input parameters, and the model was solved for the steadystate concentration profile. This approach generated a good approximation of the experimental ATP concentration, but it resulted in estimated ADP, AMP, and Ado levels Ͼ10-fold lower than the experimental values (Table 3). Altogether, these results suggested that the model was incompatible with the experimental data and was therefore missing key components.
Our current knowledge on nucleotide release is limited to ATP because the inhibitor mixture designed to prevent nucleotide metabolism (20) interferes with fluorescence HPLC. However, emerging studies suggest that exocytosis mediates the release of ATP, ADP, and AMP from mammalian cells (39 -41), including polarized epithelia (4,7,42). We recently demonstrated that human airway epithelia co-release ATP and mucin by Ca 2ϩ -regulated exocytosis (42). During their maturation, secretory vesicles accumulate ATP-and UDP-sugars in concentrations Ն50-fold above cytosolic levels, a fraction of which is converted to ADP and UDP, and eventually to AMP and UMP (40). As a result, vesicle fusion to the apical membrane is expected to release ATP, ADP, and AMP into the ASL layer. These observations motivated the addition of terms describing ADP and AMP release to the mathematical model. Constant fluxes for ADP (J ADP ) and AMP (J ADP ) were added to the linear equations (Equations 13-17), and the values were adjusted so the model would generate steady-state concentrations closest to the experimental values. Table 3 shows that adding ADP and AMP release to the model eliminated the discrepancies between the experimental and predicted steadystate concentration profiles.
Reanalysis of the Nonlinear Model-The multienzyme model with the nonlinear equations (Equations 7-11) was updated to include terms describing the release of ADP (J ADP /V in Equation 8) and AMP (J ATP /V in Equation 9). The model was then refitted using the experimental data set (Table 1, column A) and the parameters of J ADP and J AMP ( Table 1, column C). Under these conditions, the simulations closely reproduced the exper-imental transient (Fig. 3A) and steady-state (Fig. 3B) concentration profiles, and eliminated the 12-fold difference between measured and estimated J ATP . Finally, the predicted V max and K M values ( Impact of ASL Volume on Nucleotide Regulation-The model identified an important property of airway nucleotide regulation regarding ASL volume, which has been identified as a key determinant of mucus clearance (1,2). Because the rate at which ATP release raises ATP concentration (J ATP ) and the maximum reaction velocity (V max ϭ (k cat A)/(V)[E] T ), both scale inversely with ASL volume in Equation 6, the model predicts that steady-state concentrations are independent of ASL volume. Following a simulated 2-fold increase in ASL volume, which acutely dilutes ATP and ADP, nucleotide concentrations  are predicted to rapidly return to steady-state values (Fig. 4A). These simulation data were validated experimentally on HBE cultures using two luminometry assays, soluble luciferase and S. aureus protein A-luciferase attached to the epithelial surface. Fig. 4B shows that steady-state ATP and ADP levels, measured 60 min after the addition of 10 -500 l of phosphate-buffered saline, were unaffected by ASL volume changes. This result provides new insights into the biochemical processes supporting MCC. For instance, the airways of CF patients exhibit an ϳ70% reduction in ASL volume, which plays a critical role in airway obstruction (1,2). Under steady-state conditions, the model predicts that ASL dehydration in CF does not concentrate, i.e. increase the availability of ATP and Ado for P2 and P1 receptor-mediated mucus clearance as a compensatory mechanism to restore ASL volume.
Model Adaptation to Mechanical Stress-Normal tidal breathing applies compressive and shear stresses on the airway walls that stimulate epithelial ATP release and ATP-mediated MCC (8,10). This physiological response is particularly relevant to airway homeostasis as CCS mimicking tidal breathing was shown to restore ASL volume and mucus transport on the dehydrated airways of CF patients (10). Because the mathematical model was developed using experimental data acquired on static cultures, we tested its capacity to predict the ASL nucleotide and nucleoside concentrations maintained by mechanical stress mimicking tidal breathing. We recently demonstrated that CCS increases the rate of ATP release by 9.3-fold on HBE cultures (10). In this study, we conducted experiments to test the impact of CCS on the entire steady-state nucleotide concentration profile. Fig. 4C shows that CCS raised ATP, ADP, AMP, and Ado concentrations by 4 -8-fold within the ASL layer. Model simulations were conducted to determine whether this response to mechanical stress was best fitted by an increase in the rate of ATP release or by concomitant increases in the release of several nucleotides. The equations solved with a 9.3fold higher rate of ATP release (J ATP ; Table 1, column C) reproduced the effect of CCS on ATP but not on ADP or AMP (Fig.  4C). This simulation suggested that ATP release and metabolism were unable to account for the impact of CCS on ADP and AMP levels. Simultaneous 9.3-fold increases in J ATP and J ADP generated ATP and ADP levels within experimental range, whereas AMP levels remained Ͼ2-fold lower. The impact of CCS on the ASL concentration profile was best fitted by a model simulation in which J ATP , J ADP , and J AMP were all 9.3-fold higher than base line. Altogether, these data are consistent with the general proposal that mechanical stress induces the release of all three nucleotides through exocytosis in polarized epithelia (39 -42). Furthermore, they demonstrate the ability of the mathematical model to predict the complex behavior of ASL nucleotide regulation under static and dynamic conditions. Novel Approach to Biochemical Network Analysis-Although several ectonucleotidases have been characterized on the apical surface of HBE cells (13-15, 24 -29), the key enzyme(s) regulating nucleotide-mediated MCC remain unidentified because of the complexity of the multienzyme network. The mathematical model allowed us, for the first time, to determine the contribution of each enzyme and enzyme group to the regulation of both physiological and pathological ASL nucleotide concentra- tions. Several airway ectonucleotidases have been shown to hydrolyze ATP, ADP, and/or AMP as follows: NTPDase 1, NTPDase 3, NSAP, E-NPPs, and ecto-AK (14,15,28). Based on their kinetic properties, we hypothesized that specific enzymes regulate ASL nucleotide levels encountered under physiological conditions (Ͻ0.1 M) (8,12,20), whereas others regulate higher ASL concentrations reached during tissue injury (0.1-1 M) (43)(44)(45)(46) or aerosolization (1-100 M) for diagnostic (47,48) or therapeutic (49 -52) purposes. The biochemical network was first simplified by classifying the enzymes into the following two groups: 1) low affinity high capacity (K M Ն 50 M; NTPDase 3, lowNSAP), and 2) high affinity low capacity (K M Ͻ 50 M; NTPDase 1, E-NPPs, highNSAP, and ecto-AK). Model simulations were conducted by blocking the V max of each enzyme ( Table 1, column C) to reproduce the effects of selective inhibitors on HBE cultures. One of them, levamisole, inhibits both the high and low affinity activities of NSAP (13), which could introduce validation errors into group comparisons. This possibility was addressed by comparing the impact of eliminating highNSAP, lowNSAP, or high/low NSAP activities on ASL nucleotide regulation within the two concentration ranges of interest: aerosolized/pathological (10 -100 M) and physiological (0.1-1 M). For V max values corresponding to the total maximum activities of two enzymes, NTPDase 1/highNSAP (V max (1) and V max (4) ) and NTPDase 3/lowNSAP (V max (5) ), NSAP activity was blocked by reducing the estimated parameter (Table 1, column C) proportionally by its experimental contribution ( Table 2) to total enzyme activity ( Table 1, column A). Fig. 5A shows that blocking highNSAP reduced the elimination rate of 100 M ATP by Ͻ10%. In contrast, blocking lowNSAP or high/ lowNSAP reproduced the effects of 10 mM levamisole on the HBE cultures and increased the half-life (t1 ⁄ 2 ) of 100 M ATP by 1.6-fold. Conversely, the metabolism of 1 M ATP was not affected by blocking lowNSAP, whereas blocking highNSAP or high/lowNSAP mimicked the effects of levamisole with a 1.9fold increase in the t1 ⁄ 2 of ATP (Fig. 5B). These data indicate that lowNSAP and highNSAP activities exhibit negligible overlap, regulating excess and physiological ATP concentrations selectively. Consequently, model simulations in which both NSAP activities are blocked are expected to reproduce the effects of levamisole within the concentration range of interest.
We first acquired general information on the biochemical network by comparing the relative contribution of the low affinity and high affinity enzyme groups to the regulation of aerosolized/pathological nucleotide concentrations. Model simulations in which the kinetic parameters of the low affinity enzymes (high/lowNSAP and NTPDase 3) were blocked predicted a 4-fold decrease in elimination rate of 100 M ATP and significant delays in the accumulation of ADP, AMP, and Ado (Fig. 5C). Comparable profiles were obtained experimentally using HBE cultures assayed with 100 M ATP in the absence (control) or presence of an inhibitor mixture specific for the low affinity enzymes (Fig. 5C) as follows: 10 mM levamisole (NSAP) (13) ϩ 10 M NF-279 (NTPDase 3) (53). In contrast, model simulation in which all the high affinity enzymes (high/low-NSAP, NTPDase 1, E-NPPs, and ecto-AK) were blocked yielded minor modifications in the transient profile, with a 1.5fold decrease in elimination rate of 100 M ATP (Fig. 5D).
Although experimental validation of the high affinity block awaits the identification of noncompetitive inhibitors for NTPDase 1 and E-NPPs, these data suggest that high ATP levels (Ͼ10 M) are regulated mainly by the low affinity enzymes on airway epithelia.
Comparative analyses conducted over a broad range of ATP concentrations predicted a shift in the relative contribution of the high and low affinity enzyme groups. To carefully scrutinize the involvement of each group, the two NSAP activities were analyzed independently. Model simulations in which the high affinity enzymes (highNSAP, NTPDase 1, E-NPPs, and ecto-AK) or the low affinity enzymes (lowNSAP and NTPDase 3) were blocked raised the t1 ⁄ 2 of 100 M ATP by 1.2-and 2.3-fold, respectively (Fig. 5E). Although the t1 ⁄ 2 of 10 M ATP would be equally affected by either group, the high affinity group is predicted to dominate at lower ATP concentrations. These simulations were reproduced experimentally on HBE cultures with a . E, simulated half-life (t1 ⁄2 ) of 0.1-100 M ATP generated by blocking the high affinity (solid black box) or low affinity (striped box) group, and experimental data (solid red box) for the low affinity block as above (n ϭ 5; ϮS.E.). F, simulated impact of blocking the high affinity (solid black box) or low affinity (striped box) group on the steady-state profile. Validation (solid red box) was by comparing ASL nucleotide levels before/after incubating HBE cultures with the low affinity inhibitor mixture (n ϭ 5; ϮS.E.). low affinity inhibitor mixture (10 mM levamisole ϩ 10 M NF-279) (Fig. 5E). The model suggests that the higher contribution of the low affinity group measured experimentally for 0.1 M ATP reflects the inhibitory effect of levamisole on highNSAP. Together, simulated and experimental results support a major role for high affinity enzymes in the regulation of ASL ATP within the effective concentration range of P2Y 2 receptors (0.01-1 M) (16,17).
Recent studies showed that ASL nucleotide metabolism is accelerated in chronic lung diseases, including CF, chronic obstructive pulmonary disease, and primary ciliary dyskinesia (14), which could compromise nucleotide-mediated MCC (19,52). Consequently, activation of P2Y 2 receptors by aerosolized nucleotides is being tested as a therapy for CF (52,54). As an alternative strategy to raise steady-state ASL ATP concentration, we conducted simulations to test the potential of an aerosolized inhibitor mixture against specific enzymes. As expected from their kinetic properties, model simulations in which the low affinity enzymes (lowNSAP and NTPDase 3) were blocked predicted a Ͻ2-fold increase in steady-state nucleotide concentrations (Fig. 5F). This simulation was reproduced experimentally by assays conducted on HBE cultures exposed to the low affinity inhibitor mixture (10 mM levamisole ϩ 10 M NF-279) (Fig. 5F). Because levamisole does not discriminate between the lowNSAP and highNSAP activities, these assays demonstrate that NSAP does not contribute significantly to the regulation of steady-state ATP concentrations (i.e. Ͻ0.1 M) on airway epithelial surfaces. To our surprise, the model predicts that the steady-state nucleotide levels are sensitive to the high affinity enzyme group with respect to ADP but not ATP (Fig. 5F), which suggested the influence of a nonlinear component within the biochemical network.
The properties of the reversible ecto-AK reaction suggest that the net directionality depends on the relative substrate concentrations (Fig. 6A). Micromolar ATP, encountered during injury or aerosolization, would drive the forward reaction, whereas steady-state conditions would favor the reverse reaction because of the relatively high ADP concentration (Table 3). Therefore, simulations and assays were designed to test whether ecto-AK influences the ATP/ADP concentration ratio on airway epithelial surfaces in these two situations (Fig. 6A). Model simulations predicted that blocking ecto-AK would not affect the metabolism of 100 M ATP (Fig. 6B) or the steadystate concentration profile (Fig. 6E). These simulated data were reproduced experimentally for the transient (Fig. 6B) and steady-state (Fig. 6E) concentration profiles by HBE cultures exposed to the vehicle control or buffer containing the ecto-AK inhibitor, 0.5 mM Ap 5 A (28). We next tested the hypothesis that ecto-AK activity may be significant within a functional unit addressing a specific nucleotide concentration range. Fig. 6C shows that blocking ecto-AK is not expected to influence the metabolic profile of 100 M ATP regulated by the low affinity enzymes, i.e. determined in the absence of the other high affinity enzymes. In contrast, model simulations and experimental data suggest that ecto-AK contributes to the regulation of ASL nucleotides under conditions where the high affinity enzyme activities predominate. First, the model predicts that blocking ecto-AK would not affect the elimination rate of 100 M ATP by a high affinity biochemical network, i.e. determined in the absence of the low affinity enzymes (Fig. 6D). However, when ADP concentration rises above that of ATP (ϳ20 min), the absence of ecto-AK would be associated with higher rates of ADP elimination and Ado accumulation (Fig. 6D). Similar results were obtained experimentally with HBE cultures using the low affinity inhibitor mixture (10 mM levamisole and 10 M NF-279) with/without 0.5 mM Ap 5 A (Fig. 6D). These data were also consistent with the behavior of ecto-AK reported on HBE cultures exposed to an ATP concentration (1 M) (28) predicted to be regulated predominantly by high affinity enzymes (Fig. 5E). The reverse ecto-AK reaction would compete with the high affinity ADP-hydrolyzing enzymes to trap adenine nucleotides into a trans-phosphorylation cycle (28). These findings suggest that airway ecto-AK regulates ASL nucleotide levels within a high affinity biochemical network exposed to a high ADP/ATP concentration ratio.
Based on this information, we revisited the impact of ecto-AK on ASL nucleotide regulation at steady-state conditions, which is expected to favor the reverse ecto-AK reaction (2ADP 3 ATP ϩ AMP) (Fig. 6A) because of the high ADP/ATP ratio (ϳ20; see Table 3). A model simulation with a high affinity inhibitor mixture, excluding the ecto-AK inhibitor, resulted in a 6-fold increase in steady-state ATP concentration (Fig. 6E), compared with 2-fold for a complete high affinity block (Fig.  5F). The physiological relevance of this information is illustrated by the fact that a comparable increase in ASL ATP level was generated by CCS mimicking normal tidal breathing (Fig.  4C), which restored normal ASL height and mucus transport on the dehydrated airway epithelial surfaces of CF patients (10). The mathematical model of ASL nucleotide regulation therefore identified the key enzymes regulating ASL ATP and proposes a novel approach for the treatment of airway obstruction in CF, an aerosolized inhibitor mixture against high affinity NTPDase 1 and E-NPPs.
Adenosine Regulation in Airway Defenses-The model was instrumental to our investigation of airway Ado regulation, supplementing the biochemical assays with simulated information that would not be available experimentally. First, we modeled the relative contribution of the two ectonucleotidases supporting the conversion of AMP to Ado on airway epithelial surfaces, NSAP and ecto-5Ј-NT (13). Model simulations were run by blocking the parameters of the following: 1) ecto-5Ј-NT; 2) high/lowNSAP toward AMP [NSAP (AMP)], or 3) high/ lowNSAP toward ATP, ADP, and AMP [NSAP (all)]. These simulations were verified by biochemical assays conducted on HBE cultures before and after the inhibition of ecto-5Ј-NT (10 mM concanavalin A) (13,55) or NSAP (10 mM levamisole) (13,56). Although blocking all parameters of NSAP would mimic the effects of levamisole, we tested whether the concomitant reductions in ecto-ATPase and ecto-ADPase activities would influence the feed-forward inhibition of Ado accumulation (Equation 3) (13). Hence, all simulations and assays were initiated with ATP to probe the response of the entire biochemical network.
Comparative analysis conducted over a range of ATP concentrations revealed a shift in the relative contribution of ecto-5Ј-NT and NSAP (Fig. 7, A-C). Both simulated and experimental data showed that ecto-5Ј-NT dominates at physiological ATP levels (0.1-1.0 M), as indicated by larger reductions in the accumulation rate of Ado (Fig. 7, A and B). In contrast, Ado production from 10 M ATP would be equally reduced by blocking ecto-5Ј-NT or NSAP (AMP) (Fig. 7C). Interestingly, an intermediate response was predicted by blocking NSAP (all) (Fig. 7C) or measured experimentally with levamisole (Fig. 7C). These data suggest that levamisole increases the availability of ATP and ADP for the feed-forward inhibition of Ado production by the remaining enzyme, ecto-5Ј-NT.
Model simulations were designed to evaluate the physiological relevance of feed-forward inhibition for airway homeostasis by calculating the time required for 1-1000 M ATP to generate an Ado concentration equivalent to the EC 50 of A 2B receptors (0.3 M) (16,17). Both simulated and biochemical data demonstrated that the delay to reach 0.3 M Ado is influenced by the ASL ATP concentration in a biphasic manner (Fig.  7D). A relatively long delay would be expected for concentrations (i.e. 1 M ATP) below the catalytic capacity of most airway ectonucleotidases (Table 1, column A). In the 10 -100 M ATP range, both enzymes efficiently converted AMP to Ado, whereas longer delays are required for ATP levels challenging their V max .
The mathematical model was able to identify the enzyme targeted by feed-forward inhibition using simulations that cannot be reproduced experimentally: by blocking the kinetic constants (K ATP (in) and K ADP (in) ; Table 1, column C). We compared the impact of blocking each AMP-hydrolyzing enzyme on the ATP-mediated regulation of Ado production in the presence or absence of feed-forward inhibition. Model simulations in which the V max of ecto-5Ј-NT activity was blocked resulted in slightly longer Ado delays (Ͻ 30%) at all ATP concentrations, which were partially reduced to control values by blocking the inhibition constants (Fig. 7D). These simulations predict a weak control of feed-forward inhibition over the remaining enzyme, NSAP, which may be ascribed to its capacity to eliminate ATP and ADP in close vicinity to the catalytic site.
The reciprocal simulation exposed a strong relationship between ASL ATP concentration and Ado production by ecto- 5Ј-NT. The time required for Ado to reach the EC 50 of A 2B receptors (0.3 M) (16,17) increased by up to 6-fold when ecto-5Ј-NT dominated, i.e. in the absence of NSAP, with a regulatory threshold of 1 M ATP (Fig. 7E). These simulations were reproduced by HBE cultures pretreated with the NSAP inhibitor (10 mM levamisole). The model predicts that this relationship would be completely abrogated by simultaneously blocking NSAP (all) and the feed-forward inhibition (Fig. 7E). The mathematical model therefore provided new insights into the biochemical processes regulating ASL Ado, which would not have been available otherwise. The simulations predict that feedforward inhibition targets specifically Ado production by ecto-5Ј-NT on HBE surfaces. However, whereas the sensitivity of ecto-5Ј-NT activity to Ն1 M ATP is consistent with the properties of the vascular model of nucleotide regulation (57)(58)(59), airway epithelia seldom experience such high ATP concentrations, except locally during injury or aerosolization.
Finally, we determined the contribution of NSAP and ecto-5Ј-NT to the regulation of the steady-state nucleotide concentrations. The model predicts that NSAP inhibition would raise ATP, ADP, and AMP concentrations by 1.4-, 2.5-, and 2.8-fold, respectively (Fig. 7F). This inverse relationship between NSAP activity and the degree of phosphorylation of the substrate was closely reproduced experimentally by levamisole on HBE cultures (Fig. 7F) and is consistent with the properties of the purified enzyme: ATP Ͻ ADP Ͻ AMP (60,61). In contrast, model simulations in which ecto-5Ј-NT activity was eliminated only increased AMP concentration, in agreement with the substrate specificity of the enzyme (56). These results were reproduced by HBE cultures exposed to the selective inhibitor, concanavalin A (Fig. 7F). Altogether, these data demonstrate that the model captures the contributions of ecto-5Ј-NT and NSAP to the regulation of steady-state nucleotide concentrations. On the other hand, inhibiting ecto-5Ј-NT or NSAP did not influence the steady-state concentration of Ado, which suggests that one AMP-hydrolyzing enzyme may be sufficient under conditions where AMP concentration is Ͼ100-fold lower than their K M values.

DISCUSSION
This study presents the first mathematical model of extracellular nucleotide and nucleoside regulation for the respiratory system. Based on experimental data reported in the literature (12,13,20,26,28,29), and acquired in this study, for polarized HBE cultures, Equations 7-11 were designed to incorporate all known factors regulating ASL nucleotides and Ado as follows: epithelial ATP release, the reactions mediating the surface conversion of ATP into ADP, AMP, and Ado, the surface elimination of Ado by deamination and cellular uptake, the trans-phosphorylating ecto-AK activity, and the feed-forward inhibition of Ado production by ATP and ADP. The simulations closely reproduced the transient and steady-state concentration profiles generated on the apical surface following the addition of ATP. However, regardless of the fitting strategy, a 12-fold difference persisted between the simulated and experimental rates of ATP release, suggesting the model was missing key elements.
The current information on the mechanisms of nucleotide release is derived from measurements of ATP release (62,63).
Yet recent breakthroughs in our understanding of exocytosis support the notion that nucleotides, in addition to ATP, are also released via this route by a variety of cell types (4, 39 -41), including HBE cells (42). During their maturation, the vesicles accumulate ATP, ADP, and AMP, which are released during exocytosis (39 -41). Importantly, the addition of terms describing ADP and AMP release provided the best fit of the experimental data, as well as the transient and steady-state concentration profiles. Furthermore, it resolved the discrepancy between simulated and experimental rates of ATP release.
The resulting model was validated by testing the impact of mechanical stress on ASL nucleotide regulation. Airway epithelia, in vivo, are subjected to mechanical forces imparted by tidal breathing, which are known to stimulate ATP release (8,10). We recently showed that CCS mimicking normal tidal breathing causes a 9.3-fold increase in the rate of ATP release from HBE cultures (10). In this study, we showed that CCS raises the ASL concentrations of ATP, ADP, AMP, and Ado by 4 -7-fold. Model simulations in which J ATP was increased 9.3-fold reproduced the effects of CCS on ATP but not on the other nucleotides. The model therefore predicted that mechanical stressinduced ATP release and surface metabolism are unable to account for the effects of CCS on the ASL nucleotide composition. In fact, the best fit was obtained by raising simultaneously J ATP , J ADP , and J AMP by 9.3-fold. These results support the importance of including parameters for the release of all three nucleotides into the mathematical model of ASL nucleotide regulation. Furthermore, they illustrate the capacity of the model to predict the impact of physiological phenomena such as breathing on ASL nucleotide regulation.
Multiple ATP-hydrolyzing enzymes compete for the same pool of nucleotides within the ASL layer. The development of a mathematical model incorporating all the ectonucleotidases provided, for the first time, an invaluable tool to assess their individual contribution within the multienzyme system. The complex biochemical network was simplified by sorting the enzymes into kinetically distinct groups as follows: low affinity high capacity (K M Ͼ 50 M; lowNSAP and NTPDase 3) and high affinity low capacity (K M Ͻ 50 M; highNSAP, E-NTPDase 1, E-NPPs and ecto-AK) (26,64). Model simulations accurately predicted the predominant role of the low affinity group in the regulation of high ATP concentrations (Ͼ1 M) (43)(44)(45)(46), as demonstrated on HBE cultures using selective inhibitors of NSAP (10 mM levamisole) (13) and NTPDase 3 (10 M NF-279) (53). Therapies inhibiting these ectoenzymes could potentially accentuate airway epithelial apoptosis (65,66) and inflammatory responses (45) because of increased ATP-mediated activation of P2X 7 receptors, which are up-regulated in lung diseases (67). In addition, these enzymes are likely responsible for the rapid elimination of aerosolized nucleotides (1-100 M) used to diagnose asthma and chronic obstructive pulmonary disease (ATP and AMP) (47,48), for the treatment of CF and primary ciliary dyskinesia (UTP) (49 -52) and in animal studies (48,68,69).
The model predicts that high affinity enzymes play a major role in the regulation of ATP within the physiological concentration range encountered on airway epithelial surfaces (8,12,20). The recent report that CCS-mediated increase in ASL ATP level restores ASL volume and mucus transport on dehydrated HBE cultures from CF patients (10) motivated the search for an inhibitor mixture against the high affinity enzymes to enhance P2Y 2 receptor-mediated mucus clearance. However, the model predicts that blocking all four enzymes would raise ADP concentrations, instead of ATP. An explanation was provided by separating the high affinity enzymes into linear and nonlinear components. Inhibition of the enzymes supporting the linear reactions (ATP 3 ADP 3 AMP; NTPDase 1, E-NPPs and highNSAP) is expected to raise the concentrations of all the released nucleotides. In contrast, the directionality of the reversible ecto-AK reaction (ATP ϩ AMP 7 2ADP) would be influenced by the ASL nucleotide composition. Model simulations and biochemical assays suggest that ecto-AK does not influence the elimination of high ATP concentrations (Ͼ10 M). However, blocking ecto-AK significantly modified the transient metabolic profile generated by the high affinity enzyme group once ADP level exceeded that of ATP by Ͼ3-fold, as depicted by enhanced ADP, but increased Ado, accumulation. This behavior is consistent with the substrate requirements of the reversible ecto-AK reaction. The forward reaction (ATP ϩ AMP 3 2ADP) is limited on airway epithelia by the high ATP/AMP ratio maintained through efficient conversion of AMP to Ado. In contrast, the reverse reaction (2ADP 3 ATP ϩ AMP) is favored once ADP concentration exceeds that of [ATP ϩ AMP]. However, the ATP generated by ecto-AK would be rapidly dephosphorylated by the high affinity ecto-ATPases (NTPDase 1, E-NPPs, and highNSAP). Consequently, these data suggest that removing the strong competitive effect of the linear high affinity enzymes may allow ecto-AK to raise ATP concentrations within a nucleotide mixture dominated by ADP.
Based on this analysis of the high affinity group, we evaluated the role of ecto-AK under steady-state conditions, in which ASL exhibits an ADP/ATP concentration ratio of ϳ20 (Table  3). The model predicts that a high affinity inhibitor mixture, excluding the ecto-AK inhibitor (Ap 5 A), would raise the ASL ATP concentration by ϳ6-fold. The therapeutic potential of this approach is illustrated by the fact that CCS mimicking normal breathing produced a similar increase in ASL ATP level (Fig. 4C), which restored mucus transport on the airway epithelial surfaces of CF patients (10). The physiological importance of ecto-AK for the regulation of steady-state nucleotide levels was also demonstrated in human hepatocytes, where it regulates P2Y 13 (ADP) receptor-mediated high density lipoprotein endocytosis (70). The ecto-AK inhibitor stimulated high density lipoprotein endocytosis in the absence of exogenous substrates, supporting a significant role in cholesterol transport under steady-state conditions. Altogether, the model and experimental data suggest that aerosolization of an inhibitor mixture against the high affinity enzymes (excluding Ap 5 A) may stimulate airway clearance in CF patients. The fact that blocking NSAP did not affect the steady-state ATP level suggests that levamisole may be omitted. Further refinement of the high affinity enzyme inhibitor mixture for the treatment of obstructive lung diseases awaits the discovery of specific noncompetitive inhibitors for NTPDase 1 and the E-NPPs.
The regulation of ASL nucleotides through cell surface metabolism leads to the accumulation of another important signaling molecule for airway defenses, i.e. Ado. We previously demonstrated that two ectoenzymes contribute to the production of Ado on airway epithelial surfaces, ecto-5Ј-NT and NSAP (13). Although their characterization was conducted using AMP as initial substrate (13), this study documented the impact of the entire biochemical network on ASL Ado production using ATP as initial substrate. Based on model simulations and biochemical assays, we identified a negative relationship between ATP Ն1 M and the time required for Ado concentration to reach the EC 50 of A 2B receptors for CFTR activation (0.03 M) (16,17). Model simulations, in which the kinetic constants K ATP (in) and K ADP (in) were blocked, suggest that this relationship is caused by feed-forward inhibition of ATP on ecto-5Ј-NT, but not on NSAP. These data are in agreement with the reported inhibition of purified ecto-5Ј-NT by ATP and ADP (71,72). Although the inhibitory mechanism has not been identified, the insensitivity of NSAP may be partially because of its capacity to eliminate ATP and ADP. Altogether, these results emphasize the critical importance of studying the properties of ectonucleotidases using experimental designs accounting for the entire biochemical network, as the suppression of ecto-5Ј-NT activity by ATP would have been missed by protocols isolating the AMP 3 Ado reaction.
In human airways, feed-forward inhibition is not expected to play a significant role in the regulation of purinergic signaling, as epithelial surfaces are not normally exposed to Ͼ1 M ATP (73)(74)(75). However, this regulatory mechanism may become relevant during injury-induced ATP release to prevent A 2B receptor desensitization (76) and/or excessive inflammatory responses caused by chronically high Ado concentrations (77). But even under these circumstances, NSAP activity would be expected to compensate for the reduction in ecto-5Ј-NT activity caused by excess ATP. Consequently, this regulatory mechanism may be more relevant for purinergic signaling in tissues expressing ecto-5Ј-NT, with little or no NSAP. For instance, feed-forward inhibition is not detected on vascular smooth muscle cells (78), which express NSAP (79). In contrast, endothelial surfaces do not express NSAP (80) and exhibit an inhibitory effect of ADP on ecto-5Ј-NT activity (59,81). Because ADP activates platelets and Ado suppresses their activation (31,82), feed-forward inhibition was proposed to regulate thrombus formation during endothelial repair. A mathematical model of extracellular nucleotide metabolism on endothelial cells predicts an inhibition of Ado production by 1-10 M ADP (57)(58)(59). Because activated platelets release Ͼ1 mM ADP on endothelial surfaces (31), feed-forward inhibition would create a time gap between ADP-and Ado-mediated responses, allowing sufficient time to consolidate the thrombus while minimizing spread beyond the site of damage. Consequently, two factors determine the importance of feed-forward inhibition for purinergic signaling, the ecto-5Ј-NT/NSAP activity ratio and local nucleotide concentrations.
In conclusion, this study provides a robust model for the regulation of two essential signaling molecules for airway defenses, ATP and Ado. The equations capture the complex behavior of the linear and nonlinear reactions and the dynamic interplay between the ectoenzymes. Furthermore, the model reproduces ASL nucleotide regulation under both static conditions and mechanical stress mimicking normal breathing. Finally, therapeutic applications are proposed based on an interactive approach between simulations and biochemical assays. We are currently developing a model for the regulation of the ASL volume by P2Y 2 and A 2B receptor-mediated ion channel activities. In the near future, the integration of these two models should provide insightful information on the interactions between the biochemical and functional networks of purinergic signaling on airway epithelia under normal and pathological conditions, as well as novel guidelines for the development of efficient treatments of chronic obstructive lung diseases.