Although observations of endometrial-like cells growing outside of the uterus were made as early as the nineteenth century (1), endometriosis remains a significant and understudied public health challenge. Endometriosis is estimated to afflict 10% of menstruators and 20–25% of women undergoing surgery due to infertility or pelvic pain (2, 3). One challenge to estimating this prevalence is the variability in endometriosis presentation—with some only discovering endometriosis incidentally during surgery and others living with a wide range of debilitating symptoms (4). People with symptomatic endometriosis suffer an average of 7 years before diagnosis, a delay exacerbated by the lack of a non-surgical diagnostic for the disease (5). There is no cure for endometriosis. Rather, those with suspected or diagnosed endometriosis must decide how to combine interventions that primarily address symptoms (e.g., hormonal contraceptives and hysterectomy) and those that target endometriosis lesions specifically (e.g., ablation or excision surgeries). Unfortunately, all of these interventions affect a patient's ability to conceive and have 5 year symptom recurrence rates ranging from 10 to 62% (6).

Endometriosis patients are typically staged by the visual appearance of lesions and adhesions according to the American Society for Reproductive Medicine's revised system. However, this staging does not correlate with patient symptoms or treatment outcomes (7). Looking beyond visual characteristics, clinical and experimental studies suggest that the growth and survival of lesions is enabled by a combination of immune dysfunction, hormone dysregulation, and aberrant blood vessel development (8–10). Specifically, endometriosis patients have been observed to have differences in progesterone receptor expression and functioning (11, 12), in peritoneal cytokine profiles, and in immune cell functioning (13, 14)—which all have the potential to interfere with the efficacy of pharmacological and surgical interventions.

To understand how such complex systems can contribute to patient symptoms and treatment outcomes, we need to integrate quantitative and computational approaches with clinical and experimental techniques. Researchers have created mathematical models to predict patient diagnoses and treatment outcomes based on symptoms, measurements, and medical history. However, as the success (and failure) of therapies is increasingly recognized as dependent on system-wide biological differences, computational models will need to expand in order to understand the mechanisms connecting these differences to clinical presentations and treatment outcomes.

In this review, we will first summarize how mathematical models have been used and modified over the years to study, diagnose, and treat endometriosis (in section “Systems Biology and Computational Models of Endometriosis”). We will then explore three mathematical modeling approaches to endometriosis that each take advantage of increasing detail in biological mechanisms. For each modeling approach, we will investigate their design, use of experimental and clinical data, and the insight they provide. Lastly, we will discuss current limitations in mathematical modeling of endometriosis and possible future directions in the conclusion section.

Systems biology is an integrative approach to investigating how genetic, cellular, and tissue level differences can influence an organism's physiology. This could include using quantitative measurements, ranging from in vitro cell culture experiments to various clinical observations, to extensively characterize a biological system. These experimental and clinical data can then be analyzed using mathematical and computational modeling approaches to make predictions about how the biological system behaves under various conditions (Supplementary Table 1). But how do we represent this system complexity meaningfully in a model? There are several ways, with different levels of mechanistic detail.

The current “gold standard” for diagnosing endometriosis is laparoscopic surgery followed by histology to identify endometrial-like growths in the abdomen (5), but this surgery has several limitations that have led to it being commonly postponed or avoided. These include: its high cost, potential complications, and the need for a highly skilled endometriosis surgeon (32). Instead, blood tests, pelvic examinations, and ultrasound imaging are done to rule out other disorders. Of these methods, only ultrasound imaging can detect endometriosis; however, it is limited to only detecting one form of disease (ovarian endometriotic cysts) (20). As a result, researchers have turned to logistic regression to answer the question: Can a combination of clinical observations reliably predict endometriosis?

A logistic regression model estimates the probability of a binary outcome, such as having or not having endometriosis, using a set of independent observations about a patient as predictors (33). Regression models include components (Figure 1), in the form of these predictor variables, but they are not modeled as having any spatial or temporal relationships with one another; hence, these models are not considered mechanism-based (Table 1).

In medicine, logistic regression modeling is commonly applied to establish clinical scales, to identify risk factors for a disease, and to develop recommendations for treatment. For endometriosis, findings from logistic regression and related modeling have been cited as evidence for diagnosis guidelines [e.g., guidelines in (5) reference modeling in (34, 35)]. Specifically, researchers have used regression modeling to predict endometriosis from symptoms and medical history alone, blood tests alone, imaging alone, and a combination of these data sources (15).

Logistic regression does not require extensive prior knowledge of the mechanistic underpinnings of the disease, which can be difficult to ascertain. Instead, these models are entirely data-driven, using patient data that includes their known outcomes (e.g., diagnosis result) to predict the likely outcomes for other patients. Given sufficient data, logistic regression can identify the key elements that are predictive of endometriosis.

In this section, we will outline key points in creating a logistic regression model, using modeling by Nnoaham et al. (25) to illustrate these points. We discuss this study because it is one of the largest efforts so far to develop a non-surgical diagnosis for endometriosis, including more than 1,000 patients from 19 hospitals in 13 countries. The considerations detailed here will serve as a comparison point in later sections, where models are increasingly mechanism-based.

As part of the Women's Health Symptom Survey study in 2012, Nnoaham et al. (25) developed and validated symptom-based predictive models to predict the probability of a patient having endometriosis prior to any diagnostic surgery (Figure 2). The patients in their study all suffered from pelvic pain and/or infertility and answered over 200 questions, detailing their demographics, medical history, and symptoms. The effect of including any preoperative ultrasound data was also explored (25).

The authors used logistic regression to calculate the probability that a patient would visually be diagnosed with endometriosis at laparoscopy based on a combination of the patient observations. The authors also calculated the probability of finding “moderate” to “severe” endometriosis, according to the revised American Society for Reproductive Medicine classification system (r-ASRM stages III-IV) (25).

As shown here, regression models serve as valuable tools for identifying patient characteristics that can predict disease or treatment outcomes. Importantly, this modeling does not explain the “how” in this association, as in: “how do these patient characteristics contribute to endometriosis and treatment outcomes?” To answer this question, researchers must model the mechanisms by which components within the system affect each other.

Medicinal approaches for treating endometriosis primarily aim to manage symptoms but have limited efficacy, with symptoms often recurring once a patient stops treatment (38). As a first-line therapy, many patients presenting with a combination of chronic pelvic pain or pain during menstruation, sex, or urination will take medications such as NSAIDS and hormonal contraceptives (38). For those with persistent pain and confirmed endometriosis, therapeutic options can include gonadotropin-releasing hormone (GnRH) analogs and aromatase inhibitors (5). These second- and third-line therapies are effective in treating chronic pelvic pain through suppressing estrogen, thereby inhibiting the growth and survival of endometriosis lesions (8, 39). However, GnRH analogs and aromatase inhibitors can be associated with severe hypoestrogenic effects, such as decreases in bone mineral density (38). Emerging clinical trials aim to identify novel therapeutic strategies for treating endometriosis with increased safety through applying an array of pharmacokinetic (PK) and pharmacodynamic (PD) modeling approaches. Here, PK modeling is applied to answer the question: How much drug will a patient be exposed to over time? PD modeling then considers: As drug exposure varies, how much of a physiological response can be expected?

Although treatments for endometriosis are monitored in circulating blood (“centrally”) many drugs are delivered to or act throughout peripheral sites. To predict drug exposure or efficacy, we need a way to connect these sites. Pharmacokinetic modeling connects these central and peripheral sites through equations that predict the concentration of a drug as it is absorbed, distributed throughout the body, and eliminated via metabolism or excretion. Pharmacodynamic modeling uses these estimated and monitored drug levels to predict the onset, duration, and intensity of response to that drug (40). Population PK-PD models incorporate variability in select model parameters based on differences between patients (e.g., body mass, age, and genetic background), allowing for simulations of larger virtual populations to better inform recommendations (Figure 3).

The commonly applied two-compartment (central and peripheral) pharmacokinetic model incorporates all five elements of mechanism (Table 1)—modeling how a component (drug) moves through space and time in the context of drug dosing in order to predict phenomena, such as drug efficacy or toxic effects. Unlike regression modeling (Figure 1), pharmacokinetic models are composed of differential equations, where the variables are the concentration (or amount) of each component and the parameters are the rate constants representing how fast reactions and transitions between components and compartments occur (Figure 4).

Because of this structure and level of mechanistic detail, pharmacological models are ideal for simulating and comparing different dose amounts, regimens, and delivery sites for endometriosis therapies under development. By developing PK-PD models with additional mechanistic detail, researchers have been able to identify endometriosis patients with a genetic favorability for a GnRH antagonist (41), predict changes in bone mineral density following long-term GnRH antagonist treatment of endometriosis (42), and interrogate the role of chosen delivery method in the efficacy of combination progestin therapies (43).

As we discuss the unique considerations in population PK-PD models, we will use (26) as an example. This study by Reinecke et al. applied PK-PD modeling to select doses to be used in phase 2 of a clinical trial for an endometriosis therapy. In addition to modeling the distribution of the therapy throughout the body, this model predicted the influence of endogenous proteins on drug efficacy and adverse events in patients. This study is also of interest because of its application of multiscale data—ranging from in vitro experiments to animal experiments and previous phase 1 studies—to select the equations and parameters for this model.

In 2017, Reinecke et al. (26) used population PK-PD modeling to guide the development of an intravaginal ring (IVR), delivering the aromatase inhibitor anastrozole (ATZ) and the progestin levonorgestrel (LNG) for long-term, localized treatment of endometriosis and associated pain (26). This new approach to treating endometriosis targets estrogen production in endometriotic lesions through local inhibition of aromatase, thereby minimizing systemic hypoestrogenic effects. This therapy also includes a progestin to provide contraception because ATZ is a teratogen (44).

Population PK-PD modeling was used to identify doses that would achieve therapeutic levels of ATZ and LNG while minimizing the risk of ovarian cysts in a phase 2 clinical trial (EudraCT 2013-005090-53; NCT02203331) (26). These PK-PD models use ordinary differential equations (ODEs) to predict the levels of drugs in the body over time and the associated risk of developing ovarian cysts.

Each drug was modeled using a two-compartment model as a basis (Figure 4). Using data from in vitro and animal studies alongside their mechanistic understanding of the system, Reinecke et al. (26) amended the ATZ and LNG base models to more closely match the outcomes of a phase 1 clinical study in humans (EudraCT 2011-005620-18). As a result, these models included delivery via an intravaginal ring. In addition, since LNG predominantly binds to and influences the production of sex hormone binding globulin (SHBG) in serum, the LNG model also included the influence of LNG on SHBG (and vice versa) and the additional influence of circulating estradiol (E2) on SHBG (Figure 5).

As shown through the Reinecke et al. (26) example, population PK-PD modeling can be useful in deciding study treatments, simulating population heterogeneity, and predicting treatment response using information from in vitro, animal, and human studies. These insights inform the design of clinical trials and, ultimately, how a drug is used to treat disease. Differing from the logistic regression model, PK models predict changes in component concentrations over time, painting a dynamic picture of the system. Although the base two-compartment model is quick to create and often resembles typical drug exposure, this approach limits the questions researchers can address through modeling. As such, modelers often choose to include more mechanistic detail in their PK model and incorporate population variability in parameter values, as Reinecke et al. (26) has done.

Endometriosis treatment is complicated by the systemic effects of estrogens, gonadotropins, and related hormones throughout the menstrual cycle. Therefore, there is significant interest in understanding and predicting the effects of therapies that perturb the cycle, such as gonadotropin-releasing hormone (GnRH) analogs, aromatase inhibitors, and progestins, on endometriosis and subfertility. To do so, mechanism-based systems biology models have been created using differential equations to describe systemic hormone fluctuations that occur during the menstrual cycle (31). Quantitative systems pharmacology (QSP) connects pharmacokinetic (PK) models of hormone-modulating therapies with models of those hormones and of endogenous protein signaling in order to further study the effect of these drugs on the body. In contrast to pharmacodynamic (PD) modeling, which predicts the change in magnitude of a physiologic response, QSP modeling allows us to consider: What are the underlying mechanisms contributing to a physiological response and how can they be best therapeutically targeted?

Quantitative systems pharmacology (QSP) integrates systems biology approaches with both data-driven and mechanism-based computational techniques to understand and optimize therapies (28). Upon first glance, the structure of QSP models resembles that of PK-PD models, using differential equations to represent changes in proteins in the system over time (Figure 4). But while PK-PD models tend to be drug-centric (predicting the distribution and effects of exogenous compounds), QSP models also focus on processes endogenous to the body. QSP models thus allow us to answer questions that are more mechanism-focused than typical PK-PD models, because they model the influence of molecules from the sub-cellular to multi-organ levels, thereby including more components interacting over more spatial and temporal scales (Table 1).

QSP models have been created to explore the effects of therapies on protein signaling that impacts endometriosis. For example, Riggs et al. (46) expanded upon a mechanism-based model of bone remodeling (47) to study the effects of therapeutic estrogen-suppression to treat endometriosis (46). Importantly, Riggs et al. (46) combined their QSP model with a logistic regression model to assess how well patients' estrogen levels could predict their endometriosis-related pain severity—illustrating how the models discussed in this review can be used in harmony (46). In addition, QSP models have been used to predict in vivo treatment outcomes from in vitro systems, such as novel microphysiological systems that include the endometrium (48).

Röblitz et al. (27) created a QSP model of hormone cycling to aid in the development of GnRH analog therapies. GnRH analogs are critical in treating several conditions, including: cancers, uterine fibroids, and infertility (27). Although this study was not focused on endometriosis, we are discussing it because they model the GnRH antagonist, cetrorelix, which is used to treat endometriosis (38). Here, we will explore how these researchers integrated a highly mechanistic model of the menstrual cycle with pharmacokinetic models of GnRH analogs in order to compare treatments.

Röblitz et al. (27) modeled key hormones that travel and signal between the brain, ovaries, and the blood (Figure 7A). In the body, and specifically included in the model, GnRH is formed in the hypothalamus and transported to the pituitary gland where it stimulates the release of the gonadotropins, luteinizing hormone (LH) and follicle stimulating hormone (FSH), into the bloodstream. LH and FSH exert their effects on processes in the ovaries, which include follicular development, ovulation, and the development of the corpus luteum. Through these processes, the production and release of progesterone (P4), estradiol (E2), and inhibins A and B (IhA, IhB) are regulated. These circulating hormones signal back to the hypothalamus and pituitary, affecting the formation and release of GnRH, LH, and FSH.

Röblitz et al. (27) also created pharmacokinetic models of GnRH analog delivery to connect to these highly mechanistic models of the menstrual cycle. The delivery of GnRH agonist and antagonist are modeled using a one- and two- compartment PK model, respectively—similar to those previously described (Figure 4). Röblitz et al. (27) incorporated the pharmacokinetic model of the GnRH agonist, nafarelin, by modeling the drug in the central compartment (circulating blood) as being able to bind to and activate GnRH receptors, as natural GnRH does. In contrast, the GnRH antagonist, cetrorelix, is modeled as being able to bind to but not activate GnRH receptors (Figure 7B). In this way, the administered drugs either act alongside or compete with GnRH, thereby affecting the level of GnRH receptors available to activate downstream signaling.

The Röblitz et al. (27) model combines approaches from traditional PK-PD models with a highly mechanistic, QSP model to compare the effects of GnRH agonists and antagonists on people with normal menstrual cycles. This allows their model to efficiently predict clinical measures while supplying more insight into the biological processes affected by perturbations caused by disease or treatment than a PK-PD model alone could. Importantly, QSP models can be adapted to study different disease or treatment cases. This may involve applying the model to a new set of patients and/or adapting the model to include additional disease-related biological processes, such as in Riggs et al. (46). Ultimately, these highly mechanistic, systems biology models aim to expand (in both number and complexity) the biological questions researchers can explore.

In this review, we've explored three mathematical modeling techniques that have been applied to improve endometriosis diagnosis and treatment: regression, pharmaco-kinetics/dynamics (PK-PD), and quantitative systems pharmacology (QSP). Below, we'll summarize the benefits and challenges of each modeling approach and outline opportunities for future modeling of endometriosis.

Regression models represent a data-driven approach; meaning, they can be created without needing to start with a detailed mechanistic understanding of the system. As a result, regression models excel in identifying associations in data (e.g., which measured variables or combinations of variables are strong predictors of endometriosis or of clinical outcomes) without requiring advance knowledge of how these associations contribute to disease. However, these models are limited in their ability to explore the “how” in these associations.

PK±PD and QSP models both represent mechanism-based approaches that can be used to predict how biological factors will influence patient treatment. PK models are especially useful for deciding drug dosing in clinical studies. Although base compartmental PK models only predict the distribution of a drug throughout the body, researchers can add details about drug interactions within the body to the model (if that data is available). This leads to the creation of a more mechanistic PK-PD model. However, to better understand the role that endogenous pathways play on any disease and treatment, a QSP model is used.

QSP models closely resemble PK-PD models; however, QSP models add more focus on the biological mechanisms endogenous to the system. This leads to the inclusion of a wider range of experimental data to parametrize a QSP model (e.g., from molecular and cellular to tissue and multi-organ levels). As a result, QSP models can simulate changes within a biological system without any drug introduced—this is something PK models do not do. QSP models thereby become increasingly useful in interrogating the mechanisms underlying a drug response and contributing to disease. QSP models are also well-suited for comparing multiple disease and treatment scenarios. However, these models can be more time- and knowledge-intensive to create.

There are still many opportunities for the development and improvement of computational models to diagnose and treat endometriosis. Regression models need a large (many-patient) dataset across multiple clinical centers [such as in (25)] in order to have findings that can be generalized to other endometriosis patients. As discussed in previous reviews (15, 18), many studies on diagnostic indicators of endometriosis either have too few patients to be generalized or have yet to be validated with an independent patient population. Additionally, regression models predicting treatment outcomes are less common, so have not been discussed here. However, recent studies have used regression modeling to predict the efficacy of assisted reproductive technology and surgery on the fertility outcomes for endometriosis patients (50, 51).

The limitations of PK-PD studies often relate to the availability of sufficient data. How much data, and which data, is needed to model a therapy's delivery and effects will depend on the properties of that specific therapy. Models can be augmented with pre-clinical data and data from previous trials, as in Reinecke et al. (26). Additionally, there has been increased use of more mechanistic PK models, such as physiologically-based PK models, for investigating drug-drug interactions of therapies for endometriosis (42, 52). This could be due to the expanded tools for establishing, analyzing, and submitting these models for regulatory review (53, 54).

QSP models in general are a more recent approach. Several QSP models have been created to investigate the effects of hormone-modulating therapies on cell signaling in people with normal menstrual cycles and in people with polycystic ovary syndrome. So far, few of these models have directly modeled the effects of these therapies in endometriosis—with Riggs et al. (46) being one of the few. These researchers created a mechanism-based model to predict the effects of therapies on endometriosis symptoms and bone mineral density (46). As endometriosis is known to involve dysregulation in hormone, vascular, and immune signaling networks, there are several opportunities to use highly mechanistic computational modeling, such as QSP, to further our ability to understand, diagnose, and treat endometriosis.

For instance, the mechanism-based models of hormone signaling outlined in this review could be adapted to study the effects of hormones on endometrial tissue. One recent study has connected a hormone signaling model to a newly developed mechanistic model of endometrial changes during the menstrual cycle, including terms to represent growth, shedding, and blood vessel development (55). This and future studies can be used to explore the effects of endometriosis-associated hormone dysregulation on the endometrium.

Focusing on vascular and immune influences, researchers can adapt mechanism-based models of protein-signaling in blood vessel development (56) to study the impact of endometriosis lesions producing pro-angiogenic cytokines [e.g., VEGF, IL-1β, IL-6, IL-8, etc. (9)]. Additionally, mechanism-based, systems biology models can help us interrogate the interactions between endometrial and immune cells in endometriosis. Since endometriosis lesions and cancerous tumors share some immune and vascular abnormalities, cancer models may serve as a basis for this. For instance, macrophages are known to affect endometriosis lesions as they differentiate, secrete cytokines, and promote angiogenesis (13)—Mahlbacher et al. (57) have modeled these macrophage behaviors within cancerous tumors. Lastly, agent-based models (another mechanism-based approach) have been created to study signaling and development of epithelial tissues (58), such models can be adapted in order to investigate the functioning of healthy and endometriotic epithelia within organoid cultures.

Since each modeling approach yields distinct insight, data-driven and mechanism-based modeling can and have been used in concert to identify associations in biological data and interrogate the underlying mechanisms of disease, respectively. The harmony of these approaches was demonstrated as we discussed previous QSP models (27, 47). By using a multitude of computational modeling approaches, researchers can synthesize multiscale experimental and clinical data to identify predictors of endometriosis and design therapies. Furthermore, there are exciting opportunities for developing mechanism-based models to discern how disruptions in cell signaling affect immune, vascular, and hormone systems, and ultimately, contribute to endometriosis.