Individual variation in unfractionated heparin dosing after pediatric cardiac surgery

We aimed to identify attributing factors to the interindividual variabilities of the infusion rates in unfractionated heparin therapy. We included patients who required unfractionated heparin therapy to achieve the target APTT after cardiac surgery between May 2014 and February 2018. Fifty-nine patients were included, of whom 8 underwent Blalock-Taussig shunt; 27, Glenn procedure; 19, Fontan procedure; 3, mechanical valve replacement; and 2, Rastelli procedure. Previously reported variables that influenced the response to unfractionated heparin treatment were initially compared, which included age; weight; sex; type of surgery; platelet count; fibrinogen, antithrombin III, total protein, albumin, alanine transaminase, and creatinine levels; and use of fresh frozen plasma. The type of surgical procedure was found to be significantly associated with the differences in heparin infusion rate (P = 0.00073). Subsequently, the variance explained by these factors was estimated through a selection based on the minimum Akaike information criterion value; models constructed by various combinations of the surgery types were compared. The model including the Blalock-Taussig shunt, Glenn procedure, and mechanical valve replacement showed the highest summed variance explained (29.1%). More than 70% of the interindividual variability in initial heparin maintenance dosing was unexplained.

heparin therapy is more strongly associated with anti-Xa activity, measurement of anti-Xa activity as a routine practice is not necessarily the standard of care due to lack of evidence 1-4 . Intensivists frequently have difficulty predicting responses to unfractionated heparin therapy in each patient by observing contributing factors to the interindividual variabilities of the unfractionated heparin dosing required for each patient. These interindividual variabilities are well known, and previous studies reported multiple factors, including obesity, aging, hepatic or renal disease, altered production of heparin-binding proteins, general heparin resistance, antithrombin deficiency, increased heparin clearance, elevated levels of heparin-binding proteins, and increased plasma levels of factor VIII, fibrinogen, and platelet factor 4 [5][6][7][8][9][10][11] . Recently, several studies have been conducted to explain interindividual variabilities of unfractionated heparin therapy, one of which is by Al-Sallami et al. who developed a population pharmacokinetics (PK)-pharmacodynamics (PD) model in pediatric patients during cardiac angiography and showed that fat-free mass was a significant covariate for clearance 12 . The study by Delavenne et al. developed a population PK-PD model for adults patients during cardiopulmonary bypass (CPB) and reported that the inclusion of body weight in their model decreased the interindividual variabilities of clearance and central compartment volume 13 . However, existing unfractionated heparin dosing algorithms do not incorporate these factors.
Thus, in this study, we aimed to identify the attributing factors to the interindividual variabilities of the heparin infusion rate upon achieving the target APTT, which is one of the most important phenotypes in the intensive care unit and can also be used to quantify responses to unfractionated heparin therapy.

Methods
We aimed to identify the attributing factors to the interindividual variabilities of the heparin infusion rate upon achieving the target APTT. We performed a retrospective observational cohort study to examine the medical records of patients admitted to the pediatric intensive care unit (PICU) in a single tertiary care center (National Center for Child Health and Development [NCCHD], Tokyo, Japan) between May 2014 and February 2018. This study was approved by the ethics review board of the National Center for Child Health and Development, Tokyo, Japan (Receipt No. 2248). Written informed consent was waived by National Center for Child Health and Development [NCCHD], Tokyo, Japan because of the retrospective design. All methods were carried out in accordance with relevant guidelines and regulations. The inclusion criteria of this study were as follows: (1) patients who were receiving unfractionated heparin therapy; (2) patients admitted to the PICU after cardiac surgery (Blalock-Taussig shunt, Glenn procedure, Fontan procedure, biological or mechanical valve replacement, or Rastelli procedure, as we provide postoperative heparin therapy with target APTT ranges of 40 to 65 s, depending on the surgical procedure as shown in Supplementary Table 1, only for patients with postsurgical status); and (3) patients who required titration of the heparin infusion rate to achieve the target APTTs. We considered that the target therapeutic ranges of APTT had been reached when the APTT was within ± 5 s of the target APTT ranges in the Supplementary Table 1. We excluded patients who receiving ECMO or any dialysis. In addition, we excluded patients who did not reach the target APTT after 72 h of starting heparin therapy, as we usually transition to either warfarin or aspirin from unfractionated heparin for anticoagulation therapy once enteral feeding is successfully started. If patients had multiple admissions during the study period, we used the medical records from the first admission or those with full descriptions of the postoperative course in the progress notes.
Drug administration. In the PICU, pediatric intensivists started unfractionated heparin therapy at an initial rate of 10 units/(kg h). They titrated the unfractionated heparin dosage to achieve the target APTT ranges. All the doses were recorded, including the start time of the infusions.

Definitions of outcomes.
Given the fact that a steady-state drug concentration in the blood is, in general, achieved after 4 or 5 half-lives of drug elimination, for unfractionated heparin with a half-life of 0.5-2 h, the steady state would be reached 2-10 h after treatment initiation 14 . More importantly, it is well known that APTT reaches a steady state after approximately 4 h in children 10 . We excluded patients who did not reach the target APTT after 72 h of starting heparin therapy as mentioned above. Therefore, we set the primary outcome for our study to be the initial maintenance dosage of unfractionated heparin (unit/[kg h]) given to a patient between 4-72 h of starting heparin therapy. In addition, we considered the infusion rate of heparin maintenance dosage, summing all the infusion rates of heparin running in lines such as arterial lines or central venous catheters to measure blood pressure, including systemic blood and pulmonary artery pressures. In the PICU, unfractionated heparin was administered via an arterial line at a rate of 2 units/h for patients whose body weights were < 5 kg, via an arterial line at 4 units/h for those whose body weights were ≥ 5 kg, and via a central venous catheter at a rate of 2 units/h for all the patients to prevent clotting. In addition, given that the initial APTT could be largely influenced by heparin therapy during surgery, we considered the target APTT as follows: the APTT (1) reached nadir after admission to the PICU and (2) subsequently reached the target. Data collection. The baseline demographic data included age (years), weight (kg), height (cm), sex, name of performed surgical procedures, and name of cardiac diseases for which those surgical procedures were collected. Laboratory data, concomitant drugs administered, number of transfusions given, and numbers/types of lines used, such as arterial lines or central venous catheter information, were also collected.

Statistical analyses.
The Shapiro-Wilkes test was used to evaluate the normality of the heparin infusion rate. Given that the histogram of the heparin infusion rates showed an abnormal distribution (P = 3.1E−06; Fig. 1), we log-transformed the heparin infusion rate for further analysis. For correlation analysis, the Spearman ρ value was used. www.nature.com/scientificreports/ As covariates of the interindividual variabilities of the infusion rate of heparin as maintenance therapy after cardiac surgery, the factors that were previously reported to be associated with heparin response/resistance and for which data are currently available were included in our study. We set the cutoff antithrombin III activity level to ≤ 60%, platelet count to > 300,000, and increased factor VIII and fibrinogen levels [5][6][7][8][9] . We set the cutoff the fibrinogen level to > 300 mg/dL in this study 15 . In addition, we investigated body weight, types of surgical procedure, total protein level (> 5.0 or < 5.0 g/dL), albumin level (> or < 3.0 g/dL), alanine transaminase level (> 100 or < 100 IU/L) to represent liver function; serum creatinine level (> 0.8 or < 0.8 mg/dL) to represent renal function, use of fresh frozen plasma (but limited to occasions within 48 h before reaching the target APTT), use of platelet transfusion (only within 48 h before the target APTT), and use of protamine (within 48 h before the target).
Subsequently, we selected covariates with significant associations that were determined by conducting nonparametric tests, specifically the Wilcoxon rank-sum test for differences between two groups and the Kruskal-Wallis test for comparing three or more variables. We then constructed a polynomial linear regression model and evaluated the total explained variance by using the variables with significant associations in the covariate selection to investigate how much interindividual variabilities of unfractionated heparin therapy could be explained by those covariates. When the types of surgical procedure to be accounted for in the model development are required, we used a dummy variable for the Rastelli procedure.
In addition, the Akaike information criterion (AIC) for each model was calculated in order to identify the optimal model in terms of possible prediction ability 16 . We selected only the polynomial regression models where the directions of the regression coefficients and values that multiply the predictor values were estimated to be the same as the direction of the regression coefficient in the single regression model for each variable.
We also estimated the variance explained, adjusting for the APTT in each model and regressing the effect of APTT, as APTT could also be influenced by multiple factors and variabilities could be due to APTT itself 7 . All statistical analyses were two-sided, and a P value of < 0.05 was considered significant. All analyses were performed using the R version 3.5.0 statistical software 17 . Boxplots were drawn using the R package ggplot2 18 .
Given that the target APTT is slightly higher after mechanical valve replacement, as shown in Supplementary

Results
Ninety-two unique patients received heparin therapy after cardiac surgery (Blalock-Taussig shunt, Glenn procedure, Fontan procedure, biological or mechanical valve replacement, or Rastelli procedure, as we provide postoperative heparin therapy with target APTTs only for patients with postsurgical status) to achieve the target APTTs at the PICU in the NCCHD between May 2014 and February 2018. Among the patients, 30 were excluded because their target APTTs were not reached within 72 h of starting heparin therapy. In addition, we excluded 3 patients whose ages were considered to be outliers (17,20, and 22 years) as compared with the ages of the other patients that ranged from 0 to 4 years. Thus, 59 unique patients were included in the analyses, of whom 8 underwent a Blalock-Taussig shunt; 27, the Glenn procedure; 19, the Fontan procedure; 3, mechanical valve replacement; and 2, the Rastelli procedure. None of the patients with biological valve replacements remained in the study. The patients' demographic details are shown in Table 1. The heparin infusion rates upon achieving the target APTTs for the 59 patients included in this study are shown in Fig. 2. The average time to reach target APTT was 25.5 h (range, 4.1-69.1 h). Of note, none of the patients developed heparin-induced thrombocytopenia, thrombosis or apparent bleeding before achieving target APTT after starting heparin therapy.
None of the 59 patients underwent measurement for factor VIII after starting the heparin treatment, and no platelet transfusions or protamine were administered within 48 h of achieving the target APTT after starting the heparin treatment; thus, no association analyses were performed on these aspects. Table 2 shows the results of the associations between each covariate and heparin infusion rate upon achieving the target APTTs. We additionally assessed the associations of antithrombin III level as a continuous variable, given that this variable is well known www.nature.com/scientificreports/ to be involved in heparin's mechanism of action and was important in the pediatric ECMO study 19,20 , which resulted in finding no significant association (P = 0.32). The infusion rates of heparin drip by type of surgical procedure are shown in Fig. 3. The types of surgical procedure were significantly associated with the differences in heparin infusion rate (P = 0.00073; Fig. 3). Subsequently, the explained variance of these covariates was estimated by summing each explained variance. We constructed the models as shown in Supplemental Table 2. Each model considers the group of variables used for the multiple regression model. After rejecting the regression models where the directions of the regression coefficients are inconsistent with the analyses for each variable, only models 1-7, 9, and 11-13 remained, as shown in Table 3, Fig. 4, and Supplementary Table 3. The maximum variance explained among the models was from model 12, which was 29.1% (Fig. 4). The directions of the associations and AICs for each model are also shown in Table 3 and Supplementary Table 3. We also estimated the variance explained, adjusting for the APTT in each model, which resulted in no apparent changes in the results as compared with those without adjustment for the APTT (Supplementary Table 4).
In total, 56 unique patients were included in the subgroup analyses, excluding three patients with mechanical valve replacement. We assessed the associations between covariates and heparin infusion rate. Age and type of surgical procedure were significantly associated with heparin infusion rate (P = 0.038 and 0.0044, respectively; Supplementary Table 5). Subsequently, the explained variance of the covariates was estimated. The constructed models are shown in Supplementary Table 6. Models 1-4, 7, 9, 11, and 12 remained after rejecting the regression models in which the directions of the regression coefficients are inconsistent. The maximum explained variance was from model 11, which was 27.5%, and the directions of the associations and AICs for each model are also shown in Supplementary Table 7. We also estimated the explained variance, adjusting for the APTT in each model, which resulted in no apparent changes in the results compared with those without adjustment for the APTT (models 1-4, 7-9, 11, and 12 remained after rejecting the regression models where the directions of the regression coefficients were inconsistent; Supplementary Table 8).

Discussion
To the best of our knowledge, this is the first study to identify the extent of each attributing factors to the interindividual variability in unfractionated heparin therapy for patients not on ECMO who were admitted to the PICU after cardiac surgery. Our key finding is that over 70% of the interindividual variability in maintenance unfractionated heparin therapy was unexplained suggesting that responses to heparin therapy are difficult to predict. Factors/parameters from current practice which contributed to heparin dosing included the type of www.nature.com/scientificreports/ surgery in the main analysis, and age and the type of surgery in subgroup analysis. A similar study investigating variability in heparin dose response in the pediatric patients on ECMO was conducted using this approach by Moynihan et al. They showed that less than 50% of the variability was explained using a model incorporating age and antithrombin activity 19,20 . Other studies that assessed interindividual variability of other than heparin dose response were regarding capacity for motor recovery after ischemic stroke or sleeping metabolic rate. The results showed that the clinical variables explained 65-89% of the variance, which is obviously much higher than the rate in our study 21,22 . Thus, further investigation is needed to unveil the unknown factors contributing to the variabilities. Unfractionated heparin interacts with antithrombin III and inhibits activated coagulation factors involved in clotting sequence 11,23,24 . Unfractionated heparin is bound to antithrombin, fibrinogens, globulins, serum proteases, and lipoproteins 11,23,25 . The anticoagulant effect is produced by inactivating thrombin and factor Xa 26 . It is mainly cleared by the liver and reticuloendothelial cells, and clearance starts by binding to proteins, endothelial cells, and macrophages 23 . It is also eliminated via the kidneys; however, some of the elimination pathways are still unknown 27,28 . Alternate explanations for variability in dosing could be explained by differences in pentasaccharide sequence content within the vials 29-31 and the impact of developmental hemostasis [32][33][34] . As many of these factors are unmeasured in clinical practice, patient responses to unfractionated heparin therapy tend to be unpredictable.
In our study, we incorporated numerous factors from our clinical practice likely to contribute to heparin dose variability. We did not observe influences on heparin dosing from multiple factors that could affect the responses to heparin therapy as identified in previous studies [5][6][7][8][9][10][11] . The only known factors that showed any associations with the rate of heparin therapy upon achieving APTT were the type of surgery in the main analysis, and age and the type of surgery in subgroup analysis. Our inability to predict heparin dose requirements represents a clinical challenge.
Finally, we built a model with the maximum number of variables, resulting in the explained variance being < 30% (Table 3, Fig. 4, and Supplementary Table 3). In other words, when calculating the AIC for each Figure 3. Infusion rates of heparin drip by type of surgical procedure. The distribution of each type of surgical procedure was as follows: 1, Blalock-Taussig shunt; 2, Glenn procedure; 3, Fontan procedure; 5, mechanical valve replacement; and 6, Rastelli procedure. The types of surgical procedure were significantly associated with the heparin infusion rates (P = 0.00073 by the Kruskal-Wallis test). Table 3. Coefficients of the models and percentages of explained variance in the dependent variables. Only the models for which coefficients had P values of *< 0.1, **< 0.05, and ***< 0.01 are shown. www.nature.com/scientificreports/ model to assess the fitting into our data, even by the best model with the highest ability to predict, which was Model 12, accounting for Blalock-Taussig shunt plus the Glenn procedure and mechanical valve replacement, the variance explained was only 29.1%. We also ensured that other environmental factors related to the measurements of APTT did not change our results by adjusting for the APTT (Supplementary Table 3). Our results considering the type of surgery such as Blalock-Taussig shunt plus the Glenn procedure and mechanical valve replacement could contribute the most to the variabilities. Our finding that Blalock-Taussig shunt was negatively associated with the variabilities may bring new insights into the behavior of the variabilities after each surgery. A subsequent subgroup analysis was conducted to ensure that these findings were not biased by inclusion of patients with mechanical valve replacement for whom the target APTT was slightly higher than patients who underwent other surgical procedures. We then confirmed that the findings were very similar to the main analyses (Supplementary Tables 7 and 8), with the maximum explained variance being estimated by the best model, which accounted for type of surgical procedure having the lowest AIC. Our findings suggest that type of surgery may explain some variability; however, we need to interpret this with caution. No previous studies have been conducted to investigate pathophysiological mechanisms which affect responses to heparin therapy between the type of cardiac surgeries. Furthermore, these associations are not modifiable. Further investigations are required to better explain our findings, improve heparin dose prediction and ultimately enhance patient care.
Obviously, many more factors should be taken into account to fully understand the interindividual variabilities, an example of which could be genetic factors. Genetic factors have proven to be important modulators of the metabolism of medications and can influence their efficacy and toxicity. Previous studies have reported that 20%-30% of the inter-individual differences in drug metabolism and drug response were estimated to be due to genetic variations 35,36 . Until now, no pharmacogenomic studies related to heparin therapy have been conducted, except for heparin-induced thrombocytopenia [37][38][39] . Including genetic factors could improve the prediction of interindividual variabilities and will be the next step to consider.
Our study has several limitations that must be acknowledged. First, the sample size was small, given the single-center study. Second, aPTT values differ between institutions and assay methods 7,40,41 . Third, the protocols of heparin therapy after cardiac surgery vary depending on the institution, and our aPPT targets and low heparin dose requirements present unique limiting generalizability [42][43][44] . Fourth, while we attempted to identify relevant covariates for our clinical practice, but this was limited by the retrospective nature of the study. Fifth, we only measured the first heparin dose required to achieve the therapeutic range, and changes over time were not evaluated. Finally, the influence of clinician bias in titrating heparin dose at the bedside according to surgery and other patient factors is not accounted for in our study design. Future studies should be conducted to explore this individual variation across other populations and evaluate associations with clinically relevant outcomes, such as bleeding and thrombosis. Variance explained by covariates associated with the interindividual variabilities of heparin infusion rate. We plotted the explained variance using the models. We constructed the models to account for each procedure as follows: model 1 for Blalock-Taussig shunt, model 2 for the Glenn procedure, model 3 for the Fontan procedure, model 4 for mechanical valve replacement, model 5 for Blalock-Taussig shunt plus Glenn procedure, model 6 for Blalock-Taussig shunt plus Fontan procedure, model 7 for Blalock-Taussig shunt plus mechanical valve replacement, model 9 for Glenn procedure plus mechanical valve replacement, model 11 for Blalock-Taussig shunt plus Glenn procedure plus Fontan procedure, model 12 for Blalock-Taussig shunt plus Glenn procedure plus mechanical valve replacement, and model 13 for Blalock-Taussig shunt plus Fontan procedure plus mechanical valve replacement.