The effect of screening on the health burden of chlamydia: An evaluation of compartmental models based on person-days of infection

Abbreviations: ODE – Ordinary Differential Equation; STD – Sexually transmitted disease; STI – Sexually transmitted infection; PID – Pelvic inflammatory disease; HIV – Human immunodeficiency virus; DALY – Disability adjusted life-years; AIC – Akaike information criterion; RSS – Residual Sum of Squares

1.   Introduction

Sexually transmitted infections (STIs) have seen sharp climbs in incidence, with a ~30% increase in the U.S. ^[1] between 2015 and 2019. This trend will likely be further exacerbated due to the COVID-19 pandemic as evidence mounts on the social effects of lockdowns, the reduction in STI testing over the pandemic ^[2] and the diversion of health resources to more pressing matters ^[3]. While all STIs are of concern, chlamydia, in particular, represents a substantial health risk for the U.S. since it is the most common STI, at a staggering 1.8 million incidences ^[4]. In part, the reason for the high incidence is the ease with which it spreads, as transmission commonly occurs through vaginal, anal or oral intercourse and possible mother-to-child transmission during childbirth ^[5]. Chlamydia is also associated with myriad negative health outcomes including infertility ^[6], lymphogranuloma venereum ^[7], conjunctivitis ^[8], an increased risk of acquiring HIV ^[9] and social stigma ^[10], among numerous others. Consequently, action is required to stop the rising incidence of chlamydia, particularly through the effective deployment of health interventions.

To mitigate the spread of chlamydia, health authorities recommend several health policies. The simplest is recommending abstinence from sexual intercourse to younger demographics ^[5], in addition to practicing safe sex for all sexually active individuals. Health authorities also recommend that at-risk groups, namely women, gay and bisexual men younger than 25 years ^[5], annually screen for chlamydia, especially for those with multiple sexual partners ^[5]. To communicate these recommendations, awareness campaigns on STIs are periodically launched targeting these demographics to illustrate the impact of STIs on life, reduce STI-related stigma and ensure people acquire the tools and knowledge to prevent and test for infection ^[11]. Fortunately, even if these health interventions fail and chlamydia infection occurs, effective antibiotics treatments are available, such as the use of azithromycin ^[12] or doxycycline ^[13,14]. Unfortunately, due to the delays in the appearance of symptoms and seeking of treatment, an infection can negatively affect reproductive health in both men and women. Pregnant women in particular face severe risk since chlamydia infection may cause a fatal ectopic pregnancy or even permanent damage to their reproductive systems through pelvic inflammatory disease (PID) ^[5].

The recent uptick of chlamydia incidence in the U.S. calls to light an urgent need to evaluate strategies that may curtail the trend. To inform on such strategies, we evaluate the role that screening may have in reducing chlamydia incidence and the effects of symptomatic and asymptomatic durations of infection. To account for the variation in symptomatic and asymptomatic durations of infection on the trajectory of disease, we develop a mathematical model that permits non-exponential and non-Erlang distributed durations of infection. Typically, such a feature requires model formulations as systems of integral or integrodifferential equations ^[15,16] which are often regarded as beyond the capabilities of modelers in public health without specialist training ^[17]. Herein is part of the novelty of the presented work, as we further develop a mathematical framework that permits non-exponential and non-Erlang distributed durations of infection while retaining model formulation as a system of ODEs. While such a framework has been developed under the context of SIS and SIR models ^[18,19], it has yet to be cast into more elaborate compartmental models such as a SEAIR analog ^[20].

To illustrate the utility of a generalized SEAIR model (gSEAIR) that describes the time evolution of person-days of infection of chlamydia, we apply it to evaluate the effects of STD screening interventions on the sexually active population in the U.S. Unlike prior works on generalized differential equation compartmental models, we consider a model that describes the trajectory of person-days of infection based on two durations of infection, namely the durations of asymptomatic and symptomatic chlamydia infection. We use the gSEAIR model to measure how changes in the shape of the duration of infection distributions affect the quality of fit to data based on Akaike information criterion (AIC) and subsequently measure how increases in STD screening alter predictions on incidence averted and disability adjusted life-years (DALYs) saved over a five-year horizon.

2.   Materials and methods

In what follows, we detail our mathematical model of chlamydia transmission, as characterized by a system of ordinary differential equations (ODEs). The model describes the evolution of the quantity of person-days of infection ^[18,19], which is based on the multiplication of incidence and a time-varying average duration of infection. So, we also provide details on the formulation of the time-varying average durations of infection, i.e., the mean residual waiting-times of infection, in addition to model parameters, goodness of model fit to data, the calculation of incidence averted and DALYs averted for each intervention scenario.

2.1.   Mathematical model

We developed a compartmental model to predict the spread of chlamydia infection across the population of the U.S. The model considers five main compartments. Each compartment has two components: number of people and duration. Thus, we have the person-days susceptible to infection ($Sm$), person-days latently infected ($Em$), person-days asymptomatically infectious ($A{m}_{A}$), person-days symptomatically infectious ($I{m}_{S}$) and person-days removed from infection ($Rm$), where $m$ is the average duration of chlamydia infection, ${m}_{S}$ is the mean residual waiting-time of symptomatic chlamydia infection at time $t$ and ${m}_{A}$ is the mean residual waiting-time of asymptomatic chlamydia infection at time $t$. The rates of transition between each compartment are given by

Here, ${m}^{{'}}$ denotes the time derivative of $m$, $b$ is the per capita birth rate, $\beta$ is the transmission rate, $N$ is the sexually active population of the U.S., $1/\sigma$ is the incubation period of chlamydia infection and $1/\kappa$ represents the average duration of immunity to chlamydia infection. Additionally, the time-varying average duration of infection is calculated by

where $\varphi$ is the proportion of asymptomatic incidence.

For our model, we consider different functional forms of mean residual waiting-times. First, we assume the classical scenario when the duration of infection is exponentially distributed, which results in a constant mean residual waiting-time. Specifically, when $m, \; {m}_{S}$ and ${m}_{A}$ are constants, system (1) reduces to the classical SEAIR model. Next, we consider a generalization of a family of distributions with periodic hazard rates ^[21] that permits multiple troughs and peaks in the average durations of infection (Supplementary Materials).

Figure 1.  Compartmental diagram. The flow of susceptible person-days ($Sm$), to latently infected ($Em$) and either asymptomatically infectious person-days ($A{m}_{A}$) or symptomatically infectious person-days ($I{m}_{S}$) and recovered person days ($Rm$). Compartments are composed of individuals multiplied by the time-varying average duration of infection ($m$), the time-varying duration of asymptomatic infection (${m}_{A}$) or the time-varying duration of symptomatic infection (${m}_{S}$). For ease of presentation, birth and mortality rates are not included (see System 1 for details).

DownLoad: Full-Size Img PowerPoint

2.2.   Parameter estimation and the durations of infection

For our model, we estimate parameters through the literature (Table 1) and published data on chlamydia incidence ^[22]. We estimate the average duration of asymptomatic infection with chlamydia using a synthesis of data on the duration of asymptomatic chlamydia infection ^[23,24] (Table 1). Additional model parameters are estimated using a nonlinear least squares procedure, in conjunction with Matlab's ode45 and fmincon algorithms, which fit the SEAIR and gSEAIR models to weekly chlamydia incidence (Figure 2). Additional parameter details, including those for the calculation of DALYs, are available in Table 1 and the Supplementary Materials.

Table 1.  Parameters, values and sources.

Symbol	Parameter	Base value	Source
$N$	Sexually active population	$15.5$ million	Fit from data
$b$	Birth rate/death rate	$0.024/$year	[25]
$\beta$	Transmission rate	$0.050-0.062/\mathrm{y}\mathrm{e}\mathrm{a}\mathrm{r}$	Table S.1.
$1/\mathrm{\kappa }$	Avg. duration of immunity	$90$ days	[26]
$1/\sigma$	Incubation period	$14$ days	[26]
$\xi$	Proportion of new infections that are asymptomatic	$0.77$	[27]
${\mu }_{A}$	Avg. duration of asymptomatic infection	$190$ days	[23]
${\mu }_{S}$	Avg. duration of symptomatic infection	$10$ days	[28]
${\rho }_{I}$	Proportion of symptomatic incidence with no complications	$0.23$	[27]
${\rho }_{A}$	Proportion of asymptomatic incidence with no complications	$0.60$	[27]
${\rho }_{M}$	Proportion of incidence leading to PID	0.09	[28]
${\rho }_{E}$	Proportion of incidence with epididymo-orchitis	$0.0415$	[29]
${\rho }_{D}$	Proportion of PID cases classified as severe	0.175	Supplementary Materials
${\rho }_{Death}$	Proportion of PID cases resulting in death	2.9 $\times {10}^{-6}$	[30]
${\lambda }_{I}$	Avg. duration of uncomplicated symptomatic infection	$0.027$ years	[28]
${\lambda }_{A}$	Avg. duration of uncomplicated asymptomatic infection	$0.521$years	[23]
${\lambda }_{M}$	Avg. duration of moderate PID due to infection	$21.01$years	Supplementary Materials
${\lambda }_{D}$	Avg. duration of severe PID due to infection	$21.01$ years	Supplementary Materials
${\lambda }_{E}$	Average duration of epididymo-orchitis and its complications due to infection	$25.38$ years	Supplementary Materials
${\lambda }_{Death}$	Avg. duration of life lost from death due to PID	$36.5$ years	Supplementary Materials
${\gamma }_{RCF}$	Avg. time period of reproductive capability in women	$39$ years	Supplementary Materials
${\gamma }_{RCM}$	Avg. time period of optimal reproductive capability in men	$41.6$ years	Supplementary Materials
${\varrho }_{F}$	Avg. time period females are able to reproduce before contracting chlamydia infection	$8.5$ years	Supplementary Materials
${\varrho }_{M}$	Avg. time period males are able to reproduce before contracting chlamydia infection	$8.6$ years	Supplementary Materials
${\vartheta }_{FP}$	Years lost in reproductive capability in women due to PID	$9.49$ years	Supplementary Materials
${\vartheta }_{ME}$	Years lost in optimal reproductive capability in men due to epididymo-orchitis	$7.62$ years	Supplementary Materials
$\tau$	Proportion of chlamydia infections held by women	$0.6624$	[31]
${\omega }_{I}$	Disability weight of symptomatic chlamydia with no complications	$0.006$	[30]
${\omega }_{A}$	Disability weight of asymptomatic chlamydia with no complications	$0.0$	[30]
${\omega }_{M}$	Disability weight of incidence leading to PID	$0.114$	[30]
${\omega }_{E}$	Disability weight of incidence with epididymo-orchitis	$0.128$	[30]
${\omega }_{D}$	Disability weight of PID cases classified as severe	$0.324$	[30]

---

 | Show Table

DownLoad: CSV

Figure 2.  Trajectories of new chlamydia incidence and cumulative square error. a) The trajectory of chlamydia infection in the United States over 175 weeks and b) the square error of model predictions relative to reported data. Reported incidence (black curve), a classical SEAIR model fit ($K = 0,$grey curve) and the gSEAIR models with hazard rates that feature one to nine cosine terms ($K = 1$ to $K = 9$).

DownLoad: Full-Size Img PowerPoint

To provide insight into the formulation of the duration of infection distribution ${P}_{j}\left(t, {t}_{0}\right)$, where $j$ is either $S$ for symptomatic infection or $A$ for asymptomatic infection, we begin by expressing it in terms of the recovery rate ^[32]. Taking ${\theta }_{j}\left(t, t+\mathrm{\Delta }t\right)$ to be the probability of recovery over the time interval $= \left[t, t+\mathrm{\Delta }t\right]$, we have that

where ${\eta }_{j}\left(t\right)$ is a time-varying recovery rate or equivalently a hazard rate. It follows for small $\mathrm{\Delta }t$ that we can estimate ${P}_{j}\left(t+\mathrm{\Delta }t, {t}_{0}\right)$ as

Rearranging terms and taking the limit as $\mathrm{\Delta }t\to 0$, we have that

or equivalently

Traditionally, ${P}_{j}\left(t, {t}_{0}\right)$ is exponentially distributed, which is a result of assuming the infectious period is exponentially distributed ^[18]. Under such assumptions the hazard rate ${\eta }_{j}\left(t\right) = 1/{\mu }_{j}.$ Given this form of hazard rate, the survival function of the duration of infection is

and the mean residual waiting-time is

Note that because ${m}_{j}\left(t\right)$ is constant in this scenario it follows that system (1) reduces exactly to the traditional SEAIR model.

To evaluate the potential periodicity of the trajectory of chlamydia, we assume a more complicated form of ${P}_{j}(t, x)$ on both the asymptomatic and symptomatic durations of infection. We consider a generalization of the simplest probability density function, $p\left(t\right)$, with a periodic hazard rate ^[21], namely

where $C$ is a normalizing constant, $\rho > 0$, $\alpha \in \left(-\mathrm{1, 1}\right)$ and $w\in \left[\mathrm{0, 2}\pi \right]$, to that of a Fourier cosine series,

where $\rho > 0$, ${a}_{i}\in \left(-\mathrm{1, 1}\right),$ ${\sum }_{i = 1}^{n}\left|{a}_{i}\right| < 1$ and ${w}_{i}\in \left[\mathrm{0, 2}\pi \right]$.

From the probability density function $p\left(t\right)$ and its associated cumulative distribution function, we obtain the hazard rate (Supplementary Materials)

with ${\rho }_{j} > 0, {a}_{ji}\in (-\mathrm{1, 1}),$ and ${\sum }_{i = 1}^{K}\left|{a}_{ji}\right| < 1.$Here, ${\rho }_{j}$ represents the average duration of infection in the absence of periodic effects and ${a}_{ji}$ is the amplitude of variation in recovery with frequency $2\pi /{w}_{ji}$, where the subscript $j$ is used to denote either asymptomatic (A) infection or symptomatic (S) infection ^[33].

Given the hazard rate and its defined relation to the mean residual waiting-time, namely ${\eta }_{j} = {(m}_{j}{'}+1)/{m}_{j}$ ^[33] with ${m}_{j}\left(0\right) = {\mu }_{j}$, we have that

2.3.   Intervention scenarios and health metrics

To inform on the benefit of awareness campaigns for mitigating chlamydia transmission, we consider the effects of increasing STI screening rates. We assume that the proportion of people that get screened within one year ^[26] follows the distribution

Thus, the average time a person has between screenings is $\frac{1}{x}$ years, which yields the screening rate ^[26]

Imposing a baseline annual screening rate of 35% ^[34] of the population, it follows that ${c}_{baseline} = 0.0083\; \mathrm{p}\mathrm{e}\mathrm{r}\; \mathrm{w}\mathrm{e}\mathrm{e}\mathrm{k}$. Alternatively, if we consider intervention scenarios that increase the annual proportion of the population screened to 40%, 50% and 70%, we have that ${c}_{40} = 0.0098\; \mathrm{p}\mathrm{e}\mathrm{r}\; \mathrm{w}\mathrm{e}\mathrm{e}\mathrm{k}$, ${c}_{50} = 0.0133\; \mathrm{p}\mathrm{e}\mathrm{r}\; \mathrm{w}\mathrm{e}\mathrm{e}\mathrm{k}$ and ${c}_{70} = 0.0231\; \mathrm{p}\mathrm{e}\mathrm{r}\; \mathrm{w}\mathrm{e}\mathrm{e}\mathrm{k}$ respectively.

For these screening rate scenarios, we evaluate our model over 5 years. We estimate incidence under each mean residual waiting-time by subtracting predictions from 40%, 50% and 70% screening rates from the baseline. The same approach was also taken to estimate annual DALYs saved, which were discounted at the standard rate of 5% per year (see Supplementary Materials for further details).

2.4.   Goodness of fit

To evaluate the quality of model fit to incidence data, we calculate the AIC ^[35]. For ease of presentation, we define the list of variables as $X: = (S, E, A, I, R{, {m}_{A}, {m}_{S})}^{T}$ and the list of parameters as

when the duration of infection is exponentially distributed and

when the duration of infection follows the distribution with a periodic hazard rate. Thus, defining $\Psi = ({s}_{0}, {e}_{0}, {a}_{0}, {i}_{0}, {r}_{0}, {m}_{A}(0), {m}_{s}(0){)}^{T},$ we can represent the ODE system as

The new symptomatic infections are defined by

For these models the residual sum of squares (RSS) is

where ${y}_{i}$ is the observed incidence on the ${i}^{th}$ week ^[14,36] and $M = 175$ is the number of data points. Optimal parameters sets, ${\widehat{\varTheta }}_{j}$, and initial conditions $\widehat{\Psi }$ for the ${j}^{th}$ distribution types were then determined by minimizing RSS (Figure 4) through a combination of Matlab's ode45 and fmincon algorithms. Thus, given the optimal parameters, it follows that

where $K$ represents the number of model parameters to be estimated from observed data. The model with minimum $AIC$ is deemed the best fit.

From the AIC, we approximate the probability that the ${i}^{th}$ model is the best candidate among all models (in the sense of combining accurate predictions while limiting the possible number of parameters ^[37]) by calculating AIC weights. First defining $\mathrm{\Delta }AI{C}_{j} = AI{C}_{j}-\mathrm{min}\left(AIC\right)$, the Akaike weights ^[38] for each scenario are

3.   Results

We assessed the effectiveness of the gSEAIR model by estimating the health burden of chlamydia in the US and informing on the potential health benefits of increasing STI screening rates from the current annual coverage of 35% to 40%, 50% and 70%, respectively. For these scenarios, we estimated the annual incidence averted and DALYs averted per year relative to the baseline for both the traditional SEAIR and gSEAIR models. The gSEAIR models considered feature mean residual waiting-times for the duration of asymptomatic and symptomatic infection with up to 9 cosine terms ($K = 1$ to $K = 9$). To identify the most appropriate SEAIR and gSEAIR models, we used AIC in addition to AIC weights to identify the optimal model and estimate the probability it was optimal among all considered scenarios.

Our results show that increasing the number of cosine terms in gSEAIR decreases the square error of model predictions relative to the data (Figure 3). Additionally, the gSEAIR model based on a hazard rate with six cosine terms is optimal compared to the other candidate models (Figure 3) since it had the lowest AIC score (Table 2). Conversely, the SEAIR model (i.e., the gSEAIR model with $K = 0$) had the highest AIC score, with a ΔAIC of 300.6 (Table 2), which indicates it is the least effective of the modeling scenarios considered. This result is supported by the AIC weights, which indicate that the $K = 6$ scenario of the gSEAIR model is optimal with a probability of 0.54 (Table 2), where most other scenarios had ΔAIC scores of at least 4.5 and AIC weights below 0.05. The exception to this is $K = 7$, where the ΔAIC was 0.81 and the AIC weight was 0.36 (Table 2), suggesting this scenario is a viable alternative to the $K = 6$.

Figure 3.  Akaike information criterion score and square error. The AIC score (solid blue curve, left axis) and square error (dotted orange curve, right axes) for SEAIR ($K = 0$) and gSEAIR ($K = 1$ to $K = 9$) models relative to incidence data.

DownLoad: Full-Size Img PowerPoint

Table 2.  Chlamydia incidence averted and DALYS saved for screening rates that achieve 35%, 40%, 50% and 70% coverage of the population per year.

	$K$
	0	1	2	3	4	5	6	7	8	9
AIC Score (1000s)	3.333	3.319	3.318	3.157	3.131	3.124	3.033	3.034	3.037	3.03
$\mathrm{\Delta }$AIC Score	300.6	285.7	154.6	123.9	98.0	91.0	0	0.8	4.6	4.8
AIC weights	0	0	0	0	0	0	0.54	0.36	0.05	0.05
Annual incidence/1000 ppl
Baseline (${c}_{30}$)	55.5	55.9	56.0	55.9	55.6	55.9	55.9	55.9	55.9	56.1
Annual incidence averted (1000 ppl)
40% screened (${c}_{40}$)	6.3	6.1	7.2	6.8	7.9	7.9	8.1	8.2	8.3	8.2
50% screened (${c}_{50}$)	17.9	17.5	19.8	19.1	21.8	21.5	22.1	22.3	22.5	22.4
70% screened (${c}_{70}$)	35.8	35.2	37.6	36.8	39.4	39.3	39.8	40.1	40.2	40.3
DALYs saved per year (1000 ppl)
40% screened (${c}_{40}$)	1.73	1.68	1.97	1.89	2.17	2.18	2.24	2.24	2.29	2.24
50% screened (${c}_{50}$)	4.94	4.83	5.46	5.28	6.01	5.93	6.07	6.15	6.19	6.18
70% screened (${c}_{70}$)	9.89	9.74	10.40	10.18	10.90	10.86	11.01	11.09	11.10	11.14

---

 | Show Table

DownLoad: CSV

At a baseline screening rate of 35%, the SEAIR and gSEAIR models predicted 55.5–56.1 annual incidences of chlamydia per 1000 people. Increasing the screening rate to 40%, 50% or 70% averted 6.1–8.3, 17.5–22.3 and 35.2–40.3 annual incidences per 1000 people respectively, depending on the number of cosine terms included in the mean residual waiting-times (Table 2). Of these findings, the gSEAIR model generally predicted a greater benefit when increasing the screening rate to 40%–70%, with an additional 1.8–4.2 annual incidence averted per 1000 people when compared to the gSEAIR model ($K = 6$) to the SEAIR model.

Regarding the health burden of chlamydia, we predict increasing the screening rate to screening rate 40%, 50% or 70% will annually avert 1.68–2.29, 4.83–6.19 and 9.74–11.14 DALYs per 1000 people (Table 3), respectively. Typically, the gSEAIR models predicted a greater quantity of DALYs averted relative to the SEAIR model, except for the $K = 1$case (Table 2). When comparing the optimal gSEAIR model ($K = 6$) to the SEAIR model, predictions illustrate an additional annual 0.51 DALYs averted per 1000 people (Table 2). Averaging across all scenarios, annual DALYs averted per 1000 people were 2.06, 5.7 and 10.6 for screening rates of 40%, 50% and 70%, respectively (Table 2).

Table 3.  Duration of infection distribution parameters with periodic hazard rate ($K = 6$).

Symptomatic parameters				Asymptomatic parameters
Amplitude		Period (days)		Amplitude		Period (days)
${a}_{S1}$	0.19	${2\pi /w}_{S1}$	198.7	${a}_{A1}$	0.45	$2\pi /{w}_{A1}$	817.0
${a}_{S2}$	0.22	$2\pi /{w}_{S2}$	82.4	${a}_{A2}$	0.18	$2\pi /{w}_{A2}$	97.0
${a}_{S3}$	0.19	$2\pi /{w}_{S3}$	57.3	${a}_{A3}$	-0.04	$2\pi /{w}_{A3}$	56.1
${a}_{S4}$	0.21	${2\pi /w}_{S4}$	44.0	${a}_{A4}$	0.15	$2\pi /{w}_{A4}$	45.0
${a}_{S5}$	0.02	${2\pi /w}_{S5}$	36.4	${a}_{A5}$	0.01	$2\pi /{w}_{A5}$	36.5
${a}_{S6}$	0.10	${2\pi /w}_{S6}$	30.2	${a}_{A6}$	0.17	$2\pi /{w}_{A6}$	30.5

---

 | Show Table

DownLoad: CSV

Averaging all gSEAIR modeling scenarios illustrates the duration of asymptomatic infection peaked at 32.2 weeks around the 85^th week of the outbreak. The symptomatic infection peaked much earlier, specifically around week 22, with an average duration of 1.9 weeks (Figure 4). Interestingly, the only mean residual waiting-time (with $t > 0$) that was strictly greater than the constant average duration of asymptomatic infection for the SEAIR model was the gSEAIR with $K = 1$. In contrast, mean residual waiting-times were typically less than the average duration of symptomatic infection for the SEAIR model, although several cases briefly surpass this value near the end of the outbreak (Figure 4). Towards this regard, when $K = 5$, the mean residual waiting-time for symptomatic infection was at a minimum (Figure 4), although this may be a result of the associated probability density function decay rate (Figure 5). For asymptomatic infection, there is not a clear scenario where one of the mean residual waiting-times is consistently the minimum, as the majority of scenarios appear to converge to common functions (Figure 4, Figure 5).

Figure 4.  Mean residual waiting-times of the duration of infection of chlamydia. The average duration of infection over 175 weeks given a) symptomatic infection and b) asymptomatic infection. The value of $K$ corresponds to the number of cosine terms in the Fourier cosine series in the probability density function, with $K = 0$ corresponding to an exponential density function.

DownLoad: Full-Size Img PowerPoint

Figure 5.  Probability density functions. The log of probability density functions for hazard rates with $K = 0$ to $K = 9$ cosine terms for a) the duration of symptomatic infection and b) the duration of asymptomatic infection.

DownLoad: Full-Size Img PowerPoint

For the optimal gSEAIR model ($K = 6$) the cosine terms of the mean residual waiting-times that are most influential correspond to a period of 82.4 days with an amplitude of 0.22 for symptomatic infection and 817 days with an amplitude of 0.45 for asymptomatic infection (Table 3). For symptomatic infection, there were other cosine terms nearly as influential, specifically all amplitudes from ${\alpha }_{S1}$ to ${\alpha }_{S4}$where close to 0.2. For asymptomatic infection, the ${\alpha }_{A1}$term was dominant, with amplitudes ${\alpha }_{A2}$, ${\alpha }_{A4}$ and ${\alpha }_{A6}$ all about one-third of its value (Table 3).

4.   Discussion

The analysis of our gSEAIR model illustrates it is an effective approach for evaluating the health burden of chlamydia in the U.S. and in assessing the potential benefits of increasing STI screening rates. The optimal gSEAIR model, according to measures such as AIC and AIC weights, was the $K = 6$ scenario. This scenario predicted a greater reduction in chlamydia incidence and DALYs when increasing the annual screening rate from 35% to 40%, 50% or 70%, at least in comparison to the traditional SEAIR model. Also, our gSEAIR models illustrated that the inclusion of time-varying average durations of infection (i.e., the mean residual waiting-times) into model dynamics typically correlated to a greater predicted health benefit from these interventions, at least in comparison to the classical SEAIR model.

As expected, the scale-up of STI screening causes a reduction in incidence and DALYs. Our findings illustrate that this reduction is comparable with other STI interventions ^[39], with the free distribution of condoms and diaphragms serving as a notable example ^[40,41]. Our predictions on this reduction are most likely conservative, as we only account for the effects of STI screening and do not account for complementary interventions that would be deployed by health authorities including contact tracing, partner notification ^[42] and the administration of suppressive therapy ^[43].

A particular area where the work presented here could be informative is in the rollout of periodic presumptive treatment ^[44,45]. To elaborate, the basis of periodic presumptive treatment is the systematic treatment of at-risk groups periodically with a combination of drugs targeting prevalent (and curable) STIs, the results of which can cause reductions in STI prevalence up to 50% ^[45]. Thus, since our gSEAIR model provides details on multiple periods associated with transmission (i.e., the periods of the cosine terms in the mean residual waiting-times), it could help to inform on how frequent periodic presumptive treatment should be deployed to maximize the health benefit of the intervention.

In addition to its potential for informing on public health issues such as periodic presumptive treatment, our work also opens avenues of mathematical investigation. For instance, most traditional compartmental models are autonomous systems of nonlinear differential equations whose solutions feature a rich history of investigation through Jacobian, next-generation or Routh Hurwitz stability analysis techniques ^[46,47,48], the application of evolutionary invasion analysis to infer the direction of evolution ^[25,49,50] or the utilization of optimal control theory to inform on the ideal construction of public health policies and interventions ^[51]. In contrast, our model is potentially nonautonomous, which necessitates the use of Floquet theory ^[52] or Poincare maps ^[19] to understand its stability properties and implications for evolution, which are utilized far less in the realm of mathematical epidemiology.

Another realm for mathematical investigation is the conversion or comparison of gSEAIR to the framework of stochastic dynamics ^[50,53]. For instance, an advantage of many Markov chain models is the estimation of major and minor outbreak probabilities and their associated durations ^[54]. Thus, in theory, a stochastic analog of gSEAIR could inform on the role that the duration of infection distribution plays in shaping these outcomes. Stochastic differential equations also represent a potentially fruitful area for future work, as it may be possible to recast their more robust ability to inform on the evolutionary dynamics of pathogens ^[55] into a generalized differential equation compartmental model similar to the ones presented in this work.

Although our model focuses on chlamydia transmission in the U.S., it could easily be adapted to study other STD outbreaks such as syphilis ^[56], gonorrhea ^[19,57], hepatitis ^[46] and even HIV and malaria co-infection ^[58] in other population demographics. Further potential generalizations and refinements include the addition of disease states such as super-spreaders and individuals receiving treatment, the subdivision of compartments to reflect levels of at-risk behavior or age demographics and even the generalization of the transmission rate to a pair formulation ^[59].

As with the majority of compartmental models, our work has several limitations. To highlight several, first, our model assumes a well-mixed (homogenous) population, which thereby disregards the potential impact that heterogeneity may have on the transmission cycle and intervention. Naturally, it follows that incorporating more realistic individual-level characteristics and mixing patterns would enhance the accuracy of the predictions provided. Second, the calibration of our model relies on reported chlamydia incidence from health authorities and estimates on the proportion of asymptomatic cases. Implicit in this requirement are potential biases that may arise due to myriad treatment-seeking behaviors among population groups such as those that mistrust medical personnel or age demographics who experience greater social stigma from disease. Finally, the model formulation imposed a parametric form of the duration of infection distribution. While the proposed distribution is flexible, as it essentially can represent any function whose Fourier cosine series converges, further empirical work is needed truly to determine the shape of the distribution.

5.   Conclusions

In summary, our study evaluates a novel form of a compartmental model of chlamydia transmission based on the quantity person-days of infection. Through our model, we better reflect recent trends in chlamydia incidence in the U.S. based on the AIC score, relative to traditional differential equation compartmental models. Given this outcome, our model projects a greater health benefit from the upscaling of STI screening interventions than would be expected from traditional approaches. By informing on the evaluation of STI screening interventions, in addition to the time periods critical to the durations of chlamydia infection, our model has the capacity to uniquely contribute to the wealth of knowledge needed to make informed decisions and thereby may aid health officials in the construction of interventions to reverse increasing rates of STIs in the U.S.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

The authors wish to thank Dr. Emelie Kenny for constructive feedback that greatly improved the clarity of the work. SG was partially supported by the National Science Foundation Grant DMS-2052592.

Conflict of interest

All authors declare no conflict of interest in this paper.