A stochastic programming model of vaccine preparation and administration for seasonal influenza interventions

1.   Introduction

Recurrent influenza epidemics have caused loses in both health and economics. Each year the outbreak leads to three million infections and hundreds of fatalities in Taiwan ^[1]. Vaccination is one of the most effective approaches to control influenza. A success immunization program relies on a robust plan to integrate vaccine preparation and administration. Currently, vaccination rates in Taiwan are still far below other countries and the WHO recommendations ^[2]. With respect to the population group involving senior citizens, the vaccination rate during the 2005–2006 flu season approximately corresponded to 50%, and the rates then declined to 40% in the 2008–2009 flu season ^[3]. Whereas the vaccination rate for similar age groups in other countries is relatively higher than those in Taiwan. For example, the vaccination rates are 80% in Korea, 60% in the U.S., and 75% in the U.K ^[4].

Prior to the flu season, the government determines vaccine orders according to the historical coverage rate and target population. This approach ignores the random effect and is incapable to respond to the emergent demand. Furthermore, the dynamic of the transmission pattern causes mismatched vaccine distribution in terms of location and timing. Evidence can refer to the pandemic of new H1N1 viral in 2009, where multiple outbreaks were discovered in early than other years, and many infections occurred for adolescents and young adults ^[5]. Another challenge involves the uncertainty of vaccine efficacy because of the circulating viral strains and the lack of an effective approach to predict future epidemics. Poor immunization performance harms the tendency of an individual to receive vaccines. A notable example was the low vaccination rate in 2010 due to the adverse reactions of the monovalent vaccine in the previous season.

An immunization program is usually carried out through multiple vaccination activities. The first activity is known as the regular vaccination, and it commences on the first day of October to provide free-vaccines to the target population. At the end of the program, recipients may receive substitute vaccine products if the inventory is depleted ^[6]. It is feasible to share vaccine products for most of the population wherein an adult uses two doses of the pediatric vaccine, and a pediatric dosage is equal to half a vial of the adult vaccine. The remaining half vial is discarded to minimize the risks of injecting expired vaccines and contaminations. Based on the epidemic level, the government may launch supplementary vaccinations that provide a second opportunity for the population to gain protection against influenza. Such immunization activity aims to complement regular vaccination by either direct outreach to the target population or to provide an easy-to-access vaccination.

A challenge in implementing an intervention program for influenza involves vaccinating a large population in a short time frame. This study examines three potential vaccination activities with different cost settings and dose requirements: (1) administering target population using regular vaccines; (2) administering target population using sharing vaccines; (3) administering the non-target population in supplementary vaccination. We develop a stochastic programming model to capture an end-to-end immunization decision beginning with vaccine procurements to distributions under an uncertain attack rate, vaccine efficacy, and recipient arrival. The objective involves minimizing vaccination and infection costs of the target population as opposed to achieving herd immunity. This type of setting is more close to an actual application because influenza vaccines only cover partial viral strains, and achieving the herd immunity threshold does not ensure to eradicate epidemics. We proposed a solution approach based on Monte Carlo simulation and integer programming warm-start was presented to solve the stochastic integer program (SIP) with continuous random variables. Additionally, the optimal solution was analyzed for the deterministic problem by using expected values. We analyzed the risk of uncertainty when making decisions for seasonal influenza interventions.

2.   Literature review

The stochastic programming is a foundational model for the study of making decisions under uncertainties. In the model's most basic format, the decision-maker faces stochastic parameters in an optimization problem. Given the distribution of an uncertain parameter, she/he must determine decision variables to optimize the expected objective value. Birge and Louveaux provided the foundation of stochastic programming ^[7]. The applications of stochastic programming were summarized in the review paper ^[8]. Related to the vaccine supply chain, Yarmand et al. (2010) developed a two-stage stochastic programming model to determine vaccine allocation in different regions ^[9]. The study considered a particular scenario of administering the vaccine to the population in two stages. The first stage determined vaccine distribution quantities in each region at the beginning of the epidemic, and the second stage provided an opportunity to reinforce the vaccination campaign based on the vaccination outcome in the previous stage and the magnitude of the outbreak. It was assumed that the regional epidemics are uncertain and sampled from a stochastic disease model. The stochastic programming solution provided the augmented doses administered to the outbreak regions. The analysis suggested that reinforcing vaccines should only apply to the region without the onset of epidemics.

