Skip to main content

A Bayesian multivariate latent t-regression model for assessing the association between corticosteroid and cranial radiation exposures and cardiometabolic complications in survivors of childhood acute lymphoblastic leukemia: a PETALE study



Childhood acute lymphoblastic leukemia (cALL) is the most frequent pediatric cancer. Over the past decades, treatment of cALL has significantly improved, with cure rates close to 90%. However intensive chemotherapy and cranial radiotherapy (CRT) during a critical period of a child’s development have been shown to lead to significant long-term side effects including cardiometabolic complications. Using the PETALE (Prévenir les effets tardifs des traitements de la leucémie aiguë lymphoblastique chez l’enfant) cALL survivor cohort, we investigated the association between combined cumulative corticosteroids (CS) doses and CRT exposures and obesity, insulin resistance, (pre-)hypertension, and dyslipidemia jointly.


A Bayesian multivariate latent-t model which accounted for our correlated binary outcomes was used for the analyses (n = 241 survivors). CS doses were categorized as low (LD) or high (HD). Combined exposure levels investigated were: 1) LD/no CRT; 2) LD/CRT, and; 3) HD/CRT. We also performed complementary sensitivity analyses for covariate adjustment.


Prevalence of cardiometabolic complications ranged from 12.0% for (pre-)hypertension to 40.2% for dyslipidemia. The fully adjusted odds ratio (OR) for dyslipidemia associated with LD/CRT (vs. LD/No CRT) was OR = 1.98 (95% credible interval (CrI): 1.02 to 3.88). LD/CRT level also led to a 0.15 (95% CrI: 0.00 to 0.29) excess risk to develop at least one cardiometabolic complication. Except for obesity, adjusted results for the highest exposure category HD/CRT were generally similar to those for LD/CRT albeit not statistically significant. White blood cell count at diagnosis, a proxy for cALL burden at diagnosis, was found associated with insulin resistance (OR = 1.08 for a 10-unit increase (× 109/L), 95% CrI: 1.02 to 1.14).


Our results indicated that combined LD/CRT exposure is a likely determinant of dyslipidemia among cALL survivors. No evidence was found to suggest that high doses of CS lead to additional risk for obesity, insulin resistance, (pre-)hypertension, and dyslipidemia beyond that induced by CRT. The multivariate model selected for analyses was judged globally useful to assess potential exposure-related concomitance of binary outcomes.

Peer Review reports


Childhood acute lymphoblastic leukemia (cALL) is the most frequent pediatric cancer. Over the years, treatment of cALL has significantly improved such that nearly 90% of patients are cured with current therapy regimens. This success is largely due to the progressive intensification of cALL treatment and the implementation of risk-adapted protocols [1, 2]. Despite this remarkable improvement, intensive chemotherapy and cranial radiotherapy (CRT) during a vulnerable period of a child’s development can lead to significant long-term side effects. Childhood ALL survivors are known to be at higher risk of developing late treatment-associated morbidities (late adverse effects; LAEs), with a cumulative incidence of chronic health conditions exceeding 60% [3, 4]. The PETALE (Prévenir les effets tardifs des traitements de la leucémie aiguë lymphoblastique chez l’enfant) study, recently conducted at Sainte-Justine University Health Center (SJUHC), was designed to identify and characterize the most common LAEs observed in cALL survivors [5].

Levy et al. [6] investigated the prevalence of cardiometabolic complications in the PETALE cohort. A total of 88 patients (35.6%) presented a single cardiometabolic outcome among obesity, (pre)-hypertension, prediabetes (insulin resistance), and dyslipidemia, while 61 patients (24.7%) cumulated two or more. A large proportion of PETALE cALL survivors received CRT, specifically those at high risk of relapse. In multiple regression analyses, CRT was found to increase the prevalence of dyslipidemia (relative risk (RR): 1.60; 95% confidence interval (CI): 1.07 to 2.41) and high low-density lipoprotein (LDL)-cholesterol (RR: 4.78; 95% CI: 1.72 to 13.28). Moreover, male gender was found positively associated with (pre)-hypertension (RR: 5.12; 95% CI: 1.81 to 14.5). When compared to the general Canadian population, Levy et al. [6] found that cALL survivors were at higher risk of having the metabolic syndrome, dyslipidemia, (pre-)hypertension and high LDL-cholesterol. However, after stratification, only CRT recipients were found more likely than healthy controls to have cardiometabolic complications.

Building upon Levy et al. [6], we performed additional analyses to shed further light on the impact of cALL treatment (exposure) received on cardiometabolic health in the PETALE cohort of survivors. Chemotherapy dosage was not evaluated in Levy et al. [6], however both radiotherapy and chemotherapy treatments may contribute to the development of cardiometabolic complications through several mechanisms, for example by damaging endocrine organs, inducing endothelial dysfunction and perturbing adipose tissue metabolism [6,7,8]. cALL chemotherapy includes multiple agents that are believed to be highly deleterious for cardiometabolic health, such as corticosteroids (CS). We thus investigated the association between doses of CRT and CS received and cardiometabolic complications (obesity, insulin resistance, (pre-)hypertension, and dyslipidemia). Rather than performing analyses for each cardiometabolic outcome separately, our analyses considered all outcomes simultaneously to account for likely co-occurrence between them.

Multivariate statistical models are indicated when analyzing interrelated outcome variables. Grounded in the Bayesian framework, the multivariate logistic regression model of O’Brien and Dunson [9] is an appealing model to study the association between a set of covariates and several correlated binary outcomes. A distinctive feature of this model is its ability to output results with high interpretability, by allowing for the familiar log odds ratio interpretation of regression coefficients as in a standard univariate logistic model, while accounting for outcome dependencies in an unrestricted fashion [9]. Inference based on O’Brien and Dunson’s model was recently further pushed forward by Hund et al. [10] who estimated average treatment effects for types of uranium exposure on kidney disease, hypertension, and diabetes, jointly, in the Navajo Nation. Despite its attractiveness, this Bayesian multivariate model, combined with the post-hoc analysis approach of Hund et al. [10], has not received broad attention in clinical epidemiology, perhaps due to a greater complexity. In addition to the substantive objective described above, a methodological objective was therefore to detail and emphasize the application of this well thought multivariate model on our cALL survivorship data.

In this study, we used the O’Brien and Dunson / Hund et al. global modeling approach to identify combined CS and CRT treatment levels associated with cardiometabolic LAEs in cALL survivors. The insight gained on treatment-LAE associations will further guide risk-adapted follow-up of cALL survivors. From a methodological perspective, our modeling proposal (with associated computer code) will expand health practitioners’ toolkit for the analysis of correlated binary outcome variables from a Bayesian perspective.


PETALE cohort

Subjects enrolled in the PETALE study were treated for cALL at SJUHC with the Dana Farber Cancer Institute (DFCI) protocols 1987–01 to 2005–01 [11, 12]. cALL survivors less than 19 years old at diagnosis, more than 5 years post diagnosis and who had never experienced a relapse were invited to participate. Between 2012 and 2016, a total of 247 participants underwent thorough clinical, biological, and psychosocial evaluations (mean age at interview: 21.6 years; standard deviation (SD): 6.3). A detailed presentation of the PETALE study and cohort is provided in Marcoux et al. [5]

Outcome, treatment and adjustment variables

Obesity was determined by presenting at least one of two factors at interview: obese according to body mass index (BMI) or high waist circumference [13, 14]. Insulin resistance was defined by high blood fasting glucose, high glycated hemoglobin, high homeostasis model assessment (HOMA-IR) [insulin (mIU/L) x glucose (mmol/L)/22.5] values or taking medication for diabetes [15, 16]. (Pre-)hypertension was assessed according to current guidelines for arterial pressures or by being on antihypertensive medication [17, 18]. Dyslipidemia was determined by high LDL-cholesterol, high triglycerides, low high-density lipoprotein (HDL)-cholesterol, or lipid-lowering medication [19, 20]. Age- and gender-dependent cut-off values defining cardiometabolic outcomes are presented in Additional file 1: Table S1.

Cumulative doses of CS were calculated as the sum of prednisone-equivalent, body-surface area-normalized (mg/m2) dexamethasone, methylprednisolone, prednisolone, and prednisone doses received during the three phases of the chemotherapy treatment (induction, intensification, and consolidation). A third quartile split (using a 13,414 mg/m2 threshold) was used to categorize cALL survivors as having received lower (LD) or higher (HD) doses of CS. CRT exposure was defined as a binary variable (yes/no) as most survivors who received CRT had an 18 Gray (Gy) cumulative dose, with none receiving over 20 Gy. Unlike Hund et al. [10] who focused on the effect of a single binary exposure variable entered in the multivariate model each time, a joint CS/CRT exposure variable was defined using three categories: 1) LD/No CRT; 2) LD/CRT, and; 3) HD/CRT. No cALL survivors received higher doses of CS but no CRT; as such we did not consider the “HD/No CRT” category. Information concerning the three exposure categories examined was available for 241 cALL survivors (five survivors were excluded because of no information on CS and one was excluded because they were not exclusively treated by DFCI cALL protocols 1987–01 to 2005–01).

The adjustment variables considered in our analyses were sex, age at diagnosis, white blood cell (WBC) count at diagnosis, and time since diagnosis. These adjustment variables are known to be associated with cALL treatment received, the studied cardiometabolic complications, or both. Notably, age and WBC count were used as criteria to stratify patients as having standard or high risk of relapse in DFCI cALL treatment protocols 1987–01 to 2005–01 [11, 12]. WBC count was entered as a main effect term after preliminary investigations assessing other possible functional forms with one degree of freedom (e.g., log, square root). Indeed, for the sake of simplicity, and because the proposed multivariate model has a marginal logistic interpretation, we initially considered a number of multiple logistic regressions for modeling insulin resistance versus different functions of WBC. We then selected the functional form for WBC that yielded the smallest Akaike information criterion (AIC) value for the corresponding fully adjusted model.

Multivariate model

Let \( {\boldsymbol{Y}}_i=\left({Y}_i^O,{Y}_i^I,{Y}_i^H,{Y}_i^D\right) \) be the vector of binary (1/0) outcomes for the i th survivor, where O, I, H, D indexes obesity, insulin resistance, (pre-)hypertension, and dyslipidemia, respectively. Let T be a p-dimensional vector following a multivariate Student-t distribution with ν degrees of freedom, location vector μ and scale matrix Σ with density

$$ {St}_p\left(\boldsymbol{t}\left|\boldsymbol{\mu}, \sum, v\right.\right)=\frac{\Gamma \left(\left(v+p\right)/2\right)}{\Gamma \left(v/2\right){\left( v\pi \right)}^{p/2}\left|\mathbf{\sum}\left|{}^{1/2}\right.\right.}\times {\left(1+\frac{1}{v}{\left(\boldsymbol{t}-\boldsymbol{\mu} \right)}^{\hbox{'}}{\mathbf{\sum}}^{-1}\left(\boldsymbol{t}-\boldsymbol{\mu} \right)\right)}^{-\frac{p+v}{2}}. $$

We defined \( {\boldsymbol{Z}}_i=\left({Z}_i^O,{Z}_i^I,{Z}_i^H,{Z}_i^D\right) \) as a multivariate Student-t vector of latent (unobserved) variables with mean \( {\boldsymbol{\mu}}_i=\left({\mu}_i^O,{\mu}_i^I,{\mu}_i^H,{\mu}_i^D\right) \):

$$ {\mu}_i^j={\beta}_0^j+{\beta}_1^j{T}_{1,i}+{\beta}_2^j{T}_{2,i}+{\boldsymbol{\beta}}_C^j{\boldsymbol{C}}_i^{\hbox{'}},\kern1em j\in \left\{ 0,I,H,D\right\}, $$

where βj is the vector of unknown regression coefficients associated to outcome j, T1 and T2 are dummies indexing the 2nd and 3rd categories of treatment, respectively, as opposed to the baseline (T0: LD/no CRT), and C is the vector of adjustment covariates. The scale matrix was defined as Σ = δ2R, where R is an unstructured correlation matrix, δ2 = π2(ν − 2)/(3ν), and ν = 7.3. The latent variables can be viewed as underlying propensities to experience the binary cardiometabolic complications; they were linked to the outcomes by \( {Y}_i^j=\mathbb{1}\left({Z}_i^j>0\right) \), where \( \mathbb{1}\left(\bullet \right) \) is the indicator function. Let Xi = (T1,i, T2,i, Ci), then the likelihood function was written as

$$ P\left({\boldsymbol{Y}}_i={\boldsymbol{y}}_i|{\boldsymbol{X}}_i,\boldsymbol{\beta}, \boldsymbol{R}\right)={\int}_{-\infty}^{\infty}\dots {\int}_{-\infty}^{\infty}\left({\prod}_j\mathbb{1}{\left({z}_i^j>0\right)}^{y_i^j}\mathbb{1}{\left({z}_i^j\le 0\right)}^{1-{y}_i^j}\right)\times {St}_4\left({\boldsymbol{z}}_i|{\boldsymbol{\mu}}_i,{\delta}^2\boldsymbol{R},\nu \right)d{\boldsymbol{z}}_i. $$

As detailed in O’Brien and Dunson [9] and further discussed in Hund et al. [10], this latent multivariate Student-t model with δ2 = π2(ν − 2)/(3ν) and ν = 7.3 has a very close resemblance with a latent multivariate logistic model with the same location vector μi [9]. As such, each non-intercept \( {\beta}_{\bullet}^j \) in (1) can be interpreted as a conditional log odds ratio. There is a specific computational advantage of working with the Student-t model rather than the logistic model representation. Indeed, after assigning priors to the parameters β = (βO, βI, βH, βD) and R, a relatively simple data-augmented Markov chain Monte Carlo (MCMC) algorithm for sampling from these parameters’ posterior distribution can be devised. This algorithm samples from the full conditional distributions of β and Zi, and updates R using a Metropolis step.


The prior distribution was specified by π(β, R) = π(β)π(R), where π(R) is uniform on the space of correlation matrices and π(β) is multivariate normal. More precisely, \( \pi \left(\boldsymbol{\beta} \right)={\prod}_jN\left({\boldsymbol{\beta}}^j|{\boldsymbol{\beta}}_{\ast}^j,{\boldsymbol{\Sigma}}_{\ast}^j\right) \), where \( {\boldsymbol{\beta}}_{\ast}^j=\mathbf{0} \) and \( {\boldsymbol{\Sigma}}_{\ast}^j=\mathit{\operatorname{diag}}\left(1000,4,\dots, 4\right) \) for j {O, I, H, D}. Our choice of prior is very similar to a multivariate normal prior distribution proposed by O’Brien and Dunson [9]. For the non-intercept coefficients, a standard deviation of 2 (variance 4) implies that the bulk of the prior is between +/− 1.96*2 for the log odds ratio. We used a larger prior variance (1000 instead of 4) for the intercept coefficient \( {\beta}_0^j \), j {O, I, H, D}, so to only weakly inform on the log odds for outcome j at the baseline level. Crude and adjusted analyses were performed. Because WBC count at diagnosis was not taken into account in Levy et al. [6], the sensitivity of the results to the inclusion of this covariable in the model was assessed.

We used a multivariate normal as the proposal distribution in the Metropolis step for R. The parameters of this distribution were selected towards the goal to have adequate acceptance rates (near 25%). For each of the three models, we obtained 1,050,000 draws from the posterior distribution of (β, R) and discarded the first 50,000 as burn-in. The running time was about 30 min for fitting the most complex model on an Intel core i7 3930k 3.2GHz with 32.0 GB RAM. For each model, we reduced the autocorrelation between the successive draws for the elements of R by applying a thinning with a factor of 100, leaving a total of 10,000 draws available for inference. Chain convergence was assessed visually for every parameter using ten sequential boxplots of 1000 draws each. No convergence problems were further detected according to the Geweke diagnostic test [21]. For the full model, the effective sample size for each of the 28 regression coefficients (β) was 10,000 according to the R package coda function effectiveSize [22]. The effective sample sizes for the 6 unique elements of the correlation matrix R ranged between 4317 and 5325.

Standard exposure association measures were first obtained for each outcome separately for each model considered and every retained posterior draw for (β, R). Specifically, in addition to conditional odds ratios (ORs), conditional risk differences were computed for each survivor under both levels of treatment compared (Tl vs. T0, l = 1, 2) and subsequently averaged to produce population risk differences (PRDs). The multivariate model also enabled the calculation of the difference in probability of a given configuration of the outcome vector under both levels of treatment (Tl vs. T0):

$$ {\displaystyle \begin{array}{c}{P}^{T_l}\left({Y}_i^O={y}_i^O,{Y}_i^I={y}_i^I,{Y}_i^H={y}_i^H,{Y}_i^D={y}_i^D\left|{\boldsymbol{C}}_i\right.\right)\\ {}-{P}^{T_0}\left({Y}_i^O={y}_i^O,{Y}_i^I={y}_i^I,{Y}_i^H={y}_i^H,{Y}_i^D={y}_i^D\left|{\boldsymbol{C}}_i\right.\right),\end{array}} $$

for \( {y}_i^j\in \left\{0,1\right\} \), where, for example,

$$ {P}^{T_1}\left({Y}_i^O=1,{Y}_i^I=0,{Y}_i^H=0,{Y}_i^D=1|{\boldsymbol{C}}_i\right) $$
$$ ={\int}_0^{\infty }{\int}_{-\infty}^0{\int}_{-\infty}^0{\int}_0^{\infty }{St}_4\left({\boldsymbol{z}}_i|{\boldsymbol{\mu}}_i^{T_1},{\delta}^2\boldsymbol{R},\nu\ \right)d{z}_i^Od{z}_i^Id{z}_i^Hd{z}_i^D, $$

with \( {\boldsymbol{\mu}}_i^{T_1}=\left({\mu}_i^{O,{T}_1},{\mu}_i^{I,{T}_1},{\mu}_i^{H,{T}_1},{\mu}_i^{D,{T}_1}\right) \) and \( {\mu}_i^{j,{T}_1}={\beta}_0^j+{\beta}_1^j\cdot 1+{\beta}_2^j\cdot 0+{\boldsymbol{\beta}}_C^j{\boldsymbol{C}}_i^{\hbox{'}} \). For each of the 24 = 16 outcome configurations, (2) was obtained using the R function pmvt which evaluated the latent variables’ multivariate Student-t cumulative distribution function for limits of integration compatible with the configuration. Each difference (2) was then averaged across individuals to obtain a population value, for example,

$$ {P}^{T_1}\left({Y}_i^O=1,{Y}_i^I=0,{Y}_i^H=0,{Y}_i^H=1\right)-{P}^{T_0}\left({Y}_i^O=1,{Y}_i^I=0,{Y}_i^H=0,{Y}_i^H=1\right)= $$
$$ \frac{1}{n}\sum \limits_{i=1}^n\left\{\ {P}^{T_1}\left({Y}_i^O=1,{Y}_i^I=0,{Y}_i^H=0,{Y}_i^D=1|{\boldsymbol{C}}_i\right)-{P}^{T_0}\left({Y}_i^O=1,{Y}_i^I=0,{Y}_i^H=0,{Y}_i^D=1|{\boldsymbol{C}}_i\right)\right\}. $$

As in Hund et al. [10], contrasts in probabilities for the presence of multiple cardiometabolic outcomes (P(N = ∑jYj ≥ 1, 2, 3, 4)) given treatment levels compared were also obtained. For example, the population risk difference for at least three cumulative cardiometabolic outcomes for LD/CRT versus LD/no CRT (\( {P}^{T_1}\left(N\ge 3\right)-{P}^{T_0}\left(N\ge 3\right) \)) was calculated as

$$ \left\{{P}^{T_1}\left({Y}_i^O=1,{Y}_i^I=1,{Y}_i^H=1,{Y}_i^D=1\right)-{P}^{T_0}\left({Y}_i^O=1,{Y}_i^I=1,{Y}_i^H=1,{Y}_i^{\mathrm{D}}=1\right)\right\}+ $$
$$ \left\{{P}^{T_1}\left({Y}_i^O=1,{Y}_i^I=1,{Y}_i^H=1,{Y}_i^D=0\right)-{P}^{T_0}\left({Y}_i^O=1,{Y}_i^I=1,{Y}_i^H=1,{Y}_i^D=0\right)\right\}+ $$
$$ \left\{{P}^{T_1}\left({Y}_i^O=1,{Y}_i^I=1,{Y}_i^H=0,{Y}_i^D=1\right)-{P}^{T_0}\left({Y}_i^O=1,{Y}_i^I=1,{Y}_i^H=0,{Y}_i^D=1\right)\right\}+ $$
$$ \left\{{P}^{T_1}\left({Y}_i^O=1,{Y}_i^I=0,{Y}_i^H=1,{Y}_i^D=1\right)-{P}^{T_0}\left({Y}_i^O=1,{Y}_i^I=0,{Y}_i^H=1,{Y}_i^D=1\right)\right\}+ $$
$$ \left\{{P}^{T_1}\left({Y}_i^O=0,{Y}_i^I=1,{Y}_i^H=1,{Y}_i^D=1\right)-{P}^{T_0}\left({Y}_i^O=0,{Y}_i^I=1,{Y}_i^H=1,{Y}_i^D=1\right)\right\}. $$

We obtained draws from the posterior distributions of target association measures and calculated corresponding posterior means and equal-tailed 95% credible intervals (CrI). Because the posterior distributions of ORs were skewed, posterior medians were presented and favored for this association measure. Due to higher computational complexity, cumulative PRDs were computed on the supercomputer Briarée maintained by Calcul Québec.


Descriptive statistics on the 241 PETALE cALL survivors included in our analyses are presented in Table 1. Males were approximately as represented as females (118/241 or 49.0%), the mean age at cALL diagnosis was 6.6 years (SD = 4.6), and the mean time elapsed between diagnosis and interview was 15.4 years (SD = 5.1). A total of 99 (41.1%) cALL survivors received treatment level LD/No CRT, 82 (34.0%) received LD/CRT, and 60 (24.9%) HD/CRT. The prevalence of the four cardiometabolic complications at interview was relatively high, ranging from 12.0% for (pre-) hypertension to 40.2% for dyslipidemia. On average, cALL survivors who received the baseline level of treatment (LD/No CRT) were the youngest at diagnosis, had the lowest WBC count at diagnosis, and included slightly more females compared to the other two groups of exposure (LD/CRT and HD/CRT). Imbalance in the distribution of time since diagnosis was also observed across levels of treatment, with longer times associated with the highest exposure category (HD/CRT).

Table 1 Characteristics of the cALL survivor cohort

Crude and adjusted results pertaining to treatment are presented individually for each outcome in Tables 2 and 3. The median-based fully adjusted OR for dyslipidemia associated with treatment LD/CRT (as opposed to LD/No CRT) was OR = 1.98 (95% CrI: 1.02 to 3.88) and the corresponding risk difference was PRD = 0.15 (95% CrI: 0.00 to 0.29). The dyslipidemia OR and PDR for HD/CRT (as opposed to LD/No CRT) were similar to those obtained for LD/CRT but they were not statistically significant, likely due to a smaller sample size in that combined exposure category (n = 82 versus n = 60). All other fully adjusted ORs and PRDs were not significantly different from their reference values (one and zero, respectively). Association measures were globally sensitive to the inclusion or exclusion of WBC count at diagnosis in the model, suggesting that it could be a confounder for studied treatment and outcomes. In particular, association measures related to obesity and insulin resistance were attracted towards the null when adjusting for this covariable.

Table 2 Crude and adjusted odds ratios (ORs) associated with treatment for individual cardiometabolic risk factors
Table 3 Crude and adjusted population risk differences (PRDs) associated with treatment for individual cardiometabolic risk factors

Examining the elements of R in the fully adjusted model, we observed a positive residual correlation (i.e., a correlation not induced by treatment and covariates) between obesity and dyslipidemia (ρ = 0.36, 95% CrI: 0.16 to 0.54) and an even stronger one between obesity and insulin resistance (ρ = 0.69, 95% CrI: 0.51 to 0.83). Table 4 presents PRDs for accumulation of cardiometabolic complications. In the fully adjusted model, the treatment level LD/CRT led to a 0.15 increased probability (95% CrI: 0.02 to 0.27) to develop at least one cardiometabolic complication (N ≥ 1), as compared to the baseline treatment level (LD/No CRT). The magnitude of this excess risk was mainly driven by configurations (YO = 0, YI = 0, YH = 0, YD = 1) and (YO = 1, YI = 0, YH = 0, YD = 1) for which the posterior means for averaged differences (2) were 0.04 and 0.05, respectively. The PRDs for N ≥ 1, 2 were significantly different from zero in the adjusted model without WBC count, but these results should be interpreted with caution in light of the potentially confounding effect of WBC. Nonetheless, these results highlight an intrinsic advantage of multivariate modeling as opposed to standard univariate modeling. Indeed, while in that partially adjusted model exposure level LD/CRT was not found associated with any of the outcome marginally (see rightmost column of Tables 2 and 3), it was found associated with the occurrence of at least one or two metabolic complications when globally accounting for these outcomes (see rightmost column of Table 4).

Table 4 Population risk differences (PRDs) for cumulative metabolic outcomes (N = ∑jYj)

The ORs associated with the adjustment variables for the four metabolic outcomes across the two adjusted models considered (full, full without WBC count) are presented in Additional file 2: Table S2. We found that being male decreased the chance of being obese but increased the chance of being hypertensive. Moreover, time since diagnosis was positively associated with dyslipidemia in the models. These associations were also reported in Levy et al. [6], although therein male gender was found protective for obesity in unadjusted analysis only. Finally, WBC count at diagnosis was the only covariable associated with insulin resistance in the fully adjusted marginal model (1) (OR = 1.08 for a 10-unit increase (× 109/L), 95% CrI: 1.02 to 1.14).


Our study aimed at providing additional insights on the treatment-associated risk of developing cardiometabolic complications among survivors of cALL. More precisely, the present analysis of 241 children, adolescents, and young adult cALL survivors of the PETALE cohort [5] was designed to inform on the association between combined exposures to CRT and CS doses and cardiometabolic LAEs. Moreover, since obesity, insulin resistance, hypertension, and dyslipidemia form a cohesive group of outcomes as part of the metabolic syndrome, the analysis was done using a multivariate modeling strategy enabling both marginal and joint risk assessments. Indeed, for each set of adjustment covariates, we assessed the impact of exposure on the outcomes marginally, jointly and cumulatively using a single model along with corresponding estimated regression parameters. A multivariate approach such as ours provides a more thorough description of which outcomes are affected by a change in the level of exposure, better controls for Type I error rate, and generally increases statistical power as compared to a series of univariate analyses [23, 24].

Although modern cALL treatment protocols are phasing out CRT due to severe toxicities [25], a large fraction of cALL survivors have been exposed to CRT. In our study, we found that cALL survivors who had received a lower cumulative dose of CS with CRT were at increased risk of presenting dyslipidemia: there were 15 excess cases of dyslipidemia per 100 survivors in the LD/CRT group, as compared to the LD/No CRT group. For survivors who had received LD/CRT, we similarly found 15 additional cases presenting one or more cardiometabolic complications per 100 compared to the baseline treatment group. Levy et al. [6] found, using the same cohort of cALL survivors, that CRT recipients had a 60% increased risk to develop dyslipidemia. Our results are thus consistent with those, although our risk difference assessment arguably provides a more clinically useful measure of treatment-related burden.

CRT has been proposed as a contributing factor to the development of cardiometabolic complications and a higher prevalence has been reported when chemotherapy was combined with radiotherapy [26]. As a single treatment exposure, CRT has previously been found associated with dyslipidemia in different cohorts of cALL survivors [8, 27, 28]. van Waas et al. observed that cALL survivors who received CRT had higher total cholesterol levels compared with those who did not, whereas their HDL-cholesterol levels did not differ [28]. Nottage et al. revealed that CRT exposure was associated with an increased risk for abnormal values of LDL-cholesterol, triglycerides, and HDL-cholesterol using St. Jude Lifetime (SJLIFE) cALL survivors treated between 1962 and 2002 and at least 10 years from diagnosis as basis [8]. The impact of CS doses on dyslipidemia development has been less studied. Nottage et al. [8] found no clinically meaningful associations between chemotherapeutic agents, including cumulative prescribed prednisone-equivalent doses, and abnormal lipid levels or other cardiometabolic complications. In our study, the OR and PRD for dyslipidemia associated with HD/CRT treatment (versus LD/No CRT) were of similar magnitude, although not significant, as those found for LD/CRT. Taken together, reported results suggest that CRT, but not CS dose, may be the key contributing factor in the development of dyslipidemia. Notwithstanding this insight, the negative impact of high doses of CS on lipid profile has been described. Co-administration of asparaginase and steroids was shown to cause significant changes in serum lipid levels during cALL treatment [29]. A similar effect due to high-dose CS was observed after organ transplantation [30]. However, these adverse events are known to be transient and, in light of both PETALE and SJLIFE studies, results so far do not support that these acute effects translate into increased dyslipidemia risk in the long-term.

While treatment exposure was not found associated with obesity in our study, CS have been implicated in modifications of the physiology of adiposity. Significant weight gain was shown during the 2 year period of maintenance chemotherapy, when patients are exposed monthly to pulses of high-dose CS [31, 32]. Decreased total lean body mass and increased percentage of fat mass have also been strongly associated with treatment combinations that included platinum, CRT, and/or steroids [7]. Some groups of authors demonstrated that glucocorticoids, including dexamethasone, increase BMI and serum leptin levels [33, 34]. However, among survivors of cALL, dexamethasone was associated with greater BMI during therapy, but was not a contributor at follow-up [35].

Chronic low-grade inflammation is believed to be a major contributing factor for insulin resistance development [36]. ALL and its intensive treatments constitute a fertile ground for the generation of oxidative stress and inflammation. Most chemotherapy agents, particularly doxorubicin, induce pro-inflammatory reactions [37] and have pro-inflammatory effects on the adipose tissue [38]. Furthermore, literature on the impact of oxidative stress on leukemogenesis and anticancer therapies is emerging [39,40,41,42]. Interestingly, our analyses suggested that CS and CRT may not be the primary determinant in the high prevalence of insulin resistance among cALL survivors. Rather, disease burden at diagnosis, as represented by WBC count at diagnosis, might be implicated in the physiopathology of the insulin resistance development. Activation of the immune system and inflammation processes may be detected by an increase in a number of markers, including WBC. Accordingly, a link between elevated WBC count and type 2 diabetes has been found in different populations [43, 44] and in a meta-analysis based on 20 studies [45]. This meta-analysis performed in 2010 also showed that higher levels of granulocytes and lymphocytes, but not monocytes, were associated with incident type 2 diabetes. On their part, Chee-Tin et al. found that all WBC subtypes were independently associated with insulin resistance [46]. Using data from the National Health and Nutrition Examination Survey Epidemiologic Follow-up Study (NHEFS), Ford reported an increasing dose-response relationship between leukocyte concentrations and diabetes incidence (hazard ratio for concentration stratum 9.1-56 × 109/L versus stratum 2.1–5.7 × 109/L: 1.50; 95% CI: 1.12, 1.99) [47].

In our cohort, which features a large window of WBC count outside the normal (4-11 × 109/L) range (min: 0.42 × 109/L; max: 361 × 109/L; median: 9.60 × 109/L; interquartile range: 25.6 × 109/L), an increasing dose-response relationship between WBC count at diagnosis and insulin resistance was also observed. A hypothesis for that finding could be that high titers of WBC cause widespread and severe tissue infiltration, thus triggering a chronic inflammation feedback loop. This mechanism could be similar to that of type 2 diabetes, with a self-sustaining low-grade visceral inflammation loop (reviewed in Zand et al.) [48]. To our knowledge, our study is the first to have investigated and outlined a possible association between cALL burden at diagnosis through WBC count at diagnosis and cardiometabolic LAEs in cALL survivors including insulin resistance. Although further studies are needed to replicate this finding in other cALL survivor cohorts and investigate the physiological mechanisms underlying this phenomenon, WBC count at diagnosis appears to be a predictive factor of insulin resistance in cALL survivors.


In this work, we took an innovative methodological approach based on a Bayesian multivariate latent-t model for evaluating the impact of cALL treatment on the cardiometabolic health of cALL survivors. Our sophisticated analysis of the PETALE cohort data indicated that combined cranial radiotherapy and low doses of corticosteroids were a determinant of dyslipidemia among cALL survivors and were associated with presenting one or more cardiometabolic LAEs. However, we could not similarly conclude that combined exposure to cranial radiotherapy and high doses of corticosteroids were associated with the accumulation of metabolic complications. These inconclusive results should not be interpreted as a lack of global impact of treatment as many aspects of the data may have obscured findings, such as statistical power and confounding effects due to other chemotherapy agents. Notwithstanding this, our multivariate model and its post-hoc analysis were found valuable, although computationally more intensive than standard statistical approaches based on univariate logistic or log-binomial models. In particular, the use of the R function pmvt which evaluated the latent variables’ multivariate Student-t cumulative distribution function was found to be a bottleneck for the computation of the cumulative population risk differences. How to decrease the running time for these calculations would need to be investigated. One possible mitigation could be to evaluate a cumulative risk difference for a typical individual of the population (based on mean or median covariate values, for example) instead of averaging cumulative risk differences over all individuals.

Two frequentist approaches often used for modeling correlated binary data are the logistic generalized estimation equation (GEE) approach [49] and the mixed effects logistic approach, the latter often referred to as a generalized linear mixed model (GLMM) [50, 51]. Because the GEE approach does not completely specify the joint distribution of the outcomes, it would not be possible to calculate multivariate outcome probabilities for obtaining cumulative population risk differences such as was done with proposed multivariate model. While a likelihood is specified in a logistic mixed effects model, the interpretation of the regression coefficients in this model is more difficult than in ours since it is conditional on the random effect(s), and integrating out these effects would imply a loss in the logistic structure. The multivariate logistic regression approach that we considered combines a full likelihood approach with ease of interpretation of the regression coefficients by directly enabling a marginal log odds ratio interpretation. The GEE approach relies on large sample arguments which are not necessary satisfied when dealing with samples of small or moderate sizes (such as the PETALE cohort). In comparison, the proposed multivariate logistic regression approach is exact, notwithstanding the number of MCMC iterations which can be increased with relative ease. In contrast to GLMM, this Bayesian model also yields, under mild conditions that can be easily verified, a proper posterior distribution when a noninformative prior is used [9].

Exemplified by the use of an unstructured correlation matrix R in our analysis, Agresti [52] highlighted the flexibility of the multivariate logistic distribution underlying the proposed model. Copula approaches are other alternatives to flexibly model correlated binary outcomes. e.g. [53,54,55]. For example, Genest et al. [53] introduced a multivariate approach with logistic regression margins combined to a meta-elliptical copula (e.g., normal or Student-t) to model residual dependency. Similarly to our studied model, copula approaches on discrete data are prone to estimation challenges. For their model, Genest et al. [53] proposed a composite likelihood method to alleviate numerical problems for the estimation of the copula (correlation) parameters; an explicit estimator for the asymptotic variance of the correlation parameter’s estimator via linearization was further provided. Overall, we consider the Bayesian multivariate latent-t model investigated as a principled and convenient tool for obtaining a thorough assessment of the possible relationship between an exposure (possibly multicategorical) and a set of correlated binary outcomes.



(childhood) acute lymphoblastic leukemia


Akaike information criterion


Body mass index


Confidence interval


Credible interval


Cranial radiotherapy




Dana Farber Cancer Institute


Generalized estimation equation


Generalized linear mixed model




High dose


High-density lipoprotein


Homeostasis model assessment of insulin resistance


Late adverse effect


Low dose


Low-density lipoprotein


Markov chain Monte Carlo


National Health and Nutrition Examination Survey Epidemiologic Follow-up Study


Odds ratio


Prévenir les effets tardifs des traitements de la leucémie aiguë lymphoblastique chez l’enfant


Population risk difference


Relative risk


Standard deviation


St. Jude Lifetime


Sainte-Justine University Health Center


White blood cell


  1. Pui CH, Robison LL, Look AT. Acute lymphoblastic leukaemia. Lancet. 2008;371(9617):1030–43.

    Article  CAS  PubMed  Google Scholar 

  2. Pui CH, Yang JJ, Hunger SP, Pieters R, Schrappe M, Biondi A, Vora A, Baruchel A, Silverman LB, Schmiegelow K, Escherich G, Horibe K, Benoit YC, Izraeli S, Yeoh AE, Liang DC, Downing JR, Evans WE, Relling MV, Mullighan CG. Childhood acute lymphoblastic leukemia: Progress through collaboration. J Clin Oncol. 2015;33(27):2938–48.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Haddy TB, Mosher RB, Reaman GH. Late effects in long-term survivors after treatment for childhood acute leukemia. Clin Pediatr (Phila). 2009;48(6):601–8.

    Article  Google Scholar 

  4. Mody R, Li S, Dover DC, Sallan S, Leisenring W, Oeffinger KC, Yasui Y, Robison LL, Neglia JP. Twenty-five-year follow-up among survivors of childhood acute lymphoblastic leukemia: a report from the childhood Cancer survivor study. Blood. 2008;111(12):5515–23.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Marcoux S, Drouin S, Laverdiere C, Alos N, Andelfinger GU, Bertout L, Curnier D, Friedrich MG, Kritikou EA, Lefebvre G, Levy E, Lippe S, Marcil V, Raboisson MJ, Rauch F, Robaey P, Samoilenko M, Seguin C, Sultan S, Krajinovic M, Sinnett D. The PETALE study: late adverse effects and biomarkers in childhood acute lymphoblastic leukemia survivors. Pediatr Blood Cancer. 2017;64(6). Epub 2016 Dec 4

    Article  Google Scholar 

  6. Levy E, Samoilenko M, Morel S, England J, Amre D, Bertout L, Drouin S, Laverdiere C, Krajinovic M, Sinnett D, Lefebvre G, Marcil V. Cardiometabolic risk factors in childhood, adolescent and Young adult survivors of acute lymphoblastic leukemia - a Petale cohort. Sci Rep. 2017;7(1):17684.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  7. Baker KS, Chow EJ, Goodman PJ, Leisenring WM, Dietz AC, Perkins JL, Chow L, Sinaiko A, Moran A, Petryk A, Steinberger J. Impact of treatment exposures on cardiovascular risk and insulin resistance in childhood cancer survivors. Cancer Epidemiol Biomark Prev. 2013;22(11):1954–63.

    Article  CAS  Google Scholar 

  8. Nottage KA, Ness KK, Li C, Srivastava D, Robison LL, Hudson MM. Metabolic syndrome and cardiovascular risk among long-term survivors of acute lymphoblastic leukaemia - from the St. Jude lifetime cohort. Br J Haematol. 2014;165(3):364–74.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. O'Brien SM, Dunson DB. Bayesian multivariate logistic regression. Biometrics. 2004;60(3):739–46.

    Article  PubMed  Google Scholar 

  10. Hund L, Bedrick EJ, Miller C, Huerta G, Nez T, Ramone S, Shuey C, Cajero M, Lewis J. A Bayesian framework for estimating disease risk due to exposure to uranium mine and mill waste on the Navajo nation. J. R. Stat. Soc. A. Stat. Soc. 2015;178(4):1069–91.

    Article  Google Scholar 

  11. Silverman LB, Stevenson KE, O'Brien JE, Asselin BL, Barr RD, Clavell L, Cole PD, Kelly KM, Laverdiere C, Michon B, Schorin MA, Schwartz CL, O'Holleran EW, Neuberg DS, Cohen HJ, Sallan SE. Long-term results of Dana-Farber Cancer Institute ALL consortium protocols for children with newly diagnosed acute lymphoblastic leukemia (1985-2000). Leukemia. 2010;24(2):320–34.

    Article  CAS  PubMed  Google Scholar 

  12. Vrooman LM, Neuberg DS, Stevenson KE, Asselin BL, Athale UH, Clavell L, Cole PD, Kelly KM, Larsen EC, Laverdiere C, Michon B, Schorin M, Schwartz CL, Cohen HJ, Lipshultz SE, Silverman LB, Sallan SE. The low incidence of secondary acute myelogenous leukaemia in children and adolescents treated with dexrazoxane for acute lymphoblastic leukaemia: a report from the Dana-Farber Cancer Institute ALL consortium. Eur J Cancer. 2011;47(9):1373–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Katzmarzyk PT. Waist circumference percentiles for Canadian youth 11-18y of age. Eur J Clin Nutr. 2004;58(7):1011–5.

    Article  CAS  PubMed  Google Scholar 

  14. Van den Broeck J, Willie D, Younger N. The World Health Organization child growth standards: expected implications for clinical and epidemiological research. Eur J Pediatr. 2009;168(2):247–51.

    Article  PubMed  Google Scholar 

  15. Allard P, Delvin EE, Paradis G, Hanley JA, O'Loughlin J, Lavallee C, Levy E, Lambert M. Distribution of fasting plasma insulin, free fatty acids, and glucose concentrations and of homeostasis model assessment of insulin resistance in a representative sample of Quebec children and adolescents. Clin Chem. 2003;49(4):644–9.

    Article  CAS  PubMed  Google Scholar 

  16. Chen J, Wildman RP, Hamm LL, Muntner P, Reynolds K, Whelton PK, He J, Third National H, nutrition examination S. Association between inflammation and insulin resistance in U.S. nondiabetic adults: results from the third National Health and nutrition examination survey. Diabetes Care. 2004;27(12):2960–5.

    Article  PubMed  Google Scholar 

  17. Paradis G, Tremblay MS, Janssen I, Chiolero A, Bushnik T. Blood pressure in Canadian children and adolescents. Health Rep. 2010;21(2):15–22.

    PubMed  Google Scholar 

  18. Stern RH. The new hypertension guidelines. J Clin Hypertens (Greenwich). 2013;15(10):748–51.

    Google Scholar 

  19. Expert Panel on. Integrated guidelines for cardiovascular H, risk reduction in C, adolescents, National Heart L, blood I. expert panel on integrated guidelines for cardiovascular health and risk reduction in children and adolescents: summary report. Pediatrics. 2011;128(Suppl 5):S213–56.

    Google Scholar 

  20. Genest J, McPherson R, Frohlich J, Anderson T, Campbell N, Carpentier A, Couture P, Dufour R, Fodor G, Francis GA, Grover S, Gupta M, Hegele RA, Lau DC, Leiter L, Lewis GF, Lonn E, Mancini GB, Ng D, Pearson GJ, Sniderman A, Stone JA, Ur E. 2009 Canadian cardiovascular society/Canadian guidelines for the diagnosis and treatment of dyslipidemia and prevention of cardiovascular disease in the adult - 2009 recommendations. Can J Cardiol. 2009;25(10):567–79.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Geweke J. Evaluating the accuracy of sampling-based approaches to calculating posterior moments. In: Bernardo JM, Berger J, Dawid AP, Smith JFM, editors. Bayesian statistics 4: Oxford University Press; 1992. p. 169–93.

  22. Plummer M, Best N, Cowles K, Vines K. CODA: convergence diagnosis and output analysis for MCMC. R news. 2006;6(1):7–11.

    Google Scholar 

  23. Tabachnick BG, Fidell LS. Using Multivariate Statistics Pearson Education, 2013.

  24. Yang Q, Wang Y. Methods for analyzing multivariate phenotypes in genetic association studies. J Probab Stat. 2012;2012:652569.

    Article  PubMed  PubMed Central  Google Scholar 

  25. Carroll WL, Hunger SP. Therapies on the horizon for childhood acute lymphoblastic leukemia. Curr Opin Pediatr. 2016;28(1):12–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Trimis G, Moschovi M, Papassotiriou I, Chrousos G, Tzortzatou-Stathopoulou F. Early indicators of dysmetabolic syndrome in young survivors of acute lymphoblastic leukemia in childhood as a target for preventing disease. J Pediatr Hematol Oncol. 2007;29(5):309–14.

    Article  PubMed  Google Scholar 

  27. Siviero-Miachon AA, Spinola-Castro AM, Lee ML, Andreoni S, Geloneze B, Lederman H, Guerra-Junior G. Cranial radiotherapy predisposes to abdominal adiposity in survivors of childhood acute lymphocytic leukemia. Radiat Oncol. 2013;8:39.

    Article  PubMed  PubMed Central  Google Scholar 

  28. van Waas M, Neggers SJ, Pieters R, Van den Heuvel-Eibrink MM. Components of the metabolic syndrome in 500 adult long-term survivors of childhood cancer. Ann Oncol. 2010;21(5):1121–6.

    Article  PubMed  Google Scholar 

  29. Parsons SK, Skapek SX, Neufeld EJ, Kuhlman C, Young ML, Donnelly M, Brunzell JD, Otvos JD, Sallan SE, Rifai N. Asparaginase-associated lipid abnormalities in children with acute lymphoblastic leukemia. Blood. 1997;89(6):1886–95.

    Article  CAS  PubMed  Google Scholar 

  30. Mathis AS, Dave N, Knipp GT, Friedman GS. Drug-related dyslipidemia after renal transplantation. Am J Health Syst Pharm. 2004;61(6):565–85 quiz 586-7.

    Article  PubMed  Google Scholar 

  31. Esbenshade AJ, Simmons JH, Koyama T, Koehler E, Whitlock JA, Friedman DL. Body mass index and blood pressure changes over the course of treatment of pediatric acute lymphoblastic leukemia. Pediatr Blood Cancer. 2011;56(3):372–8.

    Article  PubMed  PubMed Central  Google Scholar 

  32. Withycombe JS, Post-White JE, Meza JL, Hawks RG, Smith LM, Sacks N, Seibel NL. Weight patterns in children with higher risk ALL: a report from the Children's oncology group (COG) for CCG 1961. Pediatr Blood Cancer. 2009;53(7):1249–54.

    Article  PubMed  PubMed Central  Google Scholar 

  33. Kiess W, Englaro P, Hanitsch S, Rascher W, Attanasio A, Blum WF. High leptin concentrations in serum of very obese children are further stimulated by dexamethasone. Horm Metab Res. 1996;28(12):708–10.

    Article  CAS  PubMed  Google Scholar 

  34. Wallace AM, Tucker P, Williams DM, Hughes IA, Ahmed SF. Short-term effects of prednisolone and dexamethasone on circulating concentrations of leptin and sex hormone-binding globulin in children being treated for acute lymphoblastic leukaemia. Clin Endocrinol. 2003;58(6):770–6.

    Article  CAS  Google Scholar 

  35. Lindemulder SJ, Stork LC, Bostrom B, Lu X, Devidas M, Hunger S, Neglia JP, Kadan-Lottick NS. Survivors of standard risk acute lymphoblastic leukemia do not have increased risk for overweight and obesity compared to non-cancer peers: a report from the Children's oncology group. Pediatr Blood Cancer. 2015;62(6):1035–41.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Olefsky JM, Glass CK. Macrophages, inflammation, and insulin resistance. Annu Rev Physiol. 2010;72:219–46.

    Article  CAS  PubMed  Google Scholar 

  37. Vyas D, Laput G, Vyas AK. Chemotherapy-enhanced inflammation may lead to the failure of therapy and metastasis. Onco Targets Ther. 2014;7:1015–23.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  38. Sauter KA, Wood LJ, Wong J, Iordanov M, Magun BE. Doxorubicin and daunorubicin induce processing and release of interleukin-1beta through activation of the NLRP3 inflammasome. Cancer Biol Ther. 2011;11(12):1008–16.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  39. Chen YF, Liu H, Luo XJ, Zhao Z, Zou ZY, Li J, Lin XJ, Liang Y. The roles of reactive oxygen species (ROS) and autophagy in the survival and death of leukemia cells. Crit Rev Oncol Hematol. 2017;112:21–30.

    Article  PubMed  Google Scholar 

  40. Udensi UK, Tchounwou PB. Dual effect of oxidative stress on leukemia cancer induction and treatment. J Exp Clin Cancer Res. 2014;33:106.

    Article  PubMed  PubMed Central  Google Scholar 

  41. Mazor D, Abucoider A, Meyerstein N, Kapelushnik J. Antioxidant status in pediatric acute lymphocytic leukemia (ALL) and solid tumors: the impact of oxidative stress. Pediatr Blood Cancer. 2008;51(5):613–5.

    Article  PubMed  Google Scholar 

  42. Manoharan S, Kolanjiappan K, Kayalvizhi M. Enhanced lipid peroxidation and impaired enzymic antioxidant activities in the erythrocytes of patients with cervical carcinoma. Cell Mol Biol Lett. 2004;9(4A):699–707.

    PubMed  Google Scholar 

  43. Vozarova B, Weyer C, Lindsay RS, Pratley RE, Bogardus C, Tataranni PA. High white blood cell count is associated with a worsening of insulin sensitivity and predicts the development of type 2 diabetes. Diabetes. 2002;51(2):455–61.

    Article  CAS  PubMed  Google Scholar 

  44. Yoshimura A, Ohnishi S, Orito C, Kawahara Y, Takasaki H, Takeda H, Sakamoto N, Hashino S. Association of peripheral total and differential leukocyte counts with obesity-related complications in young adults. Obes Facts. 2015;8(1):1–16.

    Article  PubMed  PubMed Central  Google Scholar 

  45. Gkrania-Klotsas E, Ye Z, Cooper AJ, Sharp SJ, Luben R, Biggs ML, Chen LK, Gokulakrishnan K, Hanefeld M, Ingelsson E, Lai WA, Lin SY, Lind L, Lohsoonthorn V, Mohan V, Muscari A, Nilsson G, Ohrvik J, Chao Qiang J, Jenny NS, Tamakoshi K, Temelkova-Kurktschiev T, Wang YY, Yajnik CS, Zoli M, Khaw KT, Forouhi NG, Wareham NJ, Langenberg C. Differential white blood cell count and type 2 diabetes: systematic review and meta-analysis of cross-sectional and prospective studies. PLoS One. 2010;5(10):e13405.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  46. Lee CT, Harris SB, Retnakaran R, Gerstein HC, Perkins BA, Zinman B, Hanley AJ. White blood cell subtypes, insulin resistance and beta-cell dysfunction in high-risk individuals--the PROMISE cohort. Clin Endocrinol. 2014;81(4):536–41.

    Article  CAS  Google Scholar 

  47. Ford ES. Leukocyte count, erythrocyte sedimentation rate, and diabetes incidence in a national sample of US adults. Am J Epidemiol. 2002;155(1):57–64.

    Article  PubMed  Google Scholar 

  48. Zand H, Morshedzadeh N, Naghashian F. Signaling pathways linking inflammation to insulin resistance. Diabetes Metab Syndr. 2017.

  49. Zeger SL, Liang K-Y, Albert PS. Models for longitudinal data: a generalized estimating equation approach. Biometrics. 1988:1049–60.

    Article  CAS  PubMed  Google Scholar 

  50. Stiratelli R, Laird N, Ware JH. Random-effects models for serial observations with binary response. Biometrics. 1984;40(4):961–71.

    Article  CAS  PubMed  Google Scholar 

  51. Breslow NE, Clayton DG. Approximate inference in generalized linear mixed models. J Am Stat Assoc. 1993;88(421):9–25.

    Google Scholar 

  52. Agresti A. Categorical data analysis. 3rd ed. Hoboken: John Wiley & Sons Inc; 2013.

    Google Scholar 

  53. Genest C, Nikoloulopoulos AK, Rivest L-P, Fortin M. Predicting dependent binary outcomes through logistic regressions and meta-elliptical copulas. Brazilian Journal of Probability and Statistics. 2013;27(3):265–84.

    Article  Google Scholar 

  54. Swihart BJ, Caffo BS, Crainiceanu CM. A unifying framework for marginalised random-intercept models of correlated binary outcomes. Int Stat Rev. 2014;82(2):275–95.

    Article  PubMed  PubMed Central  Google Scholar 

  55. Yang L, Frees EW, Zhang Z. Nonparametric estimation of copula regression models with discrete outcomes. J Am Stat Assoc. 2018(just-accepted):1–36.

Download references


The team would like to deeply thank the cALL survivors that graciously accepted and gave of their time to participate in the PETALE project. Parts of statistical analyses were made on the supercomputer Briarée from Université de Montréal, managed by Calcul Québec and Compute Canada. The operation of this supercomputer is funded by the Canada Foundation for Innovation (CFI), the ministère de l’Économie, de la science et de l’innovation du Québec (MESI) and the Fonds de recherche du Québec - Nature et technologies (FRQ-NT).


This work was supported by the Institute of Cancer Research (ICR) of the Canadian Institutes of Health Research (CIHR), in collaboration with C17 Council, Canadian Cancer Society (CCS), Cancer Research Society (CRS), Garron Family Cancer Centre at the Hospital for Sick Children, Ontario Institute for Cancer Research (OICR) and Pediatric Oncology Group of Ontario (POGO) (grant number: TCF 118694) and the Fonds de recherche Québec-Santé. Geneviève Lefebvre and Valérie Marcil are Research Scholars of the Fonds de recherche Québec-Santé. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Availability of data and materials

The dataset analyzed for this study is not publicly available due to confidentiality reasons, but is available from the corresponding author on reasonable request. The R/C++ code used to generate the results can be found in the Journal repository (Additional file 3: MCMC sampler file; Additional file 4: data generation and call to sampler file; Additional file 5: risk difference post-processing file).

Author information

Authors and Affiliations



E.L., M.K., C.L., V.M., and D.S. conceived the PETALE study. S.D., S.M., E.L. and V.M. collected the data and participated in coordination. M.C.F., S.D., M.S., V.M. and G.L. participated to the analysis plan. M.C.F. wrote the computing code (MCMC sampler) and M.C.F. and G.L. performed the analyses. G.L. and V.M. performed the literature review and G.L. wrote the manuscript. All authors have read, edited, and approved the final version of the manuscript.

Corresponding author

Correspondence to Geneviève Lefebvre.

Ethics declarations

Ethics approval and consent to participate

The study was approved by the Institutional Review Board of SJUHC and investigations were carried out in accordance with the principles of the Declaration of Helsinki. Written informed consent was obtained from study participants or parents/guardians.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1:

Table S1. Cut-off values for cardiometabolic outcomes. (DOCX 25 kb)

Additional file 2:

Table S2. Adjusted odds ratios (ORs) for individual cardiometabolic outcomes. (DOCX 18 kb)

Additional file 3:

Supplementary Material: Bayesian multivariate latent-t model MCMC sampler. Rcpp code implementing the MCMC algorithm for sampling from the posterior distribution of the model parameters (ZIP 3 kb)

Additional file 4:

Supplementary Material: Bayesian multivariate latent-t model implementation, R code example to simulate data and to fit the studied multivariate model by calling the MCMC algorithm function. (ZIP 2 kb)

Additional file 5:

Supplementary Material: population risk differences calculation. R code to calculate marginal and cumulative population risk differences using the output of the MCMC sampler. (ZIP 2 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Caubet Fernandez, M., Drouin, S., Samoilenko, M. et al. A Bayesian multivariate latent t-regression model for assessing the association between corticosteroid and cranial radiation exposures and cardiometabolic complications in survivors of childhood acute lymphoblastic leukemia: a PETALE study. BMC Med Res Methodol 19, 100 (2019).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: