Mathematical model of serine protease inhibition in the tissue factor pathway to thrombin.

A mathematical model has been developed to simulate the generation of thrombin by the tissue factor pathway. The model gives reasonable predictions of published experimental results without the adjustment of any parameter values. The model also accounts explicitly for the effects of serine protease inhibitors on thrombin generation. Simulations to define the optimum affinity profile of an inhibitor in this system indicate that for an inhibitor simultaneously potent against VIIa, IXa, and Xa, inhibition of thrombin generation decreases dramatically as the affinity for thrombin increases. Additional simulations show that the reason for this behavior is the sequestration of the inhibitor by small amounts of thrombin generated early in the reaction. This model is also useful for predicting the potency of compounds that inhibit thrombosis in rats. We believe that this is the first mathematical model of blood coagulation that considers the effects of exogenous inhibitors. Such a model, or extensions thereof, should be useful for evaluating targets for therapeutic intervention in the processes of blood coagulation.

A mathematical model has been developed to simulate the generation of thrombin by the tissue factor pathway. The model gives reasonable predictions of published experimental results without the adjustment of any parameter values. The model also accounts explicitly for the effects of serine protease inhibitors on thrombin generation. Simulations to define the optimum affinity profile of an inhibitor in this system indicate that for an inhibitor simultaneously potent against VIIa, IXa, and Xa, inhibition of thrombin generation decreases dramatically as the affinity for thrombin increases. Additional simulations show that the reason for this behavior is the sequestration of the inhibitor by small amounts of thrombin generated early in the reaction. This model is also useful for predicting the potency of compounds that inhibit thrombosis in rats. We believe that this is the first mathematical model of blood coagulation that considers the effects of exogenous inhibitors. Such a model, or extensions thereof, should be useful for evaluating targets for therapeutic intervention in the processes of blood coagulation.
The clotting of blood is an exquisitely complex process. The simultaneous requirements for the free flow of blood under normal conditions and rapid clotting to prevent blood loss in the case of injury require a delicate balance between clot formation and clot lysis. Taken together, these processes involve more than 30 proteins, at least 10 of which are serine proteases. Inhibitors of serine proteases occur as natural anticoagulants, and natural and synthetic inhibitors have been widely studied for therapeutic applications in the prevention of thrombosis (1).
When choosing candidate serine protease inhibitors for therapeutic use, the choice of a target enzyme may be less than obvious. Many of the enzymes have multiple activities, there is positive and negative feedback regulation, and there are alternative pathways for activation and inactivation. In addition, serine protease inhibitors (especially synthetic compounds) generally have a spectrum of activities due to the high degrees of homology among the blood coagulation factors. This raises the question of whether the best antithrombotic compound would be specific for a single enzyme or show affinities for several serine proteases (2).
Mathematical modeling can help us understand such complex systems, but published models of blood coagulation suffer from the following limitations: consideration of only a small part of the coagulation cascade (3)(4)(5)(6)(7), empirical description of interactions for which molecular mechanisms were known (3,8,9), determination of some (4,9,10) or all (3,11) of the parameter values from the experimental data to which the model predictions were then compared (curve fitting), and the absence of comparisons of model predictions to experimental data (6,(12)(13)(14). In the best of these studies (4,5,7), most or all of the parameter values were independently determined, model predictions were compared to experimental results that had not been used to estimate parameter values, and the model was used to understand and interpret unexpected experimental observations.
An experimental system to study the tissue factor pathway to thrombin (15) and a mathematical model of this system (9) were recently described. The human tissue factor (TF) 1 -VIIa complex was used as the initiator of coagulation; human cofactors V and VIII, serine protease zymogens IX and X, and prothrombin were used at their normal plasma concentrations.
In this report, we describe a substantial revision of the mathematical model (9) of this experimental system (15). The revised model gives reasonable predictions of the reported experimental results with no adjustment of parameter values. We then extended this revised model to include the effects of serine protease inhibitors with affinity for any or all of the factors VIIa, IXa, Xa, and thrombin. The extended model was used to predict the affinity profile of the optimum inhibitor in this experimental system and to predict the potency of compounds that inhibit thrombosis in rats. In common with the best studies described above, all the model parameter values are independently determined, model predictions are compared to experimental data, and the model is used to explain unexpected results. In addition, our model is as comprehensive in its coverage of the clotting reactions as any other published model, and it is the first to include the effects of exogenous inhibitors.

MATERIALS AND METHODS
Model Computations-The model computations were done using a FORTRAN program that numerically integrates a set of coupled, ordinary differential equations. A differential equation describes the change in concentration with time of each of 36 of the 43 distinct species in the model; concentrations of the other seven species are calculated from equations for the conservation of the initial amounts of the components added to the reaction mixture (TF-VIIa complex, V, VIII, IX, X, prothrombin, and inhibitor). The equations were integrated using backward differentiation formula routines suitable for stiff systems (Numerical Algorithms Group, Oxford, UK). Routine D02EBF was used to calculate concentrations versus time for a specified length of time; routine D02EHF was used to calculate concentrations versus time until a specified end point (80% of the final thrombin activity) was reached. All subroutine and main program calculations were done in double precision. Initial concentrations of the reaction components were the normal plasma concentrations as described (15) No interactions were assumed to be in equilibrium, but unless on-and off-rates were measured independently, the on-rate was assumed to be rapid (9). The off-rate was calculated from the dissociation constant, K d , (or the Michaelis constant, K m , if only enzyme kinetic parameters were reported) as k off ϭ k on ϫ K d (or K m ).
Inhibition of Arterial Thrombosis in the Rat-All procedures used in this study were approved by the DuPont Merck Pharmaceutical Company Institutional Animal Care and Use Committee. 400 -450-g male CD rats from Charles River (Kingston, NY) were anesthetized with ketamine (80 mg/kg intramuscular) and xylazine (10 mg/kg intramuscular). Additional anesthetic was administered as needed. The jugular vein was catheterized for compound infusion and the femoral artery was catheterized for blood sampling. A flow probe (Transonic Systems, Inc.) was placed around the carotid artery. The blood flow (ml/min) was monitored for 5 min prior to the start of infusion. After a 15-min infusion, 10 l of a 25% FeCl 3 solution was placed on the carotid artery proximal to the flow probe. The flow was recorded at 5-min intervals for 60 min unless an occlusive thrombus formed, in which case the time of occlusion was noted.

RESULTS AND DISCUSSION
Model Formulation-The reaction scheme described by the model is summarized in Fig. 1. As in the experimental system (15), the TF-VIIa complex is assumed to be fully formed at the start of the experiment, so descriptions of the binding of VII(a) to TF and the activation of VII are not included. The TF-VIIa complex activates both IX (17) and X (18). IXa forms complexes with VIII (19) or VIIIa (20); the VIIIa-IXa (tenase) complex also activates X (21,22). Loss of VIIIa activity occurs spontaneously (22), and VIII and VIIIa are degraded by IXa (23). Xa and Va form the prothrombinase complex (24, 25), which activates prothrombin (II) by way of an active intermediate, meizothrom-bin (mIIa) (26). Xa, mIIa, and thrombin (IIa) each activate the cofactors V and VIII (27)(28)(29)(30). VIII activated by Xa and VIII activated by (m)IIa are treated as different species, as it has been shown that they differ in both cofactor activity and stability (22). An inhibitor (Fig. 1, I) can bind reversibly to each of the serine proteases VIIa, IXa, Xa, and (m)IIa.
Whenever possible, parameter values determined at 37°C using human proteins were used in the model. Kinetic constants measured in the presence of phospholipid vesicles were used instead of those measured in the presence of cell surfaces, consistent with the experimental system being modeled. The parameter values used in the model are shown in Table I.
Prediction of Experimental Results in the Absence of Inhibition-The model was used to predict the time course of activation of the factors included in the experimental system (15). A comparison of the model predictions to experimental data is shown in Fig. 2. The qualitative agreement is quite good; the cofactors V and VIII are completely activated before there is any significant amount of thrombin generation, and there is little activation of IX and X throughout the course of the experiment. Quantitatively, there are several differences: the predicted thrombin activity does not show the peak arising from the accumulation of the active intermediate meizothrombin (15,31), the rapid activation of V occurs too late, the activation of VIII is too abrupt, and too much X is activated. However, given the good qualitative agreement with the experimental data (15) and the dangers inherent in adjusting parameter values (32,33), we chose to extend the model to consider the effects of serine protease inhibitors without any revisions of the model parameter values.
Prediction of the Most Potent Inhibitor-The model was extended to allow for competitive, active site inhibition of the serine proteases VIIa, IXa, Xa, and (m)IIa. The affinity of the inhibitor for each of the enzymes was defined by a separate dissociation constant (K i ). The model was then used to predict the optimum affinity profile of an inhibitor in the experimental system (15).
A concentration of 10 nM inhibitor was assumed, and the K i values for each enzyme were allowed to vary over ranges similar to those observed for compounds in our library of inhibitors. These ranges were 0.1-1000 nM for VIIa, 1-1000 nM for Xa, and 0.01-1000 nM for (m)IIa. We currently do not have a good assay for the inhibition of IXa, so K i values were assumed to range from 1-1000 nM. Simulations were run to determine the time required for the generation of 80% of the final thrombin activity. We used the ratio of the inhibited time to the uninhibited time as a measure of inhibitor potency; the higher this ratio, the more potent the inhibitor. The K i ranges were covered in half-log steps, with simulations run at all of the 4,851 K i combinations.
The results shown in Fig. 3 lead to several conclusions. First, the predicted thrombin generation times are relatively insensitive to changes in the K i values until these values are on the order of 1 nM or less (except as noted below). Second, inhibitors with affinity for only one of the serine proteases are not very potent compared to inhibitors with a spectrum of affinities. Third, compounds that simultaneously have high affinity for VIIa, Xa, and thrombin are predicted to be very potent, regardless of the IXa K i (Fig. 3, bottom row, back corners). Fourth, and most surprising, for compounds that simultaneously have high affinity for VIIa and IXa, inhibition of thrombin generation decreases as the affinity for thrombin increases from 1000 to 5 nM (Fig. 3, left column, back right planes). This effect is especially dramatic for compounds that also have high affinity for Xa (Fig. 3, bottom left, back right plane).
To explain this unexpected result, we used the model to determine where the inhibitor was bound. In Fig. 4, the amount of inhibitor bound to each of the serine proteases is shown for two inhibitors, each very potent against VIIa, IXa, and Xa; one has high affinity for thrombin, and the other has low affinity for thrombin. For the case of high affinity for thrombin (Fig.  4A), the inhibitor binds to thrombin in preference to IXa and Xa as soon as a small amount of thrombin is formed. With a reduction in the amount of inhibitor bound to IXa and Xa, the activation of prothrombin occurs very rapidly. For the case of low affinity for thrombin (Fig. 4B), the inhibitor remains bound to IXa and Xa even after a substantial amount of thrombin is formed, and the activation of prothrombin occurs only gradually.
Prediction of the Potency of Inhibitors of Thrombosis in the Rat-To determine if the model predictions had any relevance to in vivo coagulation, we predicted the potency of compounds tested as inhibitors of arterial thrombosis in rats. Inhibitor potency was characterized by the OT 30 , the dose of inhibitor required to extend to 30 min the occlusion time of an injured carotid artery (the control occlusion time was 20 min). Attempts to correlate the in vivo activity with combinations of either serine protease inhibitor K i values or inhibitor potencies in in vitro clotting assays (prothrombin time, activated partial thromboplastin time (APTT), and thrombin time) led to the conclusion that the APTT was the best predictor (Fig. 5A). However, several compounds were significantly more potent in the rat than would be predicted by the APTT.
Model predictions of the thrombin generation time ratio for 25 compounds tested in the rat were compared to the OT 30 values (Fig. 5B). The ranges of K i values for the compounds tested in rats were 0.65-1330 nM for VIIa, 0.84 -920 nM for Xa, and 0.024 -206 nM for thrombin. Because we did not have IXa K i values for these compounds, the simulations were done assuming that the IXa K i was equal to the VIIa K i . Similar results were obtained when the IXa K i was assumed to be equal to the Xa K i ; predictions were somewhat worse when all the inhibitors were assumed to be active (K i ϭ 1 nM) or inactive (K i ϭ 1000 nM) against IXa (data not shown). All of the compounds predicted to be potent by the model (ratio Ͼ 1.2) were potent in the rat (OT 30 Ͻ 10 g/kg/min); all of the compounds impotent in the rat (OT 30 Ͼ 10 g/kg/min) were predicted by the model to be impotent (ratio Ͻ 1.2). In addition, the model correctly predicted the potency of the APTT outliers.
Conclusions-This work demonstrates that a mathematical model of a significant portion of the blood coagulation system can lead to an unexpected conclusion about the desired activity of an inhibitor in that system. The model also shows utility in predicting the in vivo potency of serine protease inhibitors, despite the fact that the model was designed for humans and the in vivo data were collected in rats. These results, together with the ability of the model to provide an explanation for the FIG. 4. Distribution of inhibitor. A, an inhibitor with high affinity for thrombin. B, an inhibitor with low affinity for thrombin. The lines show the generation of thrombin activity (right vertical axes). Areas are proportional to the amount of inhibitor bound to each species (left vertical axes): free inhibitor (hatched area), bound to thrombin or meizothrombin (graphed area), bound to Xa (shaded area), bound to IXa (cross-hatched area), and bound to VIIa (white area). In each case, only 1% of the total inhibitor is shown; the rest is free or bound to thrombin or Xa.

FIG. 5. Potency of inhibitors of rat arterial thrombosis.
A, correlation with the inhibitor concentration required to double the APTT. B, correlation with model predictions. Compounds that are more potent in the rat than would be predicted by the APTT are labeled with identifying codes. Note that on the horizontal axes, potency increases from right to left for APTT and from left to right for the model predictions.
predictions, suggest that such a model could be useful for identifying targets for therapeutic intervention in the processes of blood coagulation.