Before implementing a vaccination activity, it is necessary to determine operational decisions that include vaccine composition, packaging, ordering, distribution, and vaccination. The first decision is to determine viral strains included in the seasonal flu vaccine. Kornish and Keeney (2008) first explored the vaccine composition decision with a deadline ^[10]. In such a setting, deferring the decision to select viral strains has a benefit of obtaining more information about the spread of influenzas but incurred the lost production for making vaccines. They developed a discrete-time model to capture a series of decisions under uncertain epidemics. A related study developed dynamic models to determine the optimal composition of the influenza vaccine under uncertain production yields ^[11]. Scholars investigated the relationship between viral strain prediction and strain composition decisions on immunization outcomes ^[12]. The problem is formulated as a multi-stage stochastic program to maximize societal vaccination benefits. In the aspect of vaccination strategies, a prior study developed disease models to evaluate vaccination timing and coverage population for controlling influenza epidemics ^[13].

Vaccine manufacturing decisions include package formats and doses per vial. A study developed an economic model to obtain the optimal vial size of the vaccine with the considerations of open-vial wastes and uncertain arrivals ^[14]. In the area of vaccine ordering and distribution, scholars studied vaccine-ordering decisions at the clinical level under uncertain arrivals. The study investigates the performances of ordering decisions in various settings of vaccine vial size and safety stock ^[15]. Recently, Chen et al. (2014) developed a mathematical programming model for determining vaccine distributions in low- and middle- income countries ^[16].

Sewell and Jacobson (2003) considered vaccine and injection costs for childhood immunizations ^[17]. They developed an integer programming model to determine the price of combination vaccines based on the immunization schedule for children. Chick et al. (2006) analyzed supply chain coordination between vaccine suppliers and governmental institutions and showed that both parties are benefited through the coordination of price contracts ^[18]. Recently, scholars modeled various selling strategies in the setting of a single vaccine manufacturer and single retailer under uncertain supply and demand ^[11]. The case study found that both supply and demand uncertainties can be beneficial to the retailer.

Tanner et al. (2008) presented various stochastic programming models to analyze vaccination policies against influenza. The first model minimizes the total cost spent in a vaccination campaign subject to a chance constraint of the probability of a coverage rate that exceeds the vaccination threshold ^[19]. The second model considers minimizing the probability of the disease outbreak under a limited budget with respect to vaccinations in which the objective function involves minimizing the probability of a reproduction number great than one. The third model relaxes the chance constraint by setting a penalty in the objective function, and therefore, the new objective function includes both vaccination cost and the penalty of violating the chance constraints. A related study by Tanner and Ntaimo (2010) focused on vaccine allocation decisions. The study develops a stochastic program with joint chance constraints ^[20]. The irreducibly infeasible subsystem cuts were proposed to strengthen the solution quality of linear programming (LP) relaxation. Computational results demonstrated the reductions on the runtime by posting cuts and removing redundant scenarios.

Prior studies developed models to explore different operations in a vaccine supply chain. However, to the best of our knowledge, there is a paucity of studies focusing on integrated vaccine preparation and administration. It is crucial to understand the comprehensive process of an immunization program and determine a serial of decisions at the same time. The present study focuses on this unexplored area with the consideration of a national immunization program in Taiwan. Actual vaccination activities are analyzed, and the impact of uncertainties is highlighted based on a stochastic programming framework. The following sections are organized as follows. First, the notations and the SIP formulations are presented. Second, the optimal condition in the deterministic environment is analyzed as a reference policy and compared with the stochastic solution. An approximated algorithm to solve the SIP with continuous random variables is then described, and the computational result is reported with respective to the implementation of the algorithm. Finally, the case study provides the recommendations of vaccine ordering and administration and highlights the value of using the stochastic solution.

3.   Method

3.1.   Stochastic programming

The stochastic programming model provides a framework to model decision problems with uncertain parameters. A primitive formulation involves two-stage stochastic programming that determines decisions prior to and after the disclosure of unknown information. The first-stage decisions are determined irrespective of information related to random events, and the second-stage decisions are based on the first stage decision in which full information is disclosed. In the study, vaccine ordering quantities are determined at the first stage, and subsequently, allocation quantities are determined in the second stage.

Two vaccine products are considered in our study. The adult vaccine corresponds to a volume size of 0.5 ml per dose applied for populations aged 3 years and over. The pediatric vaccine contains 0.25 ml per dose for children aged less than 3 years ^[21]. The sharing of vaccine products is considered when the inventory is depleted. Thus, when the pediatric vaccine is depleted, then a recipient less than three years old will be vaccinated with a half vial of adult vaccine. The remaining dose in the opened vial must be discarded and cannot be used for other recipients. Additionally, when the adult vaccine is in short supply, individuals aged 3 years and above will receive two doses of a pediatric vaccine as an alternative for immunization purposes ^[22].

Both regular and supplementary vaccinations with different vaccination costs are considered. The regular vaccination refers to a situation in which a recipient is injected using a matched vaccine product at the vaccination location. In this phenomenon, there is a vaccination cost v¹ for each injection. If the matched vaccine is unavailable and another product type is available, then an individual receives substituted vaccine products with a different cost v². The supplementary vaccination considers all target populations when any excess of vaccine stocks at the end of the vaccination program. Therefore, vaccination cost v³ is considered the highest, and it accounts for the lateness of immunization for a given population. There is no further sharing of vaccine products in the supplementary vaccinations. The notations used in this study are summarized in the following table.

We develop a two-stage stochastic programming model integrated ordering and allocation decision in the vaccine supply chain for influenza interventions. The recourse problem (RP) model determines the optimal decisions before and after the disclosure of the uncertain information. Detail formulations state as the following:

s.t.

The objective function is to minimize the total expected costs, where the first term is deterministic representing the vaccine costs, and the second term involves uncertain parameters for various costs spent on vaccinations. The first type corresponds to the total cost of injecting matched vaccines for pediatric patients and adults during the regular vaccination $\underset{i\in I}{\mathop \sum }\, \underset{j\in J, j = i}{\mathop \sum }\, {{v}^{1}}~{{r}_{ij}}$, the second type corresponds to the total cost of injecting substituted vaccines for population $\underset{i\in I}{\mathop \sum }\, \underset{j\in J, j\ne i}{\mathop \sum }\, {{v}^{2}}~{{r}_{ij}}$, and the last type corresponds to the cost of supplementary vaccination $\underset{j\in J}{\mathop \sum }\, {{v}^{3}}~{{s}_{j}}$. Subsequently, the infection costs $\underset{j\in J}{\mathop \sum }\, {{u}_{j}}~{{\beta }_{j}}\left(\omega \right)~\left(\varphi \left(\omega \right)~y_{j}^{+}\left(\omega \right)+y_{j}^{-}\left(\omega \right) \right)$ are described. The number of susceptible population members is equal to the sum of ineffective vaccinations $\varphi \left(\omega \right)~y_{j}^{+}\left(\omega \right)$ and unvaccinated individuals $y_{j}^{-}\left(\omega \right)$. Thus, the infection cost corresponds to the unit cost by multiplying the transmission rate and the susceptible population. The disjunctive constraints (2) and (3) are set to determine whether the ordering quantity exceeds the arrival in case of the regular vaccination. The population p_j is used as an upper bound to tighten the constraints. Constraints (4) and (5) are set to determine the allocation quantity of matched vaccine types in the regular vaccination equal to the minimum value of the ordering quantity and regular arrivals. Constraint (6) specifies the upper bound on the ordering quantity. The quantities of unmet demand and excess inventory are defined by Constraints (7) and (8), respectively. In constraint (9), the ordered vaccines subtracted from overstocks and added unmet demands are equal to the regular arrivals. Constraints (10) and (11) are used to obtain the number of shared adult vaccines with respect to pediatric demands. It is only possible to cover pediatric demand per adult vaccine. Additionally, it is not possible for the number of pediatric patients injected with matched and substituted vaccines to exceed the number of regular arrivals ${{r}_{00}}\left(\omega \right)+{{r}_{10}}\left(\omega \right)\le {{d}_{0}}\left(\omega \right)$. Similarly, constraints (12) to (13) correspond to sharing pediatric vaccines for adult demands. Each adult demand requires two doses of a pediatric vaccine, and thus the sharing quantity is less than or equal to the surplus pediatric vaccines divided by two. Additionally, the total number of injected adults should be less than or equal to the number of arrivals. In constraint (14), the sum of regular and supplementary vaccinations is equal to the total vaccinated population. Constraint (15) involves calculating vaccinated and unvaccinated populations. Constraints (16) and (17) define the non-negative continuous variables in the first and second stages, respectively. Constraint (18) corresponds to the binary variable in the second stage.

This problem corresponds to a two-stage stochastic program with a fixed recourse, while the coefficients associated with the second-stage decision variables of every constraint are fixed. The unfixed parameters are only associated with the objective function and the right-hand-side value of the second-stage constraint. The first-stage problem determines ordering quantities irrespective of the uncertain parameters. Consequently, the second-stage problem determines the allocation decision based on the first-stage solution and the realization in each scenario. The vaccine ordering quantity x_i is assumed as a non-negative continuous variable, and thus the first-stage problem can be easily solved as an LP. The second-stage problem is a mix-integer programming (MIP) problem, in which the variable ${{\pi }_{i}}\left(\omega \right)$ is binary, and the other variables are continuous. The problem size increases rapidly with the number of scenarios and posts a challenge to solve the optimal solution directly. In the later section, an algorithm is introduced to obtain an approximate solution for the two-stage program with continuous random variables.

3.2.   Stochastic programming

This subsection considers the decisions of ordering quantity and allocations when all parameters are deterministic. We apply the expected value of uncertain parameters to construct the deterministic model, known as the expected value problem (EV). We assume that the infection cost for either an adult or a pediatric case is significantly greater than the total cost of vaccination, that is c+v¹, c+v², c+v³$\ll$u_j. Furthermore, the cost of injecting a matched vaccine is less than that of using substituted vaccines v¹ < v². This assumption considers the additional costs of labor and syringes to reconstitute the substituted vaccine. Finally, the cost of the supplementary vaccination exceeds that of the regular vaccination v² < v³, and this assumption considers the additional expenditures of transportation and personnel to implement the supplementary vaccination. The sharing quantities are bound by the number of excess doses after the inventories are consumed by the regular arrivals as expressed in constraints (10) and (12). To minimize the total cost, the optimal allocation of using substituted vaccine is established in Proposition 1.

Proposition 1. The optimal sharing quantity is zero for the deterministic problem, i.e., $r_{ij}^{*} = 0, ~\forall i\in I, \forall i\in I, ~i\ne j$.

The optimal allocation quantities can be found in Theorem 1 as a consequence of the zero sharing quantity. Additionally, Theorem 1 is used to establish the optimal ordering quantity in Theorem 2.

Theorem 1. If ${{c}_{i}}+{{v}^{1}}+{{u}_{j}}~\varphi ~{{\beta }_{j}} < {{u}_{j}}~{{\beta }_{j}}$, then the optimal allocation quantity in regular vaccination is d_j; otherwise, 0. If $c+{{v}^{3}}+{{u}_{j}}~\varphi ~{{\beta }_{j}} < {{u}_{j}}~{{\beta }_{j}}$, then the optimal allocation quantity of supplementary vaccination is p_j-d_j; otherwise, 0.

Proof. By Proposition 1, only the allocation quantities using the matched vaccine in a regular vaccination (${{r}_{ij}}, ~\forall i\in I, ~\forall j\in J, j = i$) and a supplementary vaccination (s_j)) are considered. If the infected cost of the unvaccinated population exceeds the total cost of the regular vaccination, then the optimal allocation quantity involves serving all arrival recipients. In contrast, there is no regular vaccination if the unvaccinated cost is cheaper than the vaccination. Additionally, if the cost of the supplementary vaccination is less than the unvaccinated cost, then the lowest cost is associated with covering the overall population.

Theorem 2. The optimal ordering quantity for the deterministic problem is:

Proof. Using Proposition 1 and Theorem 1, the optimal ordering quantity is easy to argue, which is exactly equal to the optimal allocations in various circumstances. Thus, there are no excess orders in the optimal solution. In the first condition, the optimal ordering quantity is equal to the regular vaccinations. In the second condition, the ordering quantity is equal to the population if the total cost of the supplementary vaccination is cheaper than the infection cost of the unvaccinated population. Finally, the ordering quantity is zero if the infection cost of the unvaccinated population is less than the cost of the vaccinated population with respect to any vaccination activity.

3.3.   The sample average approximation with MIP warm start

This study assumes that the transmission rate, vaccine efficacy, and arrivals are continuous uniform distributions. The sample-average approximation (SAA) is applied to generate a finite set of scenarios from the underlying distribution for constructing the stochastic program ^[7]. Furthermore, our algorithm reuses the optimal solution of the previous iteration and then solves the stochastic integer program in an incremental manner. We utilized the mathematical programming language to construct the RP problem, where a general model was formulated beforehand, along with the sampled parameters to populate mathematical programs in each iteration. The procedure of the proposed algorithm is as follows.

4.   Results

The algorithm is developed in C# using the CPLEX 12.6 callable library. Our experiments are performed on an Intel i5 CPU with a processing speed of 2.3 GHz and a 12 GB memory. The following sections first describe the data collection and model parameters. Subsequently, the optimal solutions for deterministic and stochastic problems are presented in conjunction with the computational result of solving the recourse problem by using the proposed algorithm. A comparison of deterministic and stochastic solutions is presented as insights of the value of stochastic solutions.

4.1.   Data and parameters

We consider the following uncertain factors: vaccine efficacy, attack rate, and arrival demand. Each parameter is constructed as a continuous uniform distribution. For setting the lower and upper bounds of the uniform distribution for vaccine efficacy, we refer to the prior literature, which reports multi-year vaccine efficacies in different countries based on randomized controlled trials ^[23]. Additionally, the review study includes findings from a wide spectrum of laboratory experiments conducted in different time periods and population groups. We use the range of 50% to 60% for vaccine efficiency to cover most likelihood cases. The attack rate is the number of newly infected cases during a year divided by the population at the start of the year. Based on WHO reports, the global annual attack rates are estimated as being between 5% to 10% for adults and 20% to 30% for children ^[2]. Vaccination demands were volatile in populations. For example, in the period 2012–2015, the coverage rates ranged between 40.79% to 42.66% for the population of senior citizens, that for children aged 4 years to 6 years ranged between 39% and 49%, and that for children aged 6 months to 3 years ranged from 61% to 65% ^[24]. For our model setting, the demand is based on the actual vaccination coverage rate of each age-group population during the 2014–2015 influenza season.

The vaccine cost is based on the 2015 vaccine budget document published by the Ministry of Health and Welfare of the Republic of China ^[25]. Note that the cost presented in the following sections is in New Taiwan dollars. According to the report, the unit cost of the influenza vaccine corresponded ＄130 per dose. We set the vaccination cost to ＄300 per person to cover the expenses on the care service and syringes. The information regarding the exact cost of treating the infected individual is unavailable. We thus assume that costs are ＄20,000 for each adult and ＄10,000 for each pediatric.

4.2.   The deterministic and stochastic solutions

The expected value solution involves constructing a deterministic programming model by replacing unfixed parameters with expected values. For the problem in the present study, the average vaccine efficacy corresponds to 55%, the average adult attack rate corresponds to 8.6%, the average pediatric attack rate corresponds to 34.5%, the average arrival demand of adult cases corresponds to 2,674,181, and the average arrival demand of pediatric cases corresponds to 252,126. The RP considers each realization from the sampling. The fixed parameters for both EV and RP are common. Figure 1 illustrates the optimal objective values of EV and RP in different sample sizes. The objective value of EV (denoted by the dashed line) is approximately 8.4 million. The RP solution values vary based on the sample sizes, and the objective value tends to become stationary at 8.9 million with respect to sample sizes beyond three thousand.

A decision that uses the expected values can be careless and obtain a lessening performance in scenarios deviated from the average. The expectation of the expected value solution (EEV) is a frequently-used measure to provide insights on the robustness of the EV solution. To obtain EEV, the first-stage decision variable is fixed as the EV solution, and this is followed by solving for the second-stage variables in each realization. The EEV considers the expectation of the objective values for all realizations. In the following figure, each dot represents the optimal objective value of using the EV solution in a realization. The objective values range between 6 million and 12 million. The EEV corresponds to $8,899,142.

