Quantitative predictive approaches for Dupuytren disease: a brief review and future perspectives

1.   Introduction

The Dupuytren disease (DD) is a benign medical condition characterised by the proliferation of fibrous tissue of the palmar fascia, which eventually leads to permanent fingers contracture and impaired hand function ^[1,2]. While this disease affects primarily the hand, other parts of the body might also be affected in the Dupuytren diathesis ^[3]. The prevalence of the disease seems to depend on age, sex and genetics, as it affects predominantly older men of Caucasian/European descent and in many cases it affect members of the same family ^[4,5]. Various studies also identified possible correlations between DD and manual labour ^[6], epilepsy ^[7], diabetes ^[4], as well as alcohol abuse and smoking ^[8].

Despite the fact that this disease has been known for at least 400 years and surgeons tried to treat it for almost 200 years ^[5], currently there is still no cure and the common treatment focuses on surgical resection of the fibrous tissue ^[4]. However, even if surgical outcomes are usually good, complications are still quite common ^[9]. Moreover, post-surgery rehabilitation is relatively long (1–3 months) and the disease recurrence rate is high (towards 50–60%) ^[10]. For all these reasons, the current research focuses on a better understanding of the biology of this disease, so that new non-surgical therapies could be developed. However, as stated in ^[4], the problem faced by clinical research into DD is represented by the systemic nature of this disease and the multitude of variables that dynamically influence Dupuytren activity. This might explain the variability in disease presentation and progression. To make advances in the treatment of this disease Eaton ^[4] suggested that a shift is needed: from a surgical research (focused on anatomy and techniques, and measuring short-term outcomes) to a chronic disease research (focused on biomarkers for individualised treatments and measuring long-term outcomes). The necessity for a biomarker-based assessment is emphasised in reference ^[4], as it can represent a quantifiable measure for DD progression and recurrence (since the rate of progression of contracture is not always documented in medical records).

Various studies have investigated potential tissue biomarkers (e.g., immune cell and cytokine composition of DD tissue) ^[11,12], as well as circulating biomarkers (e.g., circulating autoantibodies, circulating fibrocytes, matrix metalloproteases) ^{[13,14,15,16]}. However, the identification of such biomarkers has not led yet to a significantly better understanding of the biology of the disease, that would ultimately generate new non-surgical therapies. Probably the most encouraging result so far was obtained with a group of matrix metalloproteases produced by the Gram-positive bacterium clostridium histolyticum ^[5,17].

Biomarkers have been successfully used to understand the biology of other diseases that have many similarities with Dupuytren's disease. Due to the increased fibroblast proliferation and extracellular matrix remodelling involved in DD, various studies have emphasised that Dupuytren's disease has many similarities with malignant processes and wound healing ^[2,18]. It should be noted that some studies also reported an increased risk factor of advanced DD for cancer mortality ^[18,19]. In terms of research, the difference between DD and diseases such as cancer and wound healing is the absence of animal models of DD. For this reason, current DD studies are usually in vitro. In vivo studies in patients are not possible, although some studies focus on ex vivo cultures to investigate short-term (i.e., up to 7 days) myofibroblasts and extracellular matrix (ECM) interactions ^[20]. In contrast, wound healing and cancer research studies are performed in vitro, in vivo, ex vivo and in silico ^[21].

In silico modelling and computational approaches have been used extensively over the past decades to understand various biological mechanisms behind wound healing and cancer ^[22,23]. These approaches allow researchers to formulate in a concise manner (e.g., via equations that describe biological mechanisms) a variety of biological hypotheses, which can then be investigated via numerical simulations to see how the biological system behaves. These in silico approaches allow the freedom to vary specific biological assumptions, to test their impact on the evolution of the system, and thus plan further experimental research.

The purpose of this review is twofold: (i) briefly summarise the pathogenesis and treatment options for Dupuytren' disease, and (ii) briefly summarise the few in silico approaches that have been already used to study the Dupuytren's disease, and to discuss new approaches that could be used in the future (since such approaches have already been applied to malignant processes and wound healing). Here, we do not focus on statistical modelling and data analysis, which are the standard approaches for investigating clinical data related to DD ^[24,25]. Rather we focus on mechanistic mathematical modelling approaches, which allow for formal representations of biological hypotheses that can be used to make new mechanistic predictions.

We start our review in Section 2 with a concise discussion of Dupuytren's disease biology, possible pathology and current treatment options. Then, in Section 3 we review briefly the use of mathematical models for the prediction of the evolution of wound healing and malignant processes, after which we discuss the very few mathematical models and approaches that have been used so far to predict the evolution of DD. We conclude in Section 4 with a discussion of some open problems around the use of these mathematical models to advance our current understanding of DD.

2.   Disease background: biology, pathology and treatments

As mentioned in the Introduction, the Dupuytren disorder is characterised by the development of fibrous tissue in the palm and/or fingers, in the form of nodules and cords; see also Figure 1. This disorder has been classified into three different stages, based on clinical presentations and biological characteristics ^[1,2,26]. Here, we follow the 1959 classification of Luck ^[26] to which we add some more recent information regarding the biological characteristics of each stage.

● Stage 1 (proliferation phase) is characterised by the formation of thickened nodules in the palm or on the fingers. During this phase, fibroblasts (around microvessels) proliferate significantly and transform into myofibroblasts. The nodule is associated with an inflammatory process, as in addition to fibroblasts and myofibroblasts one can find also M1-like (pro-inflammatory) macrophages and lymphocytes ^[12,27]. The cells secrete a variety of cytokines and growth factors (e.g., TGF-$\beta$, GM-CSF, IL-6) from both myofibroblasts as well as inflammatory immune (M1) macrophages. There is little collagen at this stage and the cells inside the nodules have a random orientation ^[26].

● Stage 2 (involutional phase) is characterised by the contraction of nodules, which become smaller in diameter ^[26]. During this second stage, the myofibroblasts are the predominant cell type, and they become aligned along the longitudinal tension lines of the hand ^[5]. As the cells become more mature, their number decreases and the proportion of collagen increases. Normal fascial bands in the hand become diseases cords, which are collagen-rich structures that lack myofibroblasts ^[28]. As time progresses, the cords shorten leading to joint and soft-tissue contracture ^[3]. Some surgically-resected cords contain both cell-dense and acellular regions, being classified as nodule-cord units ^[26,29].

● Stage 3 (residual phase) is characterised by fingers contracture. Usually this affects the metacarpophalangeal, as well as the proximal interphalangeal joints. During this phase the density of cells decreases significantly and the cords are formed mostly of collagen ^[5]. Whatever cells have deposited the collagen, it was suggested that they could not have been similar to the nodule myofibroblasts since these cells do not deposit the collagen they produce ^[30].

We emphasize here that the different hypotheses have been proposed regarding the development of cords ^[28]: (i) a nodule can develop into a cord as the disease progresses with time, and thus the phenotype of the nodule fibroblast progresses towards the phenotype of the cord fibroblast; (ii) the nodule and the cord are two different stages of the disease which arise independently, and thus the phenotypes of nodule and cord cells might arise from a different precursor cell. The in vitro study in reference ^[28] supports the first hypothesis, which is generally accepted in the literature. However, there are reports that sometimes a cord can develop without a nodule ^[3], and therefore the discussion on the origin of these two different structures is still ongoing ^[29,30]

The evolution of the disease is the result of both mechanical and chemical processes involving cell motility, cell-cell interactions via contractile and alignment forces and the cytokines/chemokines secreted by the cells themselves ^[31,32]. Over the years it has become clear that a better understanding of the molecular and cellular mechanisms involved in the evolution (and recurrence) of DD will help the development of new treatment approaches. The most important cells involved in the development and evolution of DD are the fibroblasts, myofibroblasts and macrophages (see also Table 1). A few other studies have mentioned that T lymphocytes might also play an important role in the development of this disease ^[11,33,34].

2.1.   Treatment approaches for DD

While there is no cure for this disease, there are various treatment approaches aimed at managing it. These approaches are broadly classified as surgical and non-surgical techniques. We note here that indications for treatment are based on the impact of the disease on the quality of life of patients, which can be impacted also by their professions (e.g., manual worker, musician, etc.).

Surgical techniques are usually applied to advanced disease stages, when there is an important loss of hand function, and in case of rapid DD progression ^[46]. They represent the most common treatment approaches, and they range from needle fasciotomy (a minimally-invasive approach that uses hypodermic needles to induce cord rupture) and fasciotomy through minimal skin incisions, to fasciectomy (excision of fibrous cord through larger skin incisions, which may be incomplete in the hand but must be almost complete in the fingers in order to avoid recurrence), and dermo-fasciectomy (excision of skin and tissue with skin grafts or local flaps) ^[16]. Surgical treatments can be associated with complications, which range from neuro-vascular injury, to infections, haematoma, skin necrosis, stiffness, complex regional pain syndrome, recurrence, etc. Overall recurrence rates are between 2% and 60–66% ^[3,24], and are more frequent in younger patients ^[47]. Disease recurrence depends on the specific surgical procedures. Fasciectomy has a lower recurrence rate ($\approx 15\%$) compared to fasciotomy ($\approx 43\%$) ^[3]. Dermofasciectomy has an even lower recurrence rate ^[3].

Another minimally-invasive approach, considered a surgical procedure by some studies ^[48] and a non-surgical procedure by other studies ^[5,49], is the Collagenase Clostridium Histolyticum (CCH) injection. This procedure causes collagen lysis leading to cord rupture following external traction. The in vitro study in reference ^[49] showed that CCH treatment can also lead to a reduction in cell attachment, spreading and proliferation.

Non-surgical treatments are usually applied to the early disease stages, and include irradiation, pharmacological approaches (e.g., intra-lesion steroid injection, anti-TNF-$\alpha$ and anti-TGF-$\beta$ agents), physical therapy (e.g., ultrasonic waves and heat), splinting to stretch the digits. Unfortunately, most of these non-surgical approaches are not considered very effective in reversing or delaying the disease ^[1]. This is probably because many of these non-surgical approaches are based on a less than perfect understanding of the biology of DD (e.g., intra-lesion cortisol injections induce a M1$\to$M2 macrophage re-polarisation; this re-polarisation has been investigated intensely in the context of diseases such as cancer ^[50] and wound healing ^[51], but not in DD).

For the past years, experimental studies have been focusing on understanding better the biology of DD, with the goal of identifying new therapeutic targets ^[12], as well as identifying biomarkers for disease progression and treatment outcomes ^[4,16]. While there is still significant work to do for the identification of more useful non-surgical treatment approaches and biomarkers to assess the progression of the disease as well as treatment outcomes, it is worth recalling that for other diseases (e.g., cancer, wound healing), progress has been helped also by the use of various in silico predictive approaches ^[52].

3.   Quantitative predictive approaches: mathematical modelling and simulations

Various studies have emphasised the similarities between the biology of Dupuytren disease and wound healing ^[2,53], thus stressing that DD is an exaggerated form of wound healing. Other studies have mentioned different aspects that DD and neoplastic processes have in common ^[53]: from chromosomal abnormalities ^[18], to abnormal signalling pathways ^[54], matrix remodelling via MMPs ^[16] and recurrence following surgery ^[18]. For this reason, we start our discussion on predictive mathematical approaches by reviewing briefly in Section 3.1 a few such quantitative approaches that have been applied for a very long time to understand wound healing and tumour evolution. In Section 3.2 we will review the very few mathematical studies that have been used so far to investigate some aspects of DD.

3.1.   Quantitative approaches for understanding and predicting malignant processes and wound healing

The dynamics of fibroblasts, myofibroblasts, as well as growth factors and cytokines involved in the evolution of Dupuytren disease (e.g., TGF-$\beta$, TNF-$\alpha$) have been modelled mathematically in the context of tumour evolution ^{[55,56,57,58,59,60]} and wound healing ^{[61,62,63,64,65,66,67,68,69,70,71,72,73,74,75]}. As experimental studies on cancer and wound healing have started emphasising the complexity of interactions between cells and their environment ^[76,77], the focus of mathematical models for tumour development and evolution, as well as of those models for wound healing, has shifted from single-scale modelling (e.g., ^[70]) to multi-scale modelling (e.g., ^{[60,72,78,79]}); see also Figure 2. The single-scale modelling focuses on biological processes taking place inside the cells (i.e., signalling pathways) or processes involving cell dynamics (e.g., cell activation/inactivation, proliferation, cell death, cell-cell interactions, cell movement). The multi-scale modelling focuses on the processes that occur at different spatial scales (i.e., molecular, cellular, or even tissue scale), which interact with each other to give rise to the observed biological phenomena: tumour growth and invasion, wound healing, etc.

While many of these mathematical models focus exclusively on in silico numerical simulations, some compare these in silico results with experimental in vitro and in vivo results, with the ultimate goal of using the mathematical modelling to make new predictions. Next, we briefly discuss the modelling results of two studies mentioned above: one focused on the biological mechanisms of tumour evolution through interactions with different components of the microenvironment, and the other one focused on the biological mechanisms of wound healing.

● In reference ^[55] the authors have developed a mathematical model for the interactions between tumour epithelial cells, fibroblast and myofibroblasts via TGF-$\beta$ (secreted by tumour cells, and which transforms fibroblasts into myofibroblasts) and EGF (secreted by fibroblasts and myofibroblasts, which influences the proliferation of tumour epithelial cells). They used this model to predict the effects of changing some rate parameters (e.g., the rate at which fibroblasts become myofibroblasts) on the dynamics of epithelial cells, fibroblasts and myofibroblasts. Comparison with experimental results showed that their prediction of the growth of tumour-epithelial cells as one increases the rate at which fibroblasts become myofibroblasts was correct.

● In reference ^[73] the authors developed a mathematical model for wound healing in normal and diabetic cases, and to this end they investigated the interactions between fibroblasts, macrophages and collagen, in response to TGF-$\beta$, PDGF and hyaluronic acid. The model was then used to simulate various treatment protocols (e.g., using Dermagraft or Apligraf), and showed that the numerically predicted time-to-healing results correlate with those observed in clinical trials.

To conclude this section, we note that mechanistic mathematical models have been combined with in silico simulations and wet-lab experiments to investigate the role of various other cells and cytokines in a multitude of other diseases: from the role of macrophages in the onset and progression of inflammatory bowel disease ^[80], to the role of fibroblasts, epithelial cells, TGF-$\beta 1$ and PGE$_{2}$ in tissue repair during pulmonary fibrosis ^[81], or the role of TNF-$\alpha$ on drug-induced liver injury ^[82]. Therefore, mathematical predictive approaches (and computer simulations) could be also used to investigate various biological aspects involved in the evolution of Dupuytren disease.

3.2.   Quantitative approaches for understanding Dupuytren disease

Despite the many similarities between the Dupuytren disorder and other diseases (see previous Section), there have been very few quantitative studies of DD that use mathematical modelling approaches. In the following we discuss a few such models, following a dynamic vs. static classification. A dynamic model accounts for time-dependent changes in the variables of interest (e.g., changes in the level of cytokines or growth factors as DD evolves). In contrast, a static model focuses on the equilibrium states of the system (which can be averaged over some time intervals).

Dynamic models. One of the earliest and very few mathematical studies that we could find in the literature of DD was published by Moermans ^[83]. In this study, the author combined a statistical approach (i.e., follow-up life table analysis) and a mathematical modelling approach to predict the recurrence of DD after surgery. The mathematical model was described by a simple curve for the cumulative proportion ($N$) of DD cases that develop recurrence after $t$ years since surgery: $N(t) = K(1-e^{-kt})$; see also Figure 3. Here, $e^{-kt}$ describes the probability of no recurrence, parameter $K$ describes the max proportion of cases, and $k$ can be seen as the rate of recurrence.

While we could not find any other specific dynamic model aimed at understanding the evolution of Dupuytren's disease, we mention here the studies in references ^[68,69], where the authors developed a mathematical model aimed at understanding the inflammatory growth-factor-mediated fibro-proliferative responses in wound-healing disorders. In reference ^[69] the authors mentioned Dupuytren's contracture as an example of "non-inflammatory, non-invasive (myo-) fibroplastic lesions", but they did not focus on this disease. However, since recent experimental studies ^[34] have shown that DD is characterised by an inflammatory response, the model developed in references ^[68,69] could be applicable directly also to DD. For this reason, and to exemplify how one can build deterministic mathematical equations for time-variations in cells/molecules, in the following we briefly describe this mathematical model. To this end, we start by emphasising that the mathematical model in references ^[68,69] is an example of a dynamic model described by partial differential equations, which investigates the changes with respect to time ($t$) and space ($x$) in cell densities and growth factors involved in wound healing as well as in DD. In reference ^[69] the authors focused on the biological mechanisms that lead to decreases/increases in: the densities of fibroblasts ($n(x, t)$) and myofibroblasts ($m(x, t)$), the concentrations of collagen ($\rho(x, t)$) and a generic pro-inflammatory growth factor ($c(x, t)$; e.g., platelet-derived growth factor). The equations have the following structures (see also the diagram in Figure 4, which corresponds to these equations):

It is assumed that fibroblasts can move randomly, as well as in a directed manner towards higher concentrations of the growth factor. There is also the possibility of a passive movement of cells due to moving tissue (i.e., ECM), which is described by a convection term; see also Figure 4(b). The myofibroblasts and collagen fibres are assumed to move only via passive convection due to moving tissue. A fifth equation was derived in references ^[68,69] to quantify the displacement of the tissue ($u(x, t)$) due to viscous and elastic stresses, as well as displacement in response to traction stresses exerted by the fibroblasts and myofibroblasts, and impeded by the collagen fibres:

This mathematical model was simulated numerically, to see whether it behaves as expected, i.e., the increase in the density of fibroblasts and myofibroblasts is associated with a tissue displacement. Moreover, long-term model behaviour was investigated (via steady-state analysis) to assess whether the model predicts correctly the permanence of tissue displacement (i.e., wound healing).

The model discussed above was developed using a reductionist approach, as it focused on the spatio-temporal changes in the densities of specific variables which were identified as the most important: fibroblasts, myofibroblasts, collagen fibres, and one generic pro-inflammatory growth factor. However, the microenvironments of Dupuytren contracture nodules and cords are much more complex, as discussed in Section 2. To understand them, two review studies published in 2011 ^[84] and 2012 ^[85] suggested that Dupuytren's disease could be investigated using a systems biology approach. This is a holistic approach that combines mechanistic mathematical modelling of almost all identified components of a system, with numerical simulations of the equations describing the temporal and/or spatio-temporal changes in those components. Rehman et al. ^[84,85] suggested that systems biology could be used, for example, to model and evaluate the complex networks of cell-cell and cell-matrix interactions that control the remodelling of extracellular matrix during DD progression, or to evaluate the metabolic states of various cells inside the Dupuytren nodules/cords, as a function of the contractile activity of these cells (a metabolic profiling of cell cultures from DD was performed in reference ^[25]). As discussed before, more and more experimental studies have been published on the signalling pathways that promote DD ^[12,86], and thus new data is becoming available to parametrise detailed mathematical models. However, as far as we know, no mathematical models have been developed yet to investigate these suggested pathways in the context of DD progression.

Returning to our discussion in Section 3.1, we note that a systems biology approach for DD would be consistent with the application of systems biology approaches to cancer ^[87,88] or to wound healing ^[61,89,90].

Static models. Another class of mathematical models that have application to clinical treatment approaches of various diseases, including the Dupuytren disease, is represented by the cost-effectiveness models ^{[91,92,93,94,95]}. There are various approaches to investigate the economic costs of different DD treatments and disease recurrence: from decision-tree models ^[91] (described by hierarchical one-directional diagrams for the risk of events and transitions between disease states over a certain time; see Figure 5(a)) to Markov models ^[92,95] (described by diagrams for probabilities of staying, or moving through/between different health states; see Figure 5(b)). A recent review of such models ^[94] identified a small number of articles published between 1996–2016 on this topic that satisfied the inclusion criteria (i.e., quality of reporting and model evaluation): only 4 studies used rigorous cost-effectiveness modelling approaches to investigate optimal DD treatments. It should be noted that the results of those studies were sometimes contradictory, due to differences in the models (as well as the approximation of transition probabilities) and in the pricing systems: the modelling results in reference ^[91] suggested that fasciectomy was not a cost-effective treatment, while the results in reference ^[92] suggested that fasciectomy was the most-cost effective treatment. The review in reference ^[94] concluded that rigorous cost-effectiveness models for the economic analysis of DD treatments are currently limited.

4.   Discussion and future perspectives

In this study we reviewed briefly the biology of DD and the current treatment approaches, and by drawing comparisons with research for other diseases that share some biological similarities with DD, we focused on alternative approaches that could allow us to make some progress in understanding DD.

As suggested in reference ^[25], the absence of animal models for DD (to investigate and quantify the biological mechanisms involved in the evolution of this disease, and the impact of various drug therapies) could be overcome to some extent by in silico models and simulations. Such in silico approaches are now very common in the study of biological processes involved, for example, in wound healing and cancer progression, and even if there are many differences between these diseases and DD, the mathematical and computational methods are general enough to be applied to all these different diseases. Moreover, one advantage of using mathematical modelling and computational approaches is that they can replicate and investigate systemic multi-scale biological processes (see Figures 1 and 2), which DD actually is ^[4,48].

Unfortunately, reviewing the published literature did not reveal too many mathematical/computational models for the evolution of DD. Unless we increase the quantitative predictive research in this area, progress will likely continue to be slow.

Future mathematical modelling and computational applications. Based on the previous development of mathematical models and computational simulations to address a variety of questions related to the evolution of wound healing and cancer (as well as many other diseases), next we can summarise some of the potential future applications of such predictive models to DD:

● Predictive studies of DD evolution. As discussed in Section 3.1, a variety of mathematical models have been developed to investigate various aspects of cell movement and remodelling of extracellular matrix in the context of cancer and wound healing; see references ^{[55,56,57,65,66,68,69]} and the references therein. These mathematical models can be easily adapted to investigate specific biological processes occurring in DD. For example, such models could be used to investigate the evolution of different types of cords ^[96]: e.g., peritendinous cords, natatory cords, vertical cords, central cords, spiral cords, etc. Moreover, such models could be further validated with data on various cord progression and/or change in cord size ^[37] – an aspect that was not addressed until now, and which is a necessary step in the development of these predictive approaches.

● Predictive studies of DD relapse. Mathematical models have been used to study dynamically cancer relapse following tumour excision ^[97]. Similar mathematical models could be applied to study the relapse of DD following surgical treatments.

● Computational studies for the identification of DD biomarkers. The identification of biomarkers for cancer and wound healing represents a very active area of research. Mathematical models have been used in the past to investigate biomarkers for early disease detection ^[98,99], and such mathematical studies can be easily generalised to capture the specificities of DD.

● Mechanistic studies for the impact of patients age and disease duration. Mathematical models can also contribute addressing some of the open problems identified in reference ^[4]: e.g., propose mechanistic hypotheses on the relationship between the age of onset, duration of disease, and age at time of evaluation on Dupuytren's assessment.

● Pharmacokinetic and pharmacodynamic studies of drug treatments for DD. Mathematical equations are commonly used in pharmacokinetic and pharmacodynamic studies of optimal drug treatments ^[100]. Such approaches could be used in the future, once new drugs that could keep DD under control will start to be developed.

Conflict of interest

The authors declare that they have no conflict of interest.