1 Introduction of treatment failure in the drug-sensitive TB infectious individuals and the treatment of individuals with drug-resistant TB is the most cost-effective in curbing the co-infections of drug-resistant TB and HIV/AIDS. Okosun and Makinde [20] considered a malaria-cholera co-infection model, in the presence of prevention and treatment controls for both diseases. Applying the Pontryagin's Maximum principle, they derived the necessary conditions for optimality and equally showed using simulations, that to successfully curb malaria spread, control strategies must include cholera intervention strategies as well. Tilahun et al. [32] investigated a Pneumonia-Typhoid co-infection model with cost-effective optimal control analysis, incorporating preventive and treatment controls for both diseases. They showed that the strategy that combines Pneumonia treatment and Typhoid fever prevention is the most cost-effective in reducing the burden of the co-infections of Pnemonia and Typhoid fever. The researchers in [9] analyzed an optimal control model for HIV and TB co-infection, capturing resistance of HIV to anti-retroviral therapy. The authors showed that a combination of TB treatment and anti-retroviral therapy for HIV is the most effective in reducing the burden of the co-infection of the two diseases. In addition, Tanvi and Aggarwal [31] studied an HIV-TB co-infection, incorporating detection and treatment for both diseases. They observed that the strategy that implements detection for both diseases yielded the most economic and epidemic gains infighting the co-infections of the two diseases.
To the best of the authors' knowledge, optimal control analysis has never been applied to the coinfection of Human papillomavirus and TB in the literature. Hence, it will be imperative to investigate the optimal control and cost-effectiveness analysis of the co-infections of oncogenic HPV and TB. In particular, this paper extends the work of Omame et al. [24] by assessing the impact of HPV prevention and TB treatment controls on the control and management of the co-infections of the two diseases.
The rest of the paper is organized as follows. The model is formulated in Section 2, alongside the basic properties of the model. The co-infection model without controls is analyzed qualitatively in Section 3. The optimal control model is considered in Section 4 while Section 6 gives the concluding remarks.
The total sexually active population at time t, denoted by N H (t), is divided into eleven mutually-exclusive compartments: Susceptible individuals (S(t)), individuals infected with HPV (I hp (t)), individuals who have recovered from or cleared HPV infection (R hp (t)), individuals with persistent HPV infection (P hp (t)), individuals with latent TB (E t (t)), individuals with active TB infection (A t (t)), individuals treated of TB (T t (t)), individuals dually infected with HPV and latent TB (I p he (t)), individuals dually infected with HPV and active TB I p ha , individuals dually infected with persistent HPV and latent TB, (P p he (t)) individuals dually infected with persistent HPV and active (P p ha (t)). Thus N h = S + I hp + R hp + P hp + E t + A t + T t + I p HE + I p HA + P p HE + P p Based on the above formulations and assumptions, the HPV-TB Coinfection model is given by the following deterministic system of non-linear differential equations (flow diagram of the model is shown in Figure 1 , the associated state variables and parameters are well described in Tables 1 and 2) :
where:
λ hp = β hp [I hp + η p P hp + κ t (I p ha + η t I p he + η p P p ha + η t η p P p he )]
The parameters κ t (κ t ≥ 1) is a modification parameter accounting for the increased infectiousness of dually infected individuals due to tuberculosis, η t (η t ≤ 1) is a modification parameter accounting for reduced infectiousness of dually infected individuals due to latent TB. ω p (ω p ≥ 1) is a modification parameter accounting for the increased infectiousness of dually infected individuals due to HPV while η p (η p ≤ 1) is a modification parameter accounting for the reduced infectiousness of singly infected and dually infected individuals due to persistent HPV infection. This population is further reduced by natural death (at a rate µ h natural death occurs in all epidemiological compartments at this rate). In (2) , β hp is the effective contact rate for HPV transmission of HPV infection. 2.1 Basic properties of the co-infection model (1) without controls
The basic qualitative properties of the oncogenic HPV-TB co-infection model (1) will now be considered.
Specifically, we establish the following results.
For the oncogenic HPV-TB co-infection model (1) to be meaningful in the concept of epidemiological sense, it is necessary to prove that all its state variables are non-negative over time. The approach adopted by Omame et al [23] will be used to establish the results below.
Theorem 2.1 Let the initial data be 3 Mathematical Analysis of the model without controls
The co-infection model (1) has a DFE, obtained by setting the right-hand sides of the equations in the model (1) to zero, given by ξ 0 =(S * , I * hp , P * hp , R * hp , E * T , A * t , T * t , I p* he , I p* ha , P p* he , P p* ha )
The basic reproduction number of the oncogenic HPV-TB co-infection model (1), using the approach in van den Driessche and Watmough [35] , is given by R 0 = max{R 0H , R 0T } where R 0H and R 0T are, respectively, the oncogenic HPV and TB associated reproduction numbers, given by Figure 1 : Schematic diagram of the model (1) 8 All rights reserved. No reuse allowed without permission.
perpetuity.
preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in
The copyright holder for this this version posted September 18, 2020. . where,
Using Theorem 2 of [35] , the following result is established. 3.2 Global asymptotic stability(GAS) of the disease-free equilibrium(DFE) ξ 0 of the co-infection model
We shall apply the approach illustrated in [5] to investigate the global asymptotic stability of the disease free equilibrium of the co-infection model. In this section, we list two conditions that if met, also guarantee the global asymptotic stability of the disease-free state. First, system (1) must be written in the form:
where V ∈ R m denotes (its components) the number of uninfected individuals and K ∈ R n denotes (its components) the number of infected individuals. U 0 = (V * , 0) denotes the disease-free equilibrium of this system. The conditions (W 1) and (W 2) following must be satisfied in order to guarantee local asymptotic stability:
(W 1): For dV dt = P (V, 0), V * is globally asymptotically stable (GAS),
where B = D K Q(V * , 0) is an M-matrix (the off-diagonal elements of B are nonnegative) and Ω is the region where the model makes biological sense. If System (1) satisfies the above two conditions then the following theorem holds:
is a globally asymptotic stable (GAS) equilibrium of (1) provided that R 0 < 1 (LAS) and that assumptions (W 1) and (W 2) are met
where V denotes the number of non-infectious compartments and K denotes the number of infectious
It is clear from the above, that,Q(V, K) 0. Hence the DFE may not be globally asymptotically stable, suggesting the possibility of a backward bifurcation. This supports the backward bifurcation analysis in the proceeding section.
All rights reserved. No reuse allowed without permission.
perpetuity.
preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in
The copyright holder for this this version posted September 18, 2020. . https://doi.org/10.1101/2020.09.15.20195297 doi: medRxiv preprint
We shall carry out analysis in this section to know the type of bifurcation the model (1) may undergo, using the Centre Manifold Theory as illustrated in [6] . The following result can be established.
then model (1) exhibits backward bifurcation at R 0 = 1. If a < 0, then the system (1) exhibits a forward bifurcation at R 0 = 1.
Suppose
represents any arbitrary endemic equilibrium of the model. The existence of backward bifurcation will be studied using the Centre Manifold Theory [6] . To apply this theory, it is appropriate to do the following change of variables.
Further, using the vector notation (1) can be re-written in the form
All rights reserved. No reuse allowed without permission.
perpetuity.
preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in
The copyright holder for this this version posted September 18, 2020. . https://doi.org/10.1101/2020.09.15.20195297 doi: medRxiv preprint as follows:
with:
x i Consider the case when R 0H = 1. Assume, further, that β hp is chosen as a bifurcation parameter. Solving
Evaluating the Jacobian of the system (8) at the DFE, J(ξ 0 ), and using the approach in [6] , we have that J(ξ 0 ) has a right eigenvector (associated with the simple zero eigenvalue) given by w = [ω 1 , ω 2 , ω 3 , ω 4 , ω 5 , ω 6 , ω 7 , ω 8 , ω 9 , ω 10 , ω 11 ] T All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted September 18, 2020. . where,
The components of the left eigenvector of
where,
The associated bifurcation coefficients defined by a and b, given by:
13 All rights reserved. No reuse allowed without permission.
perpetuity.
preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in
The copyright holder for this this version posted September 18, 2020. . https://doi.org/10.1101/2020.09.15.20195297 doi: medRxiv preprint
Since the bifurcation coefficient b is positive, it follows from Theorem 4.1 in [6] that the model (1), or the transformed model (8), will undergo a backward bifurcation if the backward bifurcation coefficient, a, given by (11) is positive. Setting the HPV and TB re-infection paramters ϕp = 0, σt = 0 and the TB exogenous re-infection term ε t 1 = 0, the bifurcation coefficient, 
In this section, we shall use the Pontryagin's Maximum Principle to determine the necessary conditions for the optimal control of the oncogenic HPV-TB co-infection model. We incorporate time dependent controls into the model (1) to determine the optimal strategy for curbing the co-infections of the two diseases. Thus, we have,
subject to the initial conditions
The control functions, u1(t), u2(t), and u3(t) are bounded, Lebesgue integrable functions. The control u1(t) represents the efforts at preventing incident HPV infections. The control u2(t) aims at preventing infection with HPV by TB infected individuals so as to reduce the co-infection cases. TB Treatment control for individuals dually infected with oncogenic HPV and TB is denoted by u3(t). The controls u1 and u2 satisfies 0 ≤ u1, u2 ≤ 1 while the control u3 satisfies 0 < u3 ≤ h, where h is the TB drug efficacy used for the treatment of co-infected individuals. Our optimal control problem involves a scenario where the number of HPV-infected, TB-infected, the co-infection cases and the cost of implementing preventive and treatment controls u1(t), u2(t) and u3(t) are minimized subject to the state system (12) . For this, we consider the 14 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in 
T is the final time. We seek to find an optimal control, u * 1 , u * 2 , u * 3 , such that
where
T ] is the control set. We shall now establish the existence of such an optimal solution which minimizes the objective functional J.
Theorem 4.1 Given the objective functional J, defined on the control set U , and subject to the state system (12) with non-negative initial conditions at t = 0, then there exists an optimal control triple u * = (u1, u2, u3) such that J(u * ) = min {J(u1, u2, u3)|u1, u2, u3 ∈ U }.
We shall prove the existence of the given optimal control u * , along with the corresponding solution trajectories by following the approach used in Mohammed-Awel and Numfor [18] .
Proof :
The state functions are positive and the controls are Lebesgue measurable, therefore we have that J(ui) ≥ 0 for all ui ∈ U .
As a result, inf u i ∈U J(ui) exists and is finite. Therefore, a minimizing sequnce of controls (ui) ∈ U exists such that
Let S n , I n hp , P n hp , R n hp , E n where, ξ * = (S * , I * hp , P * hp , R * hp , E * t , A * t , T * t , I p* he , I p* ha , P p* he , P p* ha )
Since ||ui||L∞ < K for some K > 0, it follows that ui ∈ L 2 ([0, T ]), such that on a sub-sequence
Applying the lower semi-continuity of L 2 norm with respect to weak convergence, we have that
Considering the convergence of the sequences
and passing to the limit in the ordinary differential equation system (12) , we have that S * , I * hp , P * hp , R * hp , E * The Pontryagin's Maximum Principle [28] gives the necessary conditions which an optimal control pair must satisfy.
This principle transforms (12) , (13) and (14) into a problem of minimizing a Hamiltonian, H, pointwisely with regards to the control functions, u1, u2, u3: 
Proof of Theorem 4.2
is an optimal control and S * , I * hp , P * hp , P * hp , Rhp * , E * t , A * t , T * t , I p* he , I p* ha , P p* he , P p* ha are the corresponding state solutions. Applying the Pontryagin's Maximum Principle [28] , there exist adjoint variables satisfying:
with transversality conditions;
We can determine the behaviour of the control by differentiating the Hamiltonian, H with respect to the controls(u1, u2, u3) All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted September 18, 2020. . https://doi.org/10.1101/2020.09.15.20195297 doi: medRxiv preprint at t. On the interior of the control set, where 0 < uj < 1 for all (j = 1, 2, 3), we obtain 
Therefore, we have that [14] u
In this section, we shall explore uncertainty and sensitivity analyses on the parameters of the model due to lack of precision in the estimates of some of the co-infection model parameters and the uncertainty which may occur when gathering data for the simulations. We shall also carry out numerical simulations of the optimal control model (12) , in order to assess the impact of different intervention strategies on the dynamics of oncogenic HPV-TB co-infections.
The oncogenic HPV-TB co-infection model (1) (1) per LHS were run. Using the total number of individuals dually infected with HPV and latent TB (I p he ) as the response function, the parameters that strongly influence the dynamics of the co-infection model (1) are the demographic parameter, µh, the effective contact rate for HPV transmission, βhp, the parameter accounting for increased susceptibility to HPV by TB infected individuals, p 1 , and the recovery rate from HPV for individuals in I p he compartment. Also, using the population of individuals dually infected with HPV and active TB (I p ha ) as the input, the five top ranked parameters are the effective contact rate for HPV transmissibility, βhp, the effective contact rate neccesary for TB transmission, βt, the modification parameter accounting for increased susceptibility to HPV infection by active TB infected individuals, p 2 , the recovery rate from HPV for individuals in I p ha compartment and the recovery rate from TB for individuals in I p ha class. Taking the population of individuals dually infected with persistent HPV and latent TB (P p he ) as the input, the three top ranked parameters are the effective contact rate for HPV transmissibility, βhp, the fraction of individuals who recover from HPV infection and do not progress to persistent HPV infection, θ p 2 , and the recovery rate from persistent HPV infection, ρ p 3 , for individuals in P p he compartment. When the population of individuals dually infected with persistent HPV and active TB is used as the response function, the two top-ranked parameters are the effctive contact rate for TB transmission, βt and the TB treatment rate, r t 3 , for individuals dually infected with persistent HPV and active TB. Considering the HPV associated basic reproduction number, R0h, as the response function, it is observed in Table 3 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted September 18, 2020. . https://doi.org/10.1101/2020.09.15.20195297 doi: medRxiv preprint that the top four PRCC-ranked parameters are the demographic parameter, µh, effective contact rate for HPV transmission, βhp, the fraction of infected individuals who recover from HPV and do not develop persistent HPV, θ p 1 and the recovery rate from HPV φ p 1 . Finally, using the TB associated reproduction number, R0T, as the response function the two key parameters that drive the dynamics of the model are the effective contact rate for TB transmission, βt and the TB treatment rate r t 1 . 
Numerical simulations of the optimal control problem (12), adjoint equations (18) and characterizations of the control (21) are implemented by the Runge Kutta method using the forward backward sweep (carried out in MATLAB). Following the choice of the weight function for HPV vaccination control used in Malik et al [17] , the balancing factor χ1 = 10 3 .
We also assume χ2 = 500 and χ3 = 400. The demographic data related to the Shanxi province in rural China is used [8] . We assume the initial conditions to be: iii. Strategy C: control against infection with HPV by TB infected individuals (u2 = 0) and TB treatment control for dually infected individuals (u3 = 0).
Applying this control strategy, we note that the total number of individuals dually infected with HPV and TB in latent and active stages, respectively, (Figures 2 (a) and 2(b) ), and the total number of individuals dually infected with persistent HPV and TB in latent and active stages, respectively (Figures 2(c) and 2(d) ) are less when this control is applied than when the control is not applied. To be specific, this control strategy averts almost 420,352 new co-infection cases, which is the greatest number of averted cases in comparison with other control strategies implemented. The control profile for this strategy given in Figure 5 , reveals that control u1 is at its upper bound for the first 3.9 years before ultimately declining to zero. Similarly, the control u2 is at the maximum value of 90% for the first 8 months before declining to zero at time, t=4.5 years.
The simulations of the total number of dually infected individuals in the presence of TB treatment controls are depicted in Figures 3(a) -3(d) . Applying this control, we observe that the total number of individuals dually infected with HPV and TB in latent and active stages of infection respectively, is less than the total population when no control is applied as expected (Figure 3(a) and 3(b) ). Similarly, it is noticed from Figures (3(c) and 3(d) , that high population level impact is observed in the total number of individuals dually infected with persistent HPV and TB in latent and active stages of infection, respectively, when this control strategy is applied. Specifically, when this control strategy is implemented, about 300,549 new co-infection cases were averted. The control profile for this strategy presented in Figure 7 shows that control u1 is at its upper bound for the first 4.0 years before steadily declining to zero at final time. Also, the control u3 is at the maximum value of 100% for the first 1.2 years and then steadily declines to zero at time, t=4.5 years. This simulation results conforms with the epidemiological report in the introduction section that gynaecologic TB is a risk factor for oncogenic HPV infection (and subsequent cancer infection) [37] . Hence, if we focus on TB treatment controls, it can significantly reduce the burden of the co-infection of oncogenic HPV and TB in a population. The simulations are also in line with the point opined by [37] , that prior TB infection was associated with persistent HPV and increased susceptibility to cervical cancer. As a result, treating TB infections in dually infected individuals will significantly curb the mixed infections of both diseases. 20 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted September 18, 2020. (Figures 2c and 2d) , respectively, the presence of optimal vaccination control against incident HPV infection (u 1 = 0) and control against infection with HPV by TB infected individuals (u 2 = 0). Here, β hp = 2.0, β t = 8.557. All other parameters as in Table  2 Table 2 All rights reserved. No reuse allowed without permission.
perpetuity.
preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in
The copyright holder for this this version posted September 18, 2020. Figure 9 reveals that control u2 is at its upper bound for the first 4.5 years before gradually falling down to zero at final time. Moreso, the control u3 is at its minimum value for the first 2.1 years and then steadily rises to its peak value of 100% at time, t=2.7 years, before gradually declining to zero at final time. It is imperative to state categorically, that the results of the simulations are based on the parameter values given in table 2 and the initial conditions and weight constants given in Section 5.2
The cost-effectiveness analysis is used to evaluate the health interventions related benefits so as to justify the costs of the strategies [4] . This is obtained by comparing the differences among the health outcomes and costs of those interventions; achieved by computing the incremental cost-effectiveness ratio (ICER), which is defined as the cost per health outcome. It is given by:
Difference in costs between strategies Difference in health effects between strategies .
We calculated the total number of co-infection cases averted and the total cost of the strategies applied in Table 4 . The total number of co-infection cases prevented is obtained by calculating the total number of individuals when controls are administered and the total number when control is not applied. Similarly, we apply the cost functions 1 2 χ1u 2 1 , 1 2 χ2u 2 2 , 1 2 χ3u 2 3 , over time, to compute the total cost for the various strategies that we implemented. We now compare the cost-effectiveness of strategy B (Optimal HPV vaccination control for sexually active susceptible individuals against incident HPV infection and TB treatment control for dually infected individuals) and strategy C (Control against HPV infection by TB infected Figure 10 : Plots of the total number of individuals dually infected with HPV and TB in latent and active stages, respectively (Figures 10a and 10b ) as well as total number of individuals dually infected with persistent HPV and TB in latent and active stages (Figures 10c and 10d) , respectively, in the presence of control against HPV infection by TB-infected individuals (u 2 = 0) and TB treatment control dually infected individuals (u 3 = 0). Here, β hp = 2.0, β t = 8.557. All other parameters as in Table 2 individuals and TB treatment controls for dually infected individuals). From ICER (B) and ICER(C), we observe a cost saving of 0.01298 observed for strategy C over strategy B. This implies that strategy B strongly dominated strategy C, showing that strategy B is more costly and less effective compared to strategy C.
Therefore, strategy B is removed from subsequent ICER computations, shown in Table 5 . We shall now compare strategies C and A. Comparing strategy C (Control against HPV infection by TB infected individuals and TB treatment controls preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in
The copyright holder for this this version posted September 18, 2020. , showing that strategy A strongly dominated strategy C and is more expensive and less effective compared to strategy C. Consequently, strategy C (the strategy that combines and implements control against HPV infection by TB infected individuals as well as TB treatment controls for dually infected individuals) has the least ICER and is the most cost-effective of all the control strategies for the control and management of oncogenic HPV and TB co-infection. This clearly agrees with the results obtained from ACER method in Table 4 that strategy C is the most cost-effective strategy. 
In this work, we have developed and presented a co-infection model for Oncogenic HPV and TB with cost-effectiveness optimal control analysis. The full co-infection model was shown to undergo the phenomenon of backward bifurcation when the associated reproduction number is less than unity. It was further shown that TB and HPV re-infection parameters (ϕp = 0 and σt = 0) as well as the exogenous re-infection term (ε1 = 0) induced the phenomenon of backward bifurcation in the oncogenic HPV-TB co-infection model. The global asymptotic stability of the disease-free equilibrium of the coinfection model was also proven not to exist, when the associated reproduction number was below unity. The necessary conditions for the existence of optimal control and the optimality system for the co-infection model was established using the Pontryagin's Maximum Principle. Uncertainty and global sensitivity analyses were also carried out to determine the top ranked parameters that influence the dynamics of the co-infection model, when the associated reproduction numbers as well as the infected populations were used as response functions. When the population of individuals dually infected with HPV and latent TB (I p he ) was used as the response function, the parameters that strongly influence the dynamics of the co-infection model (1) are the demographic parameter, µh, the effective contact rate for HPV transmission, βhp, the parameter accounting for increased susceptibility to HPV by TB infected individuals, p 1 , and the recovery rate from HPV for individuals in I p he compartment. In addition, using the population of individuals dually infected with HPV and active TB (I p ha ) as the input, the five top ranked parameters are the effective contact rate for HPV transmissibility, βhp, the effective contact rate neccesary for TB transmission, βt, the modification parameter accounting for increased susceptibility to HPV infection by active TB infected individuals, p 2 , the recovery rate from HPV for individuals in I p ha compartment and the 25 All rights reserved. No reuse allowed without permission.
perpetuity.
preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in
The copyright holder for this this version posted September 18, 2020. . https://doi.org/10.1101/2020.09.15.20195297 doi: medRxiv preprint recovery rate from TB for individuals in I p ha class. Taking the population of individuals dually infected with persistent HPV and active TB as the response function, the two top-ranked parameters are the effctive contact rate for TB transmission, βt and the TB treatment rate, r t 3 , for individuals dually infected with persistent HPV and active TB. Numerical simulations of the optimal control model showed that: i. A combination of HPV prevention control and TB treatment control has a positive population level impact in reducing the burden of oncogenic HPV and TB co-infection cases in a population.
ii. A strategy that implements control against HPV infection by TB infected individuals and TB treatment controls for dually infected individuals can significantly reduce the burden of oncogenic HPV and TB co-infections.
iii. The strategy that combines and implements control against HPV infection by TB infected individuals as well as TB treatment control for dually infected individuals has the least ICER and is the most cost-effective of all the control strategies for the control and management of the burden of oncogenic HPV and TB co-infection.
All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted September 18, 2020. . https://doi.org/10.1101/2020.09.15.20195297 doi: medRxiv preprint