Figure 3 illustrates the computational result of using the proposed algorithm to solve the RP in different sample sizes. The criterion that involves terminating the branch-and-bound procedure is set as 1%. The runtime continues to increase as the problem incorporates higher sample sizes. Most problems can be solved within 20 s. However, in a few circumstances, the runtime approximately corresponds in the range of 60–70 seconds.

Table 2 presents the optimal ordering and allocation quantities for the RP and EV solutions. The EV solution orders amount to 2,674,181 doses of adult vaccines and 252,126 doses of pediatric vaccines. All adult vaccines are administered to adults in regular vaccination. Similarly, pediatric vaccines are only administered for pediatric cases. There is no sharing of vaccines and supplementary vaccination. With respect to the RP solution, the first-stage decision variable values from the problem with the largest sample size are reported while 2,406,959 adult vaccines and 226,919 pediatric vaccines are ordered. The second stage decisions are varying with the scenarios.

Table 3 reports the objective values of RP and EV. Additionally, the expectation of objective values of using the EV solution is calculated. For general stochastic programs, the relationship between RP and EEV can be established as the following inequality RP ≤ EEV. With respect to the minimization problem, RP obtains a lower objective value because the solution considers all possible scenarios. However, the EEV only makes the decision once by referring to the average values. A detailed proof can be found in the study by Madansky in 1960 ^[26]. With respect to the problem outlined in the present study, the RP is less than EEV (＄8,884,969 < ＄8,899,142), and this result satisfies the theoretical condition. Next, the gap between EEV and RP is measured, and this is termed the value of the stochastic solution (VSS), where VSS = EEV - RP. Evidently, VSS is non-negative because RP is always less than or equal to EEV with respect to minimization problems. The gap represents the potential loss of using the expected value to make a decision. The objective values and VSS are summarized in the following table, where VSS is significantly less than then EEV. This is because the uncertain parameters have minor variations. In other cases of highly variant, RP would be much desirable give it would reduce more costs from the EEV solution.

5.   Discussion and conclusions

This study examines the important challenge in vaccine supply chains to aid public health decisions to control influenza. Ordering and allocation decisions are integrated into a vaccine supply chain by considering uncertain infection, vaccine efficacy, and recipient arrival. A two-stage stochastic programming model is proposed to determine the optimal vaccine ordering quantity and allocation. Constraints in the stochastic programming model were formulated according to the WHO guidelines on vaccine characteristics and usages. Our model can be used for not only analyzing vaccine policies in Taiwan but also for other countries against influenza pandemics.

We present an algorithm based on Monte Carlo simulation with MIP warm-start to solve the SIP with continuous random variables. The algorithm is not limited to any probability distributions. Furthermore, the optimal solution for the deterministic problem is analyzed and directly solved as opposed to using a conventional solver to obtain the optimal solution. The uncertain effect is investigated by comparing the deterministic and stochastic solutions in realizations. The result indicated that an additional cost of $14,143 incurred if policymakers ignore uncertain effects by simply determining a decision based on the expected values (Table 3).

Assumptions and limitations in this study state as follows. We assume an unlimited vaccine supply. Additionally, the vaccine ordering quantity is not bound by a budget. However, these two assumptions can be relaxed by adding additional constraints. Another assumption is that the adult and pediatric vaccines can apply to all age groups when the vaccine inventory is in short. In reality, concerns exist with respect to the substituted vaccine in populations. Additionally, the supplementary vaccination is assumed to cover all target populations in the study, while this may not be applicable in real-world scenarios. Finally, the second-dose requirement is ignored in the case of first-time vaccine recipients under eight years of age, and this may impact both ordering and allocation decisions. Both deterministic and stochastic models consider the influenza intervention by incorporating the actual vaccination activities in Taiwan. Despite the variation in vaccine usages in different countries, the proposed approach can be generalized and reused with minor modifications.

Acknowledgments

This study is supported by the Ministry of Science and Technology, R.O.C. under grants MOST 108-2221-E-009-029. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Conflict of interest

All authors declare no conflicts of interest in this paper.