Skip to main content

Studying the association between longitudinal mammographic density measurements and breast cancer risk: a joint modelling approach



Researchers have suggested that longitudinal trajectories of mammographic breast density (MD) can be used to understand changes in breast cancer (BC) risk over a woman’s lifetime. Some have suggested, based on biological arguments, that the cumulative trajectory of MD encapsulates the risk of BC across time. Others have tried to connect changes in MD to the risk of BC.


To summarize the MD–BC association, we jointly model longitudinal trajectories of MD and time to diagnosis using data from a large (\(N = 40{,}087\)) mammography cohort of Swedish women aged 40–80 years. Five hundred eighteen women were diagnosed with BC during follow-up. We fitted three joint models (JMs) with different association structures; Cumulative, current value and slope, and current value association structures.


All models showed evidence of an association between MD trajectory and BC risk (\(P < 0.001\) for current value of MD, \(P < 0.001\) and \(P =0.005\) for current value and slope of MD respectively, and \(P < 0.001\) for cumulative value of MD). Models with cumulative association structure and with current value and slope association structure had better goodness of fit than a model based only on current value. The JM with current value and slope structure suggested that a decrease in MD may be associated with an increased (instantaneous) BC risk. It is possible that this is because of increased screening sensitivity rather than being related to biology.


We argue that a JM with a cumulative association structure may be the most appropriate/biologically relevant model in this context.


A woman’s breast is a complex soft-tissue organ and its composition changes over time. Mammographic density (MD) reflects breast tissue composition; epithelial and stromal cells, collagen, and fat and varies extensively between individuals. MD is the most established image-based risk factor for BC [1,2,3,4]. MD also plays a major role in decreasing mammographic sensitivity [5,6,7], since both dense tissue and tumors appear white on mammograms.

Researchers have suggested that longitudinal trajectories of breast density can be used to understand changes in BC risk over a woman’s lifetime [8]. Several studies have tried to understand the changes of density with age using tumor free mammograms [9,10,11,12]. It has been shown that decline in MD with increasing age is most pronounced over the menopausal transition [13], and reaches a plateau around the age of 65 [9]. These trends are consistent across diverse groups of women around the globe, suggesting that they result from an intrinsic biological, most likely hormonal, mechanism [9, 10]. This pattern in density decline led researchers to draw parallels to Pike’s model for the rate of breast tissue ageing [14, 15]. This model suggests that the slowing rate of increase of age-specific BC incidence which is seen after menopause occurs because of a reduction in the rate of breast tissue ageing in postmenopausal women [16]. Boyd et al. [17] showed that age-specific cumulative MD is strongly associated with age-specific incidence of BC using cross-sectional data. They suggested that the cumulative MD profile is consistent with the accumulation of mutations with increasing time of exposure to MD, and that it may be considered as a tissue-specific marker of the biological processes underlying the rate of breast tissue ageing [9, 13]. No study in the literature has yet showed that cumulative value of MD at an individual level, i.e. using longitudinal data, is associated with BC risk.

Others studies have tried to connect changes in MD to the risk of BC via (association) analyses of longitudinal data. These have provided contrasting conclusions [8, 18,19,20,21,22,23]. These studies have used different methods for measuring MD, and have defined density change in different ways. They have used fairly standard statistical methods in their analyses, treating density measurements as fixed covariates. Recent studies [21, 22] have included density change as a time-varying covariate in a time-to-event analysis (Cox model). This approach is not ideal. Firstly, it assumes that the time-dependent covariate is exogenous, and secondly, it compensates for the fact that the covariate is not observed continuously by replacing unobserved values by the most recent observed value. This approach can provide biased parameter estimates, especially when there are irregular and infrequent measures of the marker. Moreover, it assumes that MD is measured without error. MD measurements are known to be noisy. Ignoring measurement error in a marker can lead to a severe underestimation of an association with a time-to-event/diagnosis process [24].

For observational studies with longitudinal data on MD and follow-up of cancer diagnoses it is appropriate to jointly model both processes, i.e. to fit a joint model (JM) [25]. In addition to accounting for the measurement error in the longitudinal outcome (MD), JM offers several options for modeling the association between the marker (MD) and the time-to-event (diagnosis) process [26]. JMs have been previously used to study the association between density/density change and death due to BC in BC patients [27]. As far as we are aware, there is only one publication [28] in which a JM has been used to study the association between MD and BC diagnosis/incidence. This work did not specifically model error in the measurements of density, and assumed a linear trajectory for the longitudinal marker and a simple model for the time-to-event process (current value association structure; see “Methods” section).

The objective of this paper is to explore the use of joint modelling of observational data, for investigating hypotheses that have been put forward to explain the association between the trajectory of MD and the risk of BC. We fit three JMs with different association structures to a large Swedish mammography screening cohort which has been previously used to study MD change and BC risk [22].

Materials and methods


KARMA is a population-based prospective screening cohort [29]. All women that participated in the national screening program at four mammography units in Sweden from January 2011 to March 2013, were invited to participate in the study. Women attending screening were enrolled, irrespective of how many previous screens they had attended. A total of 70,874 women were included. Informed consent was given for a continuous collection of mammograms. Using the Swedish personal identification number, KARMA has been linked to a number of registers, including the national quality register for BC. For the present study we excluded all women without informed consent and with missing information on MD, body mass index (BMI), hormone replacement therapy (MHT) and menopausal status at baseline (MP). We also excluded prevalent breast cancer cases, women with previous cancers, except non-melanoma skin cancer, women with any breast surgery, including those with breast enlargement and/or breast reduction. Moreover, we included only women aged between 40 and 74 at baseline, and ended follow-up at age 80. Because we investigate longitudinal MD, we started follow-up from the second screen after entrance to the cohort (most women - those older than the minimum screening age that had previously accepted screening invitations - had had several screens before entering the KARMA study). This selection (delayed entry) is handled appropriately in the statistical analysis (see “Methods” section). The final analysis included 40,087 women. Five hundred eighteen of these women were diagnosed with BC during follow-up. This number includes both in-situ and invasive cancers and cancers diagnosed both via screening and symptomatically. All participants signed an informed consent and the ethical review board at Karolinska Institutet approved the study.

Processed mammograms from the mediolateral oblique (MLO) view of both breasts were collected from full-field digital mammography systems. For women with BC, mammograms from the contralateral breast (i.e the breast which does not have a tumor) were used, whereas for women that were not diagnosed with BC before the end of follow-up, one side was selected randomly at baseline and follow-up measurements were based on mammograms from the same side. Dense area (cm\(^2\)) was measured using the STRATUS method which aligns mammograms from the same woman before taking density measurements [30]. We used area-based absolute density in our analysis, which is believed to be a more etiologically relevant phenotype of MD for BC risk than percent density, since it reflects the amount of tissue at risk of carcinogenesis [31].


Joint modeling is an increasingly productive area of statistical research that has developed rapidly since the 1990s [25, 32,33,34]. JMs construct a mixed-effects model to describe the evolution of the longitudinal outcome (MD) over time, and simultaneously connect that process (as a time dependent covariate) to an outcome (BC) using a time-to-event model. JM accounts for loss of follow-up (correctly handles the potential association between the longitudinal measurements and drop-out) while taking random variation into account.

Longitudinal submodel

For modeling MD over time we use a mixed-effects model [35]

$$\begin{aligned} y_i(t_{ij}) = m_i(t_{ij} )+ \epsilon _{ij},\ \epsilon _{ij} \sim N(0,\sigma ^2) \end{aligned}$$

where \(y_i(t_{ij})\) is the jth observed longitudinal response of the continuous biomarker (MD) for the ith individual taken at time \(t_{ij}\). Measurement error is incorporated through \(\epsilon _{ij}\), and \(m_i(t_{ij})\) is modeled as:

$$\begin{aligned} m_i(t_{ij}) = X_i^T (t_{ij}) \beta +Z_i^T (t_{ij}) b_i+u_i^T \delta ,\ b_i \sim N(0,D) \end{aligned}$$

where fixed effects \(\beta\), with design matrix \(X_i\) are shared across individuals (they represent the mean trajectory of the biomarker over time), and patient-specific random effects \(b_i\) with design matrix \(Z_i\) define the individual deviation relative to the mean trajectory. The random effects are assumed to have a multivariate normal distribution with mean zero and variance-covariance matrix D. It is possible to model nonlinear trajectories by including polynomials or splines of time in both \(X_i\) and \(Z_i\), and the effect of covariates on the trajectory can be modeled by including interactions with time in \(X_i\) and \(Z_i\). Baseline covariates are included using \(u_i\) with regression coefficients \(\delta\). The model naturally handles uneven spacing of repeated measurements.

Time-to-event submodel

We use a proportional hazards model to model time to diagnosis. We assume that the hazard of experiencing the event of interest (cancer diagnosis) is dependent on a subject-specific characteristic of the longitudinal trajectory and that it can be formulated as

$$\begin{aligned} h_i(t) = h_0(t)\ \text {exp}(\gamma ^T w_i + f(m_i(t),b_i,\alpha )), \end{aligned}$$

where \(h_0(t)\) describes the (baseline) hazard. The effect parameters \(\gamma\) describe how the hazard varies as a function of explanatory covariates \(w_i\). The parameter \(\alpha\) quantifies the association between a priori selected features of the longitudinal process and the hazard for the event at time t. Several options for the function f() are possible and lead to different forms of the time-to-event submodel [26]. The first considered function is:

$$\begin{aligned} f(m_i(t),b_i,\alpha ) = \alpha \ m_i(t) \end{aligned}$$

This formulation assumes that the current ‘true’ level of MD is directly associated with the instantaneous risk of BC. The second considered function includes an additional term which is the current slope (rate of change) of MD:

$$\begin{aligned} f(m_i(t),b_i,\alpha ) = \alpha _1 m_i(t)+ \alpha _2 m_i'(t) \end{aligned}$$

The coefficients \(\alpha _1\) and \(\alpha _2\) express the strength of the association between the current value and rate of change of the (true) subject trajectory (of MD) at time t and the instantaneous risk of BC diagnosis. The last considered function is based on the cumulative value of MD:

$$\begin{aligned} f(m_i(t),b_i,\alpha ) =\alpha \int _{t_0}^t\! m_i(s)\,\textrm{d}s \end{aligned}$$

This formulation accounts for the projected history of the longitudinal outcome (i.e. true/latent MD) from a user-specified initial time (age), \(t_0\), up to the current time (age), t, in the predictor of the relative risk submodel.

In the Results section we illustrate, with an example, specifically for functions (5) and (6), how choice of function impacts the relationship between the marker trajectory and (instantaneous) risk.

Specific model formulations

For the longitudinal submodel, we included an effect of MHT use (never used/former use/current use), BMI (continuous variable), and MP status (pre/post) as variables at baseline. We transformed MD by taking its square root, to have a normal distribution, prior to including it in the longitudinal submodel. To account for a nonlinear trend, we applied a natural cubic spline within the mixed-effects models (in both the fixed and random effects). We allowed the trend over time to differ according to MHT treatment and MP status. This, in an approximate manner, accounts for individual variation in the ages at which women enter menopause. We fitted the model:

$$\begin{aligned} y_i(t_{ij})&= (\beta _0 +b_{i0})+(\beta _1 +b_{i1}) B_n(t,\lambda _1)+(\beta _2 +b_{i2})B_n(t,\lambda _2)+(\beta _3 +b_{i3})B_n(t,\lambda _3) \\&\quad +\beta _4 B_n(t,\lambda _1) MHT_{1i} +\beta _5 B_n(t,\lambda _2) MHT_{1i}+\beta _6 B_n(t,\lambda _3) MHT_{1i} \\&\quad +\beta _7 B_n(t,\lambda _1) MHT_{2i} +\beta _8 B_n(t,\lambda _2) MHT_{2i}+\beta _9 B_n(t,\lambda _3) MHT_{2i} \\&\quad +\beta _{10} B_n(t,\lambda _1) MP_{1i} +\beta _{11} B_n(t,\lambda _2) MP_{1i}+\beta _{12} B_n(t,\lambda _3) MP_{1i} \\&\quad +\delta _1 BMI_i+\delta _2 MHT_{1i}+\delta _3 MHT_{2i}+ \delta _4 MP_{1i}+ \epsilon _{ij}, \ \epsilon _{ij} \sim N(0,\sigma ^2) \end{aligned}$$

where \({B_n(t,\lambda _k): k = 1,2,3}\) denotes the B-spline basis matrix for a natural cubic spline with knots \(\lambda _k\): interior knots at ages 50 and 55 years were chosen, while boundary knots at ages 43 and 65 years constrained the trends to be linear beyond these ages. We note that the multivariate distribution of the random effects naturally handles that baseline MD will be strongly associated with the level of MD reduction across time.

For the time-to-event (cancer diagnosis) process, a relative risk model with a penalized-spline-approximated baseline risk function was used in all cases. We included BMI, MHT, and family history (FH) of BC (No/yes/missing) as covariates. We used a time scale based on age since 40, as 40 is the lower limit of the screening age interval. The majority of women entered the study over the age of 40, which means that they had delayed entry times. The resulting left truncation was accounted for in our analysis. The hazard is modelled as:

$$\begin{aligned} h_i(t)&= h_0(t)\ \text {exp}(\gamma _1 BMI_i+ \gamma _2 MHT_{1i}+\gamma _3 MHT_{2i} +\gamma _4 FH_{1i}+\gamma _5 FH_{2i}\\&\quad + f(m_i(t),b_i,\alpha )) \end{aligned}$$


We used the JMbayes package in R for our analysis. A study [36] reviewing seven available JM packages in R concluded that, JMbayes [26] is the most comprehensive and expandable one. JMbayes allows for flexibility in the modeling of parametric and nonparametric baseline hazards, spline-based nonlinear longitudinal trajectories, and different association structure parameterizations. It also handles left truncation.

Model comparisons and Bayes p values

Because we use Bayesian inference software, in order to compare the different JMs in terms of goodness-of-fit, we rely on the Deviance Information Criterion (DIC), which is a hierarchical modeling generalization of the Akaike information criterion [37]. Smaller values indicate better model adjustments to the data. Reported point estimates of parameters are posterior means, and p-values are Bayes p-values, based on tail probabilities for containing the value zero [26].


Key characteristics of individuals included in our analyses are described in Table 1. The average length of follow-up was 5.44 years. The time interval between mammography rounds in this cohort was 18–24 months. The majority of women (76.3%) had completed 3 or more rounds of mammography. The maximum number of rounds of mammography was six. 518 women were diagnosed with BC during follow-up. Those women were older, more likely to be postmenopausal, to be using MHT, and to have a first-degree relative that has been diagnosed with BC, than women without a BC diagnosis by end of follow-up. Women who developed breast cancer were older at baseline and had fewer screens compared to the women who did not develop breast cancer during follow-up. This is a bi-product of the strong relationship between age and breast cancer risk. The women entering at younger ages have lower chance of being diagnosed during follow-up (and consequently have more screens) than those entering at older ages.

Fig. 1
figure 1

MD trajectory with age. Smoothed average (square root) MD values (all women), with age (the line is encapsulated within a 95% confidence band marked in blue), and individual trajectories of (square root) MD for a sample of 1000 women (black)

Table 1 Key characteristics at baseline of individuals included in the study
Table 2 Parameter estimates and 95% credibility intervals for the event process under the joint modeling analysis with different association structures

Figure 1 represents (the square root of) MD measurements (joined by lines) for 1000 randomly selected individuals. We smoothed the data points (for all KARMA women) using a natural cubic spline and added this (with 95% confidence band in blue) to the figure. Although this function represents an average, it has a pattern which is thought to be typical of an individual’s MD trajectory. MD is inversely and nonlinearly associated with age, with the largest declines observed between the ages of 45 and 60 years (during the menopausal transition) [9, 10, 13]. The large measurement error that contributes to the fluctuations in measured MD across time, within individuals, provides motivation for using a JM approach, i.e. to improve statistical efficiency.

Point estimates, with 95% credibility intervals, for all joint models—for parameters in the event and the longitudinal processes are shown in Tables 2 and 3, respectively. We refer to models using the time-to-event sub-models in Eqs. (4)–(6), as models (1)–(3), respectively. All models showed strong evidence of an association between MD trajectory and BC risk (\(P_\alpha < 0.001\) for current value of MD in model (1), \(P_{\alpha _1} < 0.001\) and \(P_{\alpha _2} = 0.005\) for current value and slope of MD respectively in model (2), and \(P_\alpha < 0.001\) for the cumulative value of MD in model (3)). According to the DIC criterion, the best model fit was obtained using the current value and slope of MD association structure (model (2), DIC = 1,074,546), followed by the cumulative association structure (model (3), DIC = 1,075,533), and then the current value of MD alone (model (1), DIC = 1,075,645). Model (1) suggested that a 1-unit increase in the value of the square root of MD corresponds to a exp(0.130) = 1.139-fold increase in the risk of BC diagnosis (2.5–97.5% credibility interval (CI): 1.084–1.190). Using model (2), we estimated that a 1-unit increase in the value of the slope is associated with a exp(− 1.160) = 0.313-fold decrease in the risk of BC diagnosis (2.5–97.5% CI 0.135–0.716)). This is counter-intuitive as it implies that reducing density increases risk. An explanation of this observation is given in the Discussion. Model (3) suggests that a unit increase in the area under the profile of the square root of MD corresponds to a exp(0.008) = 1.008-fold increase in risk (2.5–97.5% CI 1.005–1.011); Table 2.

Table 3 Parameter estimates and 95% credibility intervals for longitudinal submodels under the joint modeling analysis (2)

Table 3 shows that the estimated coefficients for the longitudinal process were similar in all JMs. In accordance with previous literature [38, 39], MD was negatively associated with BMI and MP, and was positively associated with MHT. Statistical significance of the spline interaction terms suggest that the models captured (to some degree) the influence of the timing of the menopause and MHT on the MD trajectory. For all models, the risk of being diagnosed with BC, was positively associated with current use of MHT and having a family history of BC.

Fig. 2
figure 2

Smoothed average (square root) predicted MD values. Smoothed average (square root) predicted MD values (from model (3)), by age, for women diagnosed with BC at the end of follow-up (with 95% confidence band marked by shading in blue), and for women that remained free from BC diagnosis until the end of follow-up (in black—the confidence band is not represented since it is very narrow, comparable to that based on observed MD values for all women, represented in Fig. 1)

Fig. 3
figure 3

Predicted hazard functions for the two hypothetical individuals. Trajectories of (square root) MD of two hypothetical individuals in the upper panel. Predicted hazard functions for the two individuals, based on fitted models (2) and (3), are displayed in the middle and lower panel, respectively

For each individual we obtained predicted (square root) MD values across their follow-up period based on the fitted joint model (3). We smoothed these values using splines—separately for individuals with/without a BC diagnosis during follow-up; see Fig. 2. At young ages, on average, the women with a BC diagnosis appear to have slightly higher levels of density, and their density declines later on, than is the case for women without a BC diagnosis. For the other JMs, the corresponding plots were almost identical (data not shown). We also note that when we smooth the observed MD values separately for women with and without a BC diagnosis during follow-up, we see exactly the same pattern (data not shown).

To illustrate the significance of the estimates of the parameters of the time-to-event model, we estimated hazard functions for two hypothetical individuals with trajectories of (square root) of MD; Fig. 3. These were chosen to resemble the smoothed function in Fig. 2. Both individuals were assumed to have BMI values of 23 and to not have been using HRT. We note that for model (2) because of the negative coefficient for MD slope, the corresponding plots of the hazard functions crossed each other (see Discussion for an explanation).


Using data from a large prospective mammography screening cohort, we fitted three JMs with different association structures to investigate the association between MD and BC risk. The goodness-of-fit of both models that used (additional) characteristics of the MD trajectory, improved in comparison to using only the current value of MD.

Using only one measurement of MD, it has been consistently shown that women with a high MD have a substantially higher risk of BC than women with low MD; in [1], for example, women with 75% or greater MD were estimated to have an approximately 5-fold higher risk of BC compared to women with less than 10% dense tissue. In [28] which used a (Bayesian) JM with the current value of MD considered as a categorical variable using BIRADs, the mean of the posterior distribution of the hazard ratio was around four for category d (extremely dense) versus category a (fatty). Our results using model (1) yielded a similar effect size to these studies; women with current value of MD \(>80\) cm\(^2\) have \(> 4\) fold higher risk of BC diagnosis than women with less than 25 cm\(^2\). Studies based on analysing longitudinal data on MD, have focused on specific hypotheses—namely, current change of MD having a direct effect on current BC diagnosis [21, 22]. These studies have provided conflicting conclusions and have treated density change as a categorical variable (without accounting for measurement error).

The cumulative association structure-based model (model (3)) can be motivated by Boyd et al. [17] who argue that cumulative exposure to MD may reflect cumulative exposure to hormones that stimulate cell division in the breast and hence may be an important determinant of BC incidence. Our analysis is the first to examine this hypothesis by modelling individual level information, i.e. longitudinal data.

The observation obtained from model (2), that a high (instantaneous) decrease in MD is associated with an increased rate/risk of BC (diagnosis) is contrary to a hypothesis put forward in the literature, that women who do not experience a density decrease over time have a higher risk of BC than women who experience a decrease in MD [21,22,23]. The most likely explanation of the negative coefficient for density change in model (2) is that an instantaneous decrease in density can aid tumor detection, i.e. increase screening sensitivity. Some studies have specifically demonstrated how the age-related decrease in mammographic density is related to the increase of sensitivity of mammography with age [40, 41]. An alternative explanation could be related to the fact that breast tissue undergoes a massive remodeling during menopausal transition, known as lobular involution [42]. It has been shown that women with later ages at menopause have on average shorter menopausal transition periods than women with early ages at menopause. This might indicate a higher rate of change in women with late age at menopause, which might reflect a higher risk of BC [43]; see Fig. 2. We believe that the first explanation, that MD reduction increases sensitivity, is most likely, since we model diagnosis rather than the (unobservable) onset of cancer.

Because screening is not specifically incorporated in the modelling process, the effect of density change is, inappropriately, applied across all time points and can account e.g. for the unexpected crossing of hazard curves as observed in model (2) in Fig. 3 (in practice risk of diagnosis increases at the specific times that women attend screening and may particularly increase at these time points if density has decreased). For this reason, we suggest that the cumulative association structure JM is more robust (will give parameter estimates that more closely resemble biology) than the current MD value, current MD slope structure.

Clearly there is a lag between tumor onset and detection (even more problematic is that the duration of this lag varies dramatically across individuals). A pragmatic approach to account for this could be to model lagged effects of MD characteristics on diagnosis. JMs allow for lagged effects parameterizations [26]. In the current study we did not fit such models since almost 50% of the cases did not have images more than 2.5 years before diagnosis.

Because of the nature of the study design cases had fewer screening rounds than controls (follow-up stops at diagnosis). In theory, however, a bias could be introduced if, in the general population, the disease process is associated with the visit/screening process. There are JM approaches for additionally modeling the visiting process [44], but it has been shown [45] that even when the disease process is associated with the visiting process, fitting random effects models ignoring the visiting process produces estimates with little bias.

We should note that the benefits of JM are strictly linked to the correct specification of the longitudinal marker trajectory and the baseline hazard function, indicating the need for a careful consideration of assumptions to avoid biased estimates [46]. Although there are limitations to JMs for studying the MD–BC association, JMs do provide a better framework than approaches previously used in the literature - in particular because they account for measurement error. Using longitudinal trajectories of MD that account for measurement error can be less influenced by mammography acquisition parameters, such as compressed breast thickness.

A limitation of the current study is that we lacked longitudinal information on established BC risk factors. All information on covariates was collected at baseline, or at date of the first mammogram. A recent study advocated the need for longitudinal studies assessing the impact of breast fat and body weight history on MD features and BC risk over a long follow-up among both pre-menopausal and postmenopausal women [47]. Multivariate joint modeling can offer a useful tool to determine the interrelationships between BC risk factors and MD, and their relationship with risk of BC. In addition, integrating information about MD in younger women using different imaging techniques, such as MRI [17], would be highly valuable to increase our understanding of breast density through a woman’s lifetime.

In our study we used JMs that assume that parameters measuring the strength of the association between MD and BC diagnosis are time-constant. However, it might be more natural and more, biologically, relevant to assume that the effects of characteristics of the MD trajectory changes over time. Further work is needed based on JMs with time-varying association structures [48].


Our study is the first to have used JMs to investigate the association between the longitudinal history of MD and BC risk whilst considering different association structures. We statistically investigated the hypothesis that the cumulative exposure to MD is associated with BC risk, which has been proposed earlier with biological arguments. A JM suggested that a decrease in MD may be associated with an increased (instantaneous) BC risk. It is possible that this is because of increased sensitivity rather than being related to biology. For this reason, the cumulative exposure association structure may represent a more useful model to study BC diagnosis.

Availability of data and materials

The datasets used and/or analysed during the current study are available from the corresponding author upon reasonable request and with permission of



Mammographic density


Breast cancer


Joint model


Body mass index


Menopausal hormone treatment


Menopausal status


Family history


Breast Imaging Reporting and Data System


  1. Boyd NF, Guo H, Martin LJ, Sun L, Stone J, Fishell E, Jong RA, Hislop G, Chiarelli A, Minkin S, et al. Mammographic density and the risk and detection of breast cancer. N Engl J Med. 2007;356(3):227–36.

    Article  CAS  PubMed  Google Scholar 

  2. Palomares MR, Machia JR, Lehman CD, Daling JR, McTiernan A. Mammographic density correlation with Gail model breast cancer risk estimates and component risk factors. Cancer Epidemiol Prev Biomark. 2006;15(7):1324–30.

    Article  Google Scholar 

  3. Vachon CM, Van Gils CH, Sellers TA, Ghosh K, Pruthi S, Brandt KR, Pankratz VS. Mammographic density, breast cancer risk and risk prediction. Breast Cancer Res. 2007;9(6):217.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Yaghjyan L, Colditz GA, Rosner B, Tamimi RM. Mammographic breast density and subsequent risk of breast cancer in postmenopausal women according to the time since the mammogram. Cancer Epidemiol Prev Biomark. 2013;22(6):1110–7.

    Article  Google Scholar 

  5. Lynge E, Vejborg I, Andersen Z, von Euler-Chelpin M, Napolitano G. Mammographic density and screening sensitivity, breast cancer incidence and associated risk factors in danish breast cancer screening. J Clin Med. 2019;8(11):2021.

    Article  PubMed  PubMed Central  Google Scholar 

  6. Mandelson MT, Oestreicher N, Porter PL, White D, Finder CA, Taplin SH, White E. Breast density as a predictor of mammographic detection: comparison of interval-and screen-detected cancers. J Natl Cancer Inst. 2000;92(13):1081–7.

    Article  CAS  PubMed  Google Scholar 

  7. Chiu SY-H, Duffy S, Yen AM-F, Tabár L, Smith RA, Chen H-H. Effect of baseline breast density on breast cancer incidence, stage, mortality, and screening parameters: 25-year follow-up of a swedish mammographic screening. Cancer Epidemiol Prev Biomark. 2010;19(5):1219–28.

    Article  Google Scholar 

  8. Kerlikowske K, Ichikawa L, Miglioretti DL, Buist DS, Vacek PM, Smith-Bindman R, Yankaskas B, Carney PA, Ballard-Barbash R. Longitudinal measurement of clinical mammographic breast density to improve estimation of breast cancer risk. J Natl Cancer Inst. 2007;99(5):386–95.

    Article  PubMed  Google Scholar 

  9. McCormack VA, Perry NM, Vinnicombe SJ, dos Santos Silva I. Changes and tracking of mammographic density in relation to Pike’s model of breast tissue aging: a UK longitudinal study. Int J Cancer. 2010;127(2):452–61.

    CAS  PubMed  Google Scholar 

  10. Burton A, Maskarinec G, Perez-Gomez B, Vachon C, Miao H, Lajous M, López-Ridaura R, Rice M, Pereira A, Garmendia ML, et al. Mammographic density and ageing: a collaborative pooled analysis of cross-sectional data from 22 countries worldwide. PLoS Med. 2017;14(6):e1002335.

    Article  PubMed  PubMed Central  Google Scholar 

  11. Engmann NJ, Scott C, Jensen MR, Winham SJ, Ma L, Brandt KR, Mahmoudzadeh A, Whaley DH, Hruska CB, Wu F-F, et al. Longitudinal changes in volumetric breast density in healthy women across the menopausal transition. Cancer Epidemiol Prev Biomark. 2019;28(8):1324–30.

    Article  Google Scholar 

  12. Maskarinec G, Pagano I, Lurie G, Kolonel LN. A longitudinal investigation of mammographic density: the multiethnic cohort. Cancer Epidemiol Prev Biomark. 2006;15(4):732–9.

    Article  Google Scholar 

  13. Kelemen LE, Pankratz VS, Sellers TA, Brandt KR, Wang A, Janney C, Fredericksen ZS, Cerhan JR, Vachon CM. Age-specific trends in mammographic density: the Minnesota breast cancer family study. Am J Epidemiol. 2008;167(9):1027–36.

    Article  PubMed  Google Scholar 

  14. Boyd NF, Martin LJ, Bronskill M, Yaffe MJ, Duric N, Minkin S. Breast tissue composition and susceptibility to breast cancer. J Natl Cancer Inst. 2010;102(16):1224–37.

    Article  PubMed  PubMed Central  Google Scholar 

  15. McCormack VA, dos Santos Silva I. Breast density and parenchymal patterns as markers of breast cancer risk: a meta-analysis. Cancer Epidemiol Prev Biomark. 2006;15(6):1159–69.

    Article  Google Scholar 

  16. Pike MC, Krailo M, Henderson B, Casagrande J, Hoel D. ‘Hormonal’ risk factors, ‘breast tissue age’ and the age-incidence of breast cancer. Nature. 1983;303(5920):767–70.

    Article  CAS  PubMed  Google Scholar 

  17. Boyd N, Berman H, Zhu J, Martin LJ, Yaffe MJ, Chavez S, Stanisz G, Hislop G, Chiarelli AM, Minkin S, et al. The origins of breast cancer associated with mammographic density: a testable biological hypothesis. Breast Cancer Res. 2018;20(1):17.

    Article  PubMed  PubMed Central  Google Scholar 

  18. Lokate M, Stellato RK, Veldhuis WB, Peeters PH, van Gils CH. Age-related changes in mammographic density and breast cancer risk. Am J Epidemiol. 2013;178(1):101–9.

    Article  PubMed  Google Scholar 

  19. Work M, Reimers L, Quante A, Crew K, Whiffen A, Terry M. Differences in mammographic density decline over time in breast cancer cases and women at high risk for breast cancer. Cancer Epidemiol Prev Biomark. 2013;22(3):476–476.

    Article  Google Scholar 

  20. Vachon CM, Pankratz VS, Scott CG, Maloney SD, Ghosh K, Brandt KR, Milanese T, Carston MJ, Sellers TA. Longitudinal trends in mammographic percent density and breast cancer risk. Cancer Epidemiol Prev Biomark. 2007;16(5):921–8.

    Article  Google Scholar 

  21. Kim EY, Chang Y, Ahn J, Yun J-S, Park YL, Park CH, Shin H, Ryu S. Mammographic breast density, its changes, and breast cancer risk in premenopausal and postmenopausal women. Cancer. 2020;126(21):4687–96.

    Article  PubMed  Google Scholar 

  22. Azam S, Eriksson M, Sjölander A, Hellgren R, Gabrielson M, Czene K, Hall P. Mammographic density change and risk of breast cancer. JNCI J Natl Cancer Inst. 2020;112(4):391–9.

    Article  PubMed  Google Scholar 

  23. Sartor H, Kontos D, Ullén S, Förnvik H, Förnvik D. Changes in breast density over serial mammograms: a case-control study. Eur J Radiol. 2020;127:108980.

    Article  PubMed  Google Scholar 

  24. Sweeting MJ, Thompson SG. Joint modelling of longitudinal and time-to-event data with application to predicting abdominal aortic aneurysm growth and rupture. Biom J. 2011;53(5):750–63.

    Article  PubMed  PubMed Central  Google Scholar 

  25. Wulfsohn MS, Tsiatis AA. A joint model for survival and longitudinal data measured with error. Biometrics. 1997;330–9.

  26. Rizopoulos, D.: The R package Jmbayes for fitting joint models for longitudinal and time-to-event data using MCMC. J Stat Softw. 2016;72(i07).

  27. Andersson TM-L, Crowther MJ, Czene K, Hall P, Humphreys K. Mammographic density reduction as a prognostic marker for postmenopausal breast cancer: results using a joint longitudinal-survival modeling approach. Am J Epidemiol. 2017;186(9):1065–73.

    Article  PubMed  PubMed Central  Google Scholar 

  28. Armero C, Forné C, Rué M, Forte A, Perpiñán H, Gómez G, Baré M. Bayesian joint ordinal and survival modeling for breast cancer risk assessment. Stat Med. 2016;35(28):5267–82.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Gabrielson M, Hammarström M, Borgquist S, Leifland K, Czene K, Hall P. Cohort profile: the Karolinska mammography project for risk prediction of breast cancer (KARMA). Int J Epidemiol. 2017;46(6):1740–1.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Eriksson M, Li J, Leifland K, Czene K, Hall P. A comprehensive tool for measuring mammographic density changes over time. Breast Cancer Res Treat. 2018;169(2):371–9.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Haars G, van Noord PA, van Gils CH, Grobbee DE, Peeters PH. Measurements of breast density: no ratio for a ratio. Cancer Epidemiol Prevent Biomark. 2005;14(11):2634–40.

    Article  Google Scholar 

  32. Rizopoulos D. Joint models for longitudinal and time-to-event data: with applications in R. Boca Raton: CRC Press; 2012.

    Book  Google Scholar 

  33. Tsiatis AA, Degruttola V, Wulfsohn MS. Modeling the relationship of survival to longitudinal data measured with error. Applications to survival and CD4 counts in patients with AIDS. J Am Stat Assoc. 1995;90(429):27–37.

    Article  Google Scholar 

  34. Faucett CL, Thomas DC. Simultaneously modelling censored survival data and repeatedly measured covariates: a Gibbs sampling approach. Stat Med. 1996;15(15):1663–85.

    Article  CAS  PubMed  Google Scholar 

  35. Laird NM, Ware JH. Random-effects models for longitudinal data. Biometrics 1982;963–74.

  36. Cekic S, Aichele S, Brandmaier AM, Köhncke Y, Ghisletta P. A tutorial for joint modeling of longitudinal and time-to-event data in R. 2019.

  37. Spiegelhalter DJ, Best NG, Carlin BP, Van der Linde A. The deviance information criterion: 12 years on. J R Stat Soc Ser B Stat Methodol. 2014;87:485–93.

    Article  Google Scholar 

  38. Titus-Ernstoff L, Tosteson AN, Kasales C, Weiss J, Goodrich M, Hatch EE, Carney PA. Breast cancer risk factors in relation to breast density (United States). Cancer Causes Control. 2006;17(10):1281–90.

    Article  PubMed  PubMed Central  Google Scholar 

  39. Sung H, Ren J, Li J, Pfeiffer RM, Wang Y, Guida JL, Fang Y, Shi J, Zhang K, Li N, et al. Breast cancer risk factors and mammographic density among high-risk women in Urban China. NPJ Breast Cancer. 2018;4(1):1–12.

    Article  Google Scholar 

  40. Keen JD, Keen JE. How does age affect baseline screening mammography performance measures? A decision model. BMC Med Inform Decis Mak. 2008;8(1):1–16.

    Article  Google Scholar 

  41. Buist DS, Porter PL, Lehman C, Taplin SH, White E. Factors contributing to mammography failure in women aged 40–49 years. J Natl Cancer Inst. 2004;96(19):1432–40.

    Article  PubMed  Google Scholar 

  42. Martinson HA, Lyons TR, Giles ED, Borges VF, Schedin P. Developmental windows of breast cancer risk provide opportunities for targeted chemoprevention. Exp Cell Res. 2013;319(11):1671–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Paramsothy P, Harlow SD, Nan B, Greendale GA, Santoro N, Crawford SL, Gold EB, Tepper PG, Randolph JF Jr. Duration of the menopausal transition is longer in women with young age at onset: the multi-ethnic study of women’s health across the nation. Menopause (New York, NY). 2017;24(2):142.

    Article  Google Scholar 

  44. Gasparini A, Abrams KR, Barrett JK, Major RW, Sweeting MJ, Brunskill NJ, Crowther MJ. Mixed-effects models for health care longitudinal data with an informative visiting process: a Monte Carlo simulation study. Stat Neerl. 2020;74(1):5–23.

    PubMed  Google Scholar 

  45. Neuhaus JM, McCulloch CE, Boylan RD. Analysis of longitudinal data from outcome-dependent visit processes: failure of proposed methods in realistic settings and potential improvements. Stat Med. 2018;37(29):4457–71.

    Article  PubMed  Google Scholar 

  46. Arisido MW, Antolini L, Bernasconi DP, Valsecchi MG, Rebora P. Joint model robustness compared with the time-varying covariate cox model to evaluate the association between a longitudinal marker and a time-to-event endpoint. BMC Med Res Methodol. 2019;19(1):222.

    Article  PubMed  PubMed Central  Google Scholar 

  47. Soguel L, Durocher F, Tchernof A, Diorio C. Adiposity, breast density, and breast cancer risk: epidemiological and biological considerations. Eur J Cancer Prev. 2017;26(6):511.

    Article  PubMed  PubMed Central  Google Scholar 

  48. Andrinopoulou E-R, Eilers PH, Takkenberg JJ, Rizopoulos D. Improved dynamic predictions from joint models of longitudinal and survival data with time-varying effects using p-splines. Biometrics. 2018;74(2):685–93.

    Article  PubMed  Google Scholar 

Download references


The authors thank all the participants in the KARMA study and personnel for their devoted work during data collection.


Open access funding provided by Karolinska Institute. This work was supported by the Swedish Cancer Society (Grant Number CAN 2020-0714); the Swedish Research Council (Grant Number 2020-01302) and the Swedish e-Science Research Centre.

Author information

Authors and Affiliations



MI and KH participated in the study design, performing the statistical analyses, interpretation of the results, drafting, and revising the manuscript. KC and PH participated in the study design, interpretation of the results, and helped in manuscript revision. All authors have read and approved the final manuscript.

Corresponding author

Correspondence to Keith Humphreys.

Ethics declarations

Ethical approval and consent to participate

All participants in the studies had provided written informed consent, and the studies had the approval of the ethics review board at Karolinska Institutet, Stockholm, Sweden.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher's Note

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

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Illipse, M., Czene, K., Hall, P. et al. Studying the association between longitudinal mammographic density measurements and breast cancer risk: a joint modelling approach. Breast Cancer Res 25, 64 (2023).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: