Introduction

Recent experimental studies have shown that the plasma levels of blood clotting proteins vary considerably, both with age, as well as between healthy individuals of the same age1. The haemostatic system in neonates and children also shows considerably different functional characteristics compared to adults2,3,4. Therefore, clinical studies of the adult haemostatic system may not be transferrable to the treatment of thrombotic and haemorrhagic disorders in young patients.

The high variability of blood clotting protein abundances is an intrinsic property of the haemostatic system rather than being caused by flaws of the available experimental methods. The complexity of the haemostatic system in combination with the strong fluctuations of protein concentrations between individuals and with age makes determining the impact of experimentally measured differences in plasma levels of individual haemostatic proteins extremely challenging. Mathematical modelling is a useful tool for investigating the relative importance of a variety of complicated interdependencies in complex systems, and has also been successfully applied to the haemostatic system. Early qualitative models focused on investigating the mechanism of haemostatic response in rapidly forming a blood clot shortly following an injury5,6. More recently, this behaviour known as excitability was investigated in mathematical models by Jesty et al.7 and Beltrami and Jesty8.

The next generation of haemostasis models was based on more detailed representations of the biochemical reaction network of haemostatic proteins. Representing chemical reactions by mass-action kinetics allows for a more data-driven approach by measuring rate constants experimentally. Thus, given that sufficient experimental data is available, these models not only capture qualitative aspects such as excitability but may also be used to study haemostasis quantitatively. The model by Hockin, et al.9, henceforth referred to as the Hockin-Mann model, consists of a system of 34 ordinary differential equations that represents the dynamics of 44 biochemical reactions of the haemostatic network by mass action kinetics. Each equation represents the time-dependent dynamics of one of the haemostatic proteins, or a complex of proteins. For a given set of initial concentrations, the Hockin-Mann model enables us to calculate how these haemostatic factors change over time.

The Hockin-Mann model continues to be used as the starting point or a building block in new models that are being developed. Wajima et al. extend a previous version of the Hockin-Mann model in order to study several experimental tests of the haemostatic system, its response to different treatments (warfarin, heparin, vitamin K) and to perturbations by Taipan snake venom10. Models based on the Hockin-Mann model have been developed for studying the effect of anti-coagulants such as e.g. the novel anti-coagulant (NOAC) rivaroxaban11,12. Most recently Mitrophanov et al. used the Hockin-Mann model for investigating the effects of acidosis on thrombin generation13.

Gaining systematic insight into the behaviour of highly detailed models within the order of tens or even hundreds of equations is a difficult problem. Danforth et al. follow a computational approach for investigating the parameter sensitivity of the Hockin-Mann model14. In this study they examine the sensitivity of thrombin generation as well as the behaviour of the full model on variations of individual rate constants in a range from 10 to 1000% of the reference values given in the original publication9. In Danforth, et al.15 the statistical approach from14 is used to assess the sensitivity of thrombin generation to variations in blood haemostatic factors, both utilising synthetically generated combinations of concentrations, as well as data from healthy adults (patients with haemophilia A and a population receiving warfarin). The study presented here follows a similar approach as Danforth, et al.15 but we investigate healthy individuals from different age groups rather than different adult populations. Most recently, Dunster and King16 demonstrated, via a thorough mathematical analysis of one of the earlier models of haemostasis by Willems, et al.17, how, taking advantage of distinct time scales, simplified models for different phases of thrombin generation can be derived. These phases can be related to the initiation, propagation and termination phase of thrombin generation. Because the simplified models are much easier to analyse than the full model, the dominant mechanisms in each of six different phases can be clearly identified.

This study aims to investigate age-dependent changes of the blood clotting system. In the absence of knowledge regarding changes to the network of biochemical reactions itself, we use the Hockin-Mann model as a representation of the haemostatic system of all age classes. In order to test the hypothesis that age-dependent changes are due to differences in the plasma levels of blood clotting proteins we use age-stratified data by Attard, et al.1. In contrast to previous studies we explicitly account for variability between individuals by parametrising the initial concentrations of haemostatic factors with data from individual patients – in fact, our results show that choosing the levels of haemostatic factors based on data aggregated for different age groups would be misleading.

Different to earlier experimental studies of the age-dependent blood clotting system which were based on functional assays2,3, Attard et al. have presented the first comprehensive data set that provides quantitative age-stratified protein levels. Quantitative assays are clearly more suitable than functional data for obtaining blood clotting factor concentrations necessary for parametrising models of haemostasis because functional assays only provide a proxy for concentrations. For this reason, caution is required when comparing with older data sets based on functional assays.

To the best of our knowledge, this is the first time that a model of the haemostatic system is parameterised with age-stratified haemostatic protein abundance data. This allows us to link observed age-dependent differences in the concentrations of individual haemostatic factors to functional implications for the activity of the haemostatic system.

Results

In order to account for the strong inter-individual variability in different age groups we parametrised the initial protein concentrations in the Hockin-Mann model with levels of prothrombin, FV, FVII, FVIII, FIX and FX measured in individual patients1. For all age groups, we compared TGCs calculated from the mean levels of these haemostatic proteins, see Supplementary Material, Fig. S1, with TGCs calculated for individual patients (Fig. S2). In analogy to experimental studies, four indices were used for quantifying different aspects of a given thrombin curve (Fig. 1): The lag time (LAG) – the time it takes until thrombin exceeds a certain fraction of its maximum - and the time it takes until thrombin reaches its maximum, the time to thrombin peak (TTP) indicate how quickly thrombin is generated in response to stimulation by TF. The other two indices, the maximum total thrombin generation (MAX) and the area under the thrombin curve (AUC), are measures for the strength of this response.

Figure 1
figure 1

Thrombin generation curve simulated using the Hockin-Mann model with parameters corresponding to an individual adult subject. Four indices that were calculated for comparison of TGCs calculated for different subjects are shown; these are lag time (LAG), time to peak (TTP), maximum thrombin concentration (MAX) and area under thrombin curve (AUC). For the TGC plotted here we have LAG = 554 s, TTP = 742 s, MAX = 570 nM and AUC = 2191 nM · min.

The results show that a TGC calculated for the mean plasma levels of haemostatic proteins of a particular age class fails to accurately represent some of the properties of the TGCs calculated for individual patients. For example, whereas the TGCs calculated from mean plasma levels indicate that the TTP is approximately t = 300 s, especially the results for teenagers and adults suggest that no such pattern exists because the TTP strongly varies within age groups. It was therefore considered essential to carry out further analysis based on TGCs calculated for individual patients. The results for individuals in each of the age groups for LAG, TTP, MAX and AUC are presented in Fig. 2. The levels of LAG and TTP show no significant differences between different age groups. In contrast, both thrombin maximum (MAX) and area under thrombin curve (AUC) exhibit a significant increase from a low level observed for neonates at day 1 and day 3 compared to other age groups. Also, the higher value of MAX and AUC for children less than one year old until adults is statistically similar (Fig. 2).

Figure 2
figure 2

Age-dependent changes of thrombin generation curve. Whereas lag time (LAG) and time to thrombin peak (TTP) seem not to vary with age, thrombin maximum (MAX) and the area under the thrombin curve (AUC) increase roughly two-fold with age.

We then investigated the relationship between MAX and AUC and the initial concentrations of the blood clotting factors reported by Attard, et al.1. Figure 3 shows that MAX appeared to increase linearly with prothrombin concentration. Figure 4 demonstrates that AUC was solely determined by the concentration of prothrombin and was not dependent of the concentrations of all other haemostatic factors.

Figure 3
figure 3

The thrombin maximum MAX increases approximately linearly with prothrombin. Similar relationships are seen for the other factors but the variation for the predicted thrombin maxima is much higher than the dependency on prothrombin shown here. The prothrombin level is both shown on the original scale in units of IU/mL from Attard, et al.1 as well as in units of nM obtained by conversion according to Table 1.

Figure 4
figure 4

The area under the thrombin curve (AUC) is solely determined by the concentration of prothrombin. The prothrombin level is both shown on the original scale in units of IU/ml from Attard, et al.1 as well as in units of nM obtained by conversion according to Table 1.

Age-related differences in the concentrations of haemostatic factors, based on Attard, et al.1 are shown in Fig. 5. The concentrations of factor VII and factor IX clearly increased with age. The concentration of factor VIII did not change with age whereas for concentrations of factors V and X the dependence with age is less obvious. However, for the concentration of prothrombin we observed a clear jump from a low level of approximately 0.5 IU/ml to more than 1.0 IU/ml in the other age groups. Hence there were essentially two different levels of prothrombin concentrations for neonates compared with the other age groups, in contrast to factors VII and IX whose concentrations increased continuously with age.

Figure 5
figure 5

Age-dependent variability of blood clotting factors (see Attard, et al.1). Haemostatic factor abundances are both shown on the original scales in units of IU/mL or percentages (%) as well as in units of nM obtained by conversion according to Table 1.

Discussion

Developmental haemostasis is a concept that describes multiple and complex age-specific differences in the haemostatic system. The overall impact of these differences is believed to provide protection for the neonate and child in terms of response to bleeding and clotting stimuli. However, predicting the response to diseases that affect haemostasis, or to the multitude of new anticoagulant drugs that are becoming available for clinical care in children is challenging. Conducting large-scale trials of these drugs in children, such as those performed in adults, has proven to be very difficult. Thus, the development of a mathematical model that simulates the age-appropriate haemostatic system would be very advantageous.

We presented a novel mathematical modelling approach that enables us to link age-related differences of haemostatic factors in neonates, children and adults to functional differences in their haemostatic system. The key contribution of this study is that data from individual patients has been used for parametrising the Hockin-Mann model. This enables us to appropriately account for the strong inter-individual variability present in different age groups, in contrast to the more common approach of obtaining aggregated parameters from population statistics—compare Figs S1 and S2 in the Supplementary Material. Because we calculate a thrombin curve for each individual patient, unlike previous studies, we are in a unique and novel position to draw conclusions regarding the variability of haemostatic activity based on realistic distributions of blood clotting factors from individual patients. We note that this is different to simply investigating the effect of upper and lower bounds of individual blood clotting factors because this ignores correlations between different factor levels—the strong variations of the correlation structure between various haemostatic factors is shown in Fig. S3 of the Supplementary Material.

Our approach is based on the hypothesis that the biochemical reactions of the haemostatic system are similar between different subjects and remain unchanged with age so that differences in haemostatic factor concentrations are the main source of observed age-related differences in haemostatic activity. Under these considerations, a model parameterised with age-stratified abundances of haemostatic factors (while leaving the reaction rate constants unchanged) fully accounts for age-related differences in the haemostatic system. The implications of these assumptions, in particular, the role of inhibitors, will be discussed in more detail below. Moreover, we will present a preliminary data set collected by two of the co-authors that confirms the qualitative differences between neonates and older age groups predicted by the model.

Our main finding is that for neonates at day 1 and day 3 post-birth, MAX and AUC were significantly lower than in older age groups indicating a lower activity of the haemostatic system. Importantly, our model enables us to identify a significantly lower prothrombin level as the most likely cause for this observation. This simple relationship between prothrombin and indicators of blood clotting activity bypasses the complex correlation structure between haemostatic factors apparent in Fig. S3 of the Supplementary Material.

First, we found that whereas LAG and TTP are unchanged with age, MAX and AUC are significantly different between neonates and all other age groups. Both are elevated from a low level of approximately 200 nM or 500 nM · min, respectively, for neonates at day 1 and day 3 to significantly higher levels of 400 nM or 1000 nM · min for all other age groups. In Fig. 6 we demonstrate that the increase of AUC for age groups older than day 3 predicted by the model is consistent with age-stratified measurements of AUC (Ignjatovic and Monagle, unpublished). This confirms that although some possible age-dependent changes in the haemostatic system are not accounted for, the model behaviour related to the indices considered in this study is nevertheless qualitatively correct. We also remark that the fact that both LAG and TTP are unchanged between different age classes indicates that using an different threshold concentration for defining the LAG time will not affect our results.

Figure 6
figure 6

Age-stratified experimental data for the area under the thrombin curve (AUC) is compared with the model. Both data and model simulations are scaled by the mean adult level to demonstrate age-related differences. The data (white) confirms the approximate two-fold relative increase of AUC to the adult level after day 3.

Second, we observed that prothrombin concentration is a strong predictor for MAX and AUC. This has been reported in both the modelling15 as well as in the experimental literature18 albeit in different settings as our study. Danforth, et al.15 showed that both MAX and AUC depend most sensitively on prothrombin in their simulation study of the Hockin-Mann model. For the Willems, et al.17 model, Dunster and King16 demonstrate that both MAX as well as AUC increase linearly with prothrombin concentration and give Duchemin, et al.18,19 as a reference which experimentally confirms this observation from the mathematical model. None of the above studies refers to different age groups.

Third, by combining our observations, we identify the lower prothrombin concentration as the cause for the lower values of MAX and AUC in neonates. Whereas our model shows that MAX and AUC are lower in neonates at day 1 and day 3 post-birth than in the older age groups, we know from the Attard, et al.1 data that neonates also have a decreased prothrombin concentration compared to adults. Taking into account that MAX and AUC are largely determined by prothrombin we conclude that the lower thrombin production in neonates is primarily due to lower prothrombin availability. Although a strong influence of prothrombin and indices related to thrombin abundance is not unexpected, it is interesting that the model suggests that other blood clotting factors hardly play any role in determining MAX and AUC. Given that the haemostatic network is in general characterised by a high level of complexity this observation is potentially very useful because it suggests that targeting prothrombin may to a great extent be sufficient for regulating the amount of thrombin generated.

Because our study relies strongly both on data from individuals as well as quantitative rather than functional assays, the influence of several blood clotting proteins such as the inhibitors Tissue Factor Pathway Inhibitor (TFPI) and Antithrombin as well as Fibrinogen and thrombomodulin could not be accounted for because they were not measured by Attard, et al.1. To the best of our knowledge, comparable age-stratified data sets for these blood clotting factors based on quantitative assays are currently not available. Moreover, although Attard et al. measured Factor XI (FXI), FXII as well as protein C and protein S, these blood clotting proteins could not be included because they are not accounted for by the Hockin, et al.9 model. Before we discuss the possible roles of different inhibitors, we emphasise that our approach to age-dependent modelling can be easily applied to an extended model as more comprehensive data becomes available.

Andrew, et al.20 observe that the inhibitor Antithrombin was decreased in neonates by 60% and it has been suggested that this enabled normal haemostatic function. But the role of inhibitors in the neonate haemostatic system remains controversial because whereas Antithrombin is decreased, the inhibitor α2-Macroglobulin (α2M) is increased. In fact, our results are consistent with a recent article by Kremers, et al.21 which attributed the decreased thrombin generation in young patients to decreased prothrombin conversion and to a lesser extent to elevated levels of α2M. Nevertheless, the issue of clarifying the relative importance of the different inhibitors in neonates clearly deserves further investigation. Because this requires extending the Hockin-Mann model by the inhibitor α2M we leave the study of this interesting question to future work.

In this study, we aimed to represent the experimental assay for determining the activity of the haemostatic system via measurement of the thrombin generation curve. For representing the in vivo system more detailed models of haemostasis have been developed that account for the fluid dynamics of blood within the constraints of the blood geometry as well as the interaction with platelets, see, e.g. Cito et al.22 or the relevant chapters in Ambrosi, et al.23.

In summary, we have demonstrated that parametrising a mathematical model with data of individual patients from different age groups is a useful tool for investigating age-dependent differences in the haemostatic system. As additional data become available, our approach can be easily transferred to more comprehensive models of the in vitro or in vivo haemostatic system which will enable us to explore the complex system underlying developmental haemostasis in more detail. Considering the considerable uncertainty regarding treatment of young patients with coagulation disorders, an extremely useful direction will be to incorporate the interactions with haemostatic drugs in the age-stratified model presented here.

Methods

In order to account for the strong variability between individuals apparent in the data by Attard et al., we simulated the Hockin-Mann model with different initial conditions for each patient sample. This approach has been used earlier for investigating inter-individual variability in adults15. Based on these simulation results we then quantified the activity of the haemostatic system by four commonly used indices that are derived from the simulated thrombin concentration over time, the so-called thrombin generation curve (TGC): the lag time (LAG), the time to thrombin peak (TTP), the maximum thrombin concentration (MAX) and the area under the thrombin curve (AUC). Finally, we investigated which individual haemostatic factors had the strongest influence on each of the four metrics and related our results to differences in the availability of haemostatic factors between different age groups.

Parametrisation

We parameterised the Hockin-Mann model using age-stratified data from the study by Attard, et al.1, that consists of measurements for factors II and V given as international units (IU/ml) and factors VII, VIII, IX and X reported as percentages (%). Data for these six haemostatic factors were available for patients from seven age groups, 10 samples each for neonates at day 1 and day 3 of age and 20 samples each for the remaining age groups (younger than 1 year, 1–5 years, 6–10 years, 11–16 years, adults). The distributions over the individual age groups are summarised in Fig. 5 but note that in this study we used individual samples—i.e. the raw data. Measurements were converted from the original units to concentrations by choosing the average adult results measured by Attard et al. as a reference parameter set, see Table 1. We required that these average adult results correspond to the initial concentrations of the Hockin-Mann model, see Table 3 in Hockin, et al.9 or Table 1, which were chosen as average values in human plasma. Thus, measurements from an individual subject were converted to concentrations by scaling with respect to the average result for adults, as reported in Attard, et al.1. For all other parameters such as rate constants we kept the original values reported in Tables 2 and 3 in Hockin et al.9. TFPI, VIIa and ATIII were set to the plasma concentrations stated in the main text of Hockin et al.9. The equations of the Hockin-Mann model can be derived from Table 19. We have stated them explicitly in Section 3 in the supplementary material of this article.

Table 1 Average, minimum and maximum levels of haemostatic factors for the data from Attard et al.1.

Simulation

Simulations were run with our own implementation of the Hockin-Mann model using a numerical solver of the GNU Scientific Library (GSL)24. For t = 0 initial conditions for factors II, V, VII, VIII, IX and X were chosen according to a sample from the data by Attard, et al.1 as described previously. Thrombin generation was initiated by initialising Tissue Factor (TF) at [TF] = 5 M which is a concentration commonly used in experimental studies of thrombin generation. The model was simulated until tend = 5000 to ensure that thrombin generation proceeded through all phases from initiation via propagation to termination. Total thrombin is represented by two fractions in the Hockin-Mann model, meizothrombin (mIIa) and thrombin (IIa). Due to the higher activity of meizothrombin, the total amount of thrombin over time, the thrombin generation curve (TGC), is calculated by [IIa] (t) + 1.2[mIIa] (t) – for details, see Hockin, et al.9. A representative example (where termination occurred approximately at t = 1400 s) is shown in Fig. 1.

Indices of thrombin generation

In analogy to experimental studies, thrombin generation was characterized quantitatively by four metrics: The lag time (LAG) and the time to thrombin peak (TTP) indicate how quickly thrombin is generated in response to stimulation by TF. In contrast, the remaining two indices, the maximum total thrombin generation (MAX) and the area under the thrombin curve (AUC), are measures for the strength of this response. The maximum total thrombin concentration is defined as

$${\rm{MAX}}\,:\,=\mathop{{\rm{\max }}}\limits_{t\in [0,{t}_{{\rm{end}}}]}\,([{\rm{IIa}}]\,(t)+1.2[{\rm{mIIa}}]\,(t))$$
(1)

where [IIa] (t) and [mIIa] (t) are thrombin and meizothrombin concentration over time, respectively, and tend is the end time of the simulation. The time to thrombin peak (TTP) is then defined as the elapsed time after stimulating the system with tissue factor until the time the concentration MAX is attained. Similarly, the lag time LAG is defined as the time until 1/6 of the maximum thrombin concentration MAX is reached (Note that the threshold concentration used for defining varies in the literature. We explain in the Discussion why the exact choice does not influence our results). The area under the thrombin curve AUC is given by the integral over the thrombin curve

$${\rm{AUC}}:\,={\int }_{0}^{{t}_{{\rm{end}}}}\,[{\rm{IIa}}]\,(t)+1.2[{\rm{mIIa}}]\,(t)dt$$
(2)

All four indices are visualised in Fig. 1.

Analysis of age-stratified simulation results

All statistical analyses of the simulation results were performed in the R software25. From the simulations based on the age-stratified data from Attard et al. we obtained age-stratified results for the four indices LAG, TTP, MAX and AUC. Differences between age groups were determined by comparing the distributions of each index. For indices that showed age-dependent differences we compared the dependency on measurements of the six haemostatic factors II, V, VII, VIII, IX, X by Attard, et al.1. Finally, we inspected the data from Attard et al. for differences in the measurements of the six haemostatic factors II, V, VII, VIII, IX, X for different age groups.