# Design of Targeted Community-Based Resource Allocation in the Presence of
Vaccine Hesitancy via a Data-Driven Compartmental Stochastic Optimization
Model

Hieu Bui
Department of Industrial Engineering
University of Arkansas
Fayetteville, AR 72701
&Sandra Ekşioğlu
Department of Industrial Engineering
University of Arkansas
Fayetteville, AR 72701
&Rubén Proaño
Department of Industrial and Systems Engineering
Rochester Institute of Technology
Rochester, NY 14623
&Haoming Shen
Department of Industrial Engineering
University of Arkansas
Fayetteville, AR 72701
Corresponding author

###### Abstract

Vaccines have proven effective in mitigating the threat of severe infections
and deaths during outbreaks of infectious diseases. However, the reluctance or
refusal to get vaccinated, known as vaccine hesitancy (VH), presents a
significant challenge in predicting the spread of the disease and assessing
the need for healthcare resources in different regions and population groups.
We propose a modeling framework that integrates an epidemiological
compartmental model that captures the spread of an infectious disease within a
multi-stage stochastic program (MSP) that determines the allocation of
critical resources under uncertainty. The proposed compartmental MSP model
adaptively manages the allocation of resources to account for changes in
population behavior toward vaccines (i.e., variability in VH), the unique
patterns of disease spread, and the availability of healthcare resources over
time and space. The compartmental MSP model allowed us to analyze the price of
fairness in resource allocation. We developed a case study using real-life
data about COVID-19 vaccination uptake from January to May 2021 and data about
healthcare resources in Arkansas, U.S. Using our numerical analysis, we
observed that ($i$) delaying the initial deployment of additional ventilators
(critical resource) by one month could lead to an average increase in the
expected number of deaths by 285.41/month, highlighting the importance of
prompt action; ($ii$) each additional ventilator in the initial stockpile and
in supply lead to a decrease in the expected number of deaths by 1.09/month
and 0.962/month, respectively, highlighting the importance of maintaining a
large stockpile and scalable production response; ($iii$) the price of
ensuring equitable allocation of resources changes over time and space. This
price is highest during the peak of a disease outbreak and in highly populated
regions. In summary, this study highlights the importance of flexible and
informed decision-making in public health crises, and the crucial role of
advanced planning and preparedness for future outbreaks, offering a framework
for effective resource allocation adaptable to various public health
emergencies.

_Keywords_ compartmental model in epidemiology, compartmental multi-stage
stochastic program, decision making under uncertainty, resource allocation,
vaccine hesitancy

##  1 Introduction

Infectious disease outbreaks with high reproduction numbers ($R_{O}$), such as
Measles ($R_{O}=16$), COVID-19 ($R_{0}\in[1.35,2.11]$), etc., can affect a
large portion of the population in a short period of time (Sy et al., 2021).
Therefore, the implementation and effective execution of public health
strategies, such as vaccination campaigns, mask mandates, and social
distancing, can greatly reduce the spread of the disease and save lives
(Moghadas et al., 2021; Brooks and Butler, 2021; Chu et al., 2020). However,
the success of vaccination campaigns, which are crucial to controlling disease
outbreaks, is often impeded by vaccine hesitancy (VH). This reluctance or
refusal to get vaccinated leads to a disparity in vaccination rates and
prolongs the repercussions of the disease, as in the case of the COVID-19
pandemic. Furthermore, VH complicates the prediction of vaccine demand and the
allocation of critical resources (Blasioli et al., 2023).

Making well-informed and timely decisions is essential to effectively manage
disease outbreaks. Compartmentalized epidemiological models, like the
susceptible-infectious-recovered (SIR) and susceptible-exposed-infectious-
recovered (SEIR), play a critical role in forecasting the progression of
infectious diseases, thus aiding in the design of containment strategies. The
SIR model divides the population into three groups: Susceptible, Infectious,
and Recovered, offering a basic framework for the dynamics of the disease
(Kermack et al., 1927). The SEIR model expands on this by including an
‘Exposed’ group, which is crucial for understanding the progression of
diseases that have a latency period, such as COVID-19. These models provide
insight into disease transmission patterns that are vital to planning and
implementing public health responses (Kong et al., 2022; Xiang et al., 2021).
The traditional SIR and SEIR models typically assume constant transition rates
between compartments, such as infection and recovery rates, yet this
simplification does not always hold true. Advanced models address this by
incorporating the uncertainty of these model parameters. However, even these
adaptations of the SEIR model that utilize stochastic compartments to inform
COVID-19 strategies often neglect critical aspects of healthcare capacity and
VH (Buonomo et al., 2022; Mesa et al., 2023). Several researchers use
operations research (OR) models to determine resource allocation strategies,
which are critical decisions for healthcare management during disease
outbreaks. These models are used to optimize the distribution of limited
resources, such as ventilators during the COVID-19 pandemic, ensuring
effective utilization (Mehrotra et al., 2020; Bhavani et al., 2021; Koonin et
al., 2020; Bertsimas et al., 2021). From the wide range of OR models,
stochastic programming (SP) and multi-stage stochastic program (MSP) provide a
robust framework for adaptive decision-making under uncertainty (Birge and
Louveaux, 2011).

Incorporating dynamic disease progression data from compartmental models would
significantly enhance the effectiveness of both the SP and MSP models.
However, the integration of these models has been impeded by the computational
challenges associated with solving each of them. Only a few studies have
attempted to integrate these models to enhance the quality of resource
allocation decisions (Yin and Buyuktahtakin, 2021a; Yin et al., 2023). By
incorporating the predictive power of compartmental models into SP and MSP
models, we can identify and evaluate resource allocation strategies that are
informed by the specifics of disease spread and progression and the
availability of healthcare resources. This realization leads to our first
research question (RQ1): “Can we develop an effective resource allocation
model of critical resources needed during disease outbreaks?” To demonstrate
how we can tackle this, we develop an MSP compartmental model for the COVID-19
pandemic. We first develop a compartmental model that is an extension of the
SEIR model to capture the specific dynamics of the spread of the COVID-19
virus, the impact of VH on the number of people exposed, and the availability
of healthcare resources on the expected number of people recovered and of
deaths. The outcomes of the modified SEIR model are incorporated in an MSP
model that determines resource allocation.

During disease outbreaks and other health emergencies, resource allocation
becomes a critical concern. For example, the COVID-19 pandemic caused an
urgent demand for essential medical equipment, such as ventilators for people
in critical condition and personal protective equipment. Recognizing the
gravity of the situation, the U.S. government invoked the Defense Production
Act (DPA), urging car manufacturers to pivot their production lines to
manufacture ventilators (Albergotti and Siddiqui, 2020). This proactive
measure showcased the influence of government interventions in crisis
management. However, it also prompts reflection on its execution and raises
questions such as: Could different results have been achieved with an earlier
or more intensive intervention? As we anticipate and prepare for potential
future outbreaks, it becomes imperative to assess the timing and intensity of
interventions on healthcare outcomes. This leads to the following research
question (RQ2): “What should be effective resource allocation strategies at
the start of an outbreak?” Specifically, how does the timing and scope of
interventions, such as DPA, impact the spread of the disease? How does VH
impact the effectiveness of resource allocation? By exploring this research
question, our aim is to provide guidance on optimizing resource allocation
strategies to ensure that they are timely and effective.

Ensuring access to healthcare resources is crucial, especially during disease
outbreaks, since strategies used for resource allocation can have significant
impacts on the affected population. Several OR models take a utilitarian
approach by focusing on minimizing costs or maximizing the benefits from the
allocation of limited healthcare resources. Such an approach allocates
resources to highly populated areas and those that have access to healthcare
care, which may lead to inequity. In the U.S., access to healthcare varies
greatly, and rural and marginalized communities often face more challenges in
accessing the necessary services and resources (Mishra et al., 2021). These
disparities can affect the effectiveness and success of public health
intervention strategies. This leads to our next research question (RQ3): “What
are the trade-offs between equity, equality, and the effectiveness of resource
allocation?” Exploring this question will help to understand how resource
allocation strategies influence disease outbreak management, both in the short
and long term.

##  2 Literature review

We identified three streams of literature related to our research:
compartmental models in epidemiology, stochastic optimization models for
resource allocation, and modeling equity in healthcare. The first stream of
literature focuses on compartmental models, such as the SIR and SEIR. These
models categorize populations into different subgroups according to their
health status. These models provide a structured approach to assess disease
progression and assess the impact of interventions on the affected population.
The second stream of literature focuses on multistage stochastic programs
(MSPs), models that aid in making decisions under uncertainty. Our study
integrates compartmental models within an MSP framework to optimize resource
allocation to manage and contain infectious diseases. The third stream
addresses equity in resource allocation. The goal is to ensure a fair and
accessible distribution of healthcare resources, which is vital during
infectious disease outbreaks. Each of the following subsections provides a
comprehensive overview of the literature. We highlight their relevance to our
proposed research model.

###  2.1 Compartmental Models in Epidemiology

Epidemiological models are key to understanding and predicting the course of
infectious disease outbreaks. With the emergence of COVID-19, the field of
mathematical epidemiology quickly attracted broader interest. Taking advantage
of existing approaches, researchers adapted and applied various methods
specifically to study COVID-19, including compartmental models, structured
metapopulations, agent-based networks, deep learning, and complex networks (Wu
et al., 2020; Adiga et al., 2020; Rodríguez et al., 2021; Zhou et al., 2020a;
Chang et al., 2021a; Kerr et al., 2021; Chang et al., 2021b). Among these,
compartmental models are essential in infectious disease research due to their
effective representation of the dynamics of disease transmission. In these
models, the population is divided into several subgroups or compartments, each
representing a specific stage in the progression of the disease. Consider, for
example, the SIR model that divides the population into three compartments:
Susceptible, Infectious, and Recovered. This model provides a basic framework
for tracking disease spread by monitoring changes in population proportions
within these compartments. Building on the SIR model, the SEIR model adds an
Exposed compartment, highlighting a critical phase where individuals have
contracted the disease but are not yet infectious. This addition is
particularly important for diseases that have an incubation period, such as
COVID-19. Several researchers extended the SEIR model to capture the unique
transmission dynamics of the COVID-19 virus, as well as the impacts of
vaccination and public health mandates on the spread of the disease. Kong et
al. (2022) provide a detailed overview of these extensions, which consider
various stages of disease progression, intervention strategies, and
demographic factors (Kimathi et al., 2021; Tuite et al., 2020; Foy et al.,
2021).

Classic epidemiological models assume that model parameters, such as the
transmission rate, the recovery rate, and the contact rate, are constant
throughout the course of the epidemic being modeled. Most recently,
epidemiological models have been updated to incorporate the variability and
uncertainty of these parameters. Models that take into account the randomness
that exists in the system can lead to better decision-making. For example,
(Lekone and Finkenstädt, 2006; Gray et al., 2011; Faranda and Alberti, 2020;
Ganyani et al., 2021; Maltsev and Stern, 2021; Gatto et al., 2020; Giordano et
al., 2020; Kretzschmar et al., 2020; Kucharski et al., 2020) propose
stochastic epidemiological models. These models use historical data to
estimate how parameters change over time, and thus, provide a more realistic
representation of the dynamics of infectious diseases.

VH has consistently presented a significant challenge to the healthcare
system, a concern that has been recognized from the World Health Organization
(WHO) as a global health threat ((WHO, 2019). The advent of COVID-19 vaccines
once again highlighted this issue. Several compartmental models capture the
impact of vaccine availability and hesitancy on transmission rate. Buonomo et
al. (2022) demonstrate that voluntary vaccination alone cannot stop the spread
of the disease. Choi and Shim (2020) develop a game-theoretic epidemiological
model to explore vaccination and social distancing as strategies to control
the spread of the disease. They use a group behavioral model to determine
optimal strategies for individuals. Awad et al. (2021) propose a compartmental
model that uses geospatial and VH information to develop strategies to
distribute vaccines. Jing et al. (2023) examine the role of VH in the
emergence of new SARS-CoV-2 variants. Bui et al. (2024) present a
compartmental model that captures the impact of VH on disease outcomes, such
as the expected number of hospitalizations and deaths. The model runs for
different values of VH (which are taken from a certain distribution) to
estimate the impacts of uncertain VH on model outcomes. The authors performed
a sensitivity analysis to assess the impact of the availability of healthcare
resources and vaccines on the expected number of deaths. Despite the variety
of compartmental models available, none examines the impacts of the temporal
and spatial variations of VH on demand for healthcare resources. This gap
underscores the need for models that use predictions of temporal and spatial
disease spread to determine the effective allocation of critical resources.

Based on our review of the literature, epidemiological models are typically
used to support decision making via “What-if" analysis. This approach explores
how changes in model parameters or assumptions affect the spread of the
disease and the need for resources. Some researchers integrate epidemiological
models into optimization models to support optimal decision making during
disease outbreaks. For example, Abdin et al. (2023) propose a non-linear
programming formulation of a compartmental model to present the complex
dynamics of infectious diseases, offering a nuanced perspective on disease
transmission. Yin and Buyuktahtakin (2021a) propose an MSP compartmental model
to support resource allocation during the Ebola outbreak in West Africa. The
MSP uses information on stochastic disease progression to determine resource
allocation to control the infectious disease outbreak. The model determines
the fair allocation of critical resources.

Our study investigates the impact of vaccination and VH on disease
transmission dynamics. The model we propose enables decision makers to assess
the impact that the timing, scale, and scope of resource allocation have on
the spread of the disease and disease outcomes. While previous research has
addressed aspects of resource management during disease outbreaks, none has
presented a comprehensive, data-driven model, as we propose. Our model
highlights the advantages of merging the predictive capabilities of
compartmental models with optimization techniques, providing valuable insights
into effective disease control strategies.

###  2.2 Stochastic Optimization Models for Resource Allocation

MSP models are used to support sequential decision making under uncertainty.
The simplest MSP model, the two-stage stochastic program (2-SP), is commonly
used to minimize expected costs (maximize expected benefits) in decisions made
today and tomorrow in the face of uncertainties faced tomorrow. MSPs are an
extension of 2-SPs, which involves making decisions sequentially based on new
information that emerges over time. For an overview of SP and MSP models and
their applications, we refer the reader to Pereira and Pinto (1991); Goh et
al. (2007); Zhang and Cardin (2017); Li and Grossmann (2021).

MSP models pose significant computational challenges that stem from the
inherent uncertainty in the data and the intricate nested structure
characteristic of multistage decision-making problems. A common approach to
handling uncertainty is to approximate the underlying stochastic process using
a scenario tree. A scenario tree branches at each stage to consider a finite
number of possible scenarios. When the number of scenarios is small, the
deterministic equivalent of the MSP problem can be solved using commercial
optimization solvers, such as CPLEX or Gurobi. Lagrangian and Benders
decomposition methods are often used to solve problems with a large number of
scenarios. These methods break down the problem into smaller subproblems that
are easier to solve (Laporte and Louveaux, 1993; Guignard, 2003). Rolling-
horizon and heuristic approaches are also used to solve larger-scale MSP
models (Christian and Cremaschi, 2015).

MSP models are often used for resource allocation in the healthcare sector, as
the size and timing of the demand for these resources are often uncertain, and
the availability of resources is limited. For example, Yin and Buyuktahtakin
(2021b) develop an MSP model to effectively manage resource allocation during
the 2018-2020 Ebola outbreak. In response to the COVID-19 pandemic, Yin et al.
(2023) propose an MSP model for allocating ventilators across counties in New
York and New Jersey. Their MSP incorporates a compartmental model to account
for disease dynamics under varying levels of uncertainty. Hosseini-Motlagh et
al. (2023) introduce a multi-stage fuzzy SP approach for patient allocation
across healthcare facilities during COVID-19. Tanner et al. (2008) and Yarmand
et al. (2014) propose SP models to determine optimal vaccination strategies.
Mehrotra et al. (2020) propose an SP model for ventilator allocation during
COVID-19, focusing on cost minimization. These publications highlight the
relevance of MSPs to support sequential decisions related to resource
allocation under uncertainty.

Integrating compartmental models within MSP resource allocation models
provides a robust framework for planning responses to disease outbreaks. This
approach enables high-quality decision making over time, taking into account
the unpredictable nature of disease spread and resource availability. While
some researchers are investigating this modeling approach, there is more to be
done. Our proposed data-driven compartmental MSP model addresses VH and
fairness in resource allocation within this modeling framework. Our proposed
model supports effective vaccination strategies and resource allocations that
ensure greater immunity and mitigate the impact of outbreaks.

###  2.3 Equity in Healthcare Resource Allocation

A WHO report highlights the importance of following ethical principles in the
allocation of healthcare resources by analyzing the utility, efficiency, and
fairness approaches taken by decision makers (WHO, 2004). Utilitarian decision
makers implement strategies that yield the highest overall benefit or the
lowest overall cost. Efficient decision makers advocate using the least amount
of resources to achieve the maximum possible impact. Fair decision makers aim
to ensure equitable allocation and reduce disparities between groups and
individuals. Aligning resource allocation strategies with these principles can
be challenging as each strategy, while meeting one principle, falls short of
meeting others. The work of Lane et al. (2017) highlights the complexity of
ensuring equity in allocating healthcare resources, emphasizing the need for
clear definitions that reflect societal values to ensure equal access to
healthcare, address patient needs effectively, and distribute potential
benefits fairly. Equitable strategies can build public trust and encourage
cooperation, which are essential during health crises.

Several OR models for resource allocation take a utilitarian approach to
optimize logistics and operations to control disease spread (Sun et al., 2014;
Liu et al., 2015; Liu and Zhang, 2016; Zaric and Brandeau, 2001; Day et al.,
2020). Other studies introduce fairness metrics and propose models that
optimize fairness in resource allocation. For example, Orgut et al. (2016)
proposes a model that allocates and distributes food donations proportionally
to demand, aiming for an equitable distribution while minimizing waste.
Similarly, Davis et al. (2015) propose a multi-step optimization approach to
improve fairness in the allocation of organs (kidneys) across the U.S.,
thereby improving equity of organ transplant. Dönmez et al. (2022) presents a
multi-objective, multi-period, nonlinear model for equitable personal
protective equipment (PPE) allocation to health centers in pandemics, with the
aim of minimizing infections and ensuring fair distribution amidst resource
scarcity. The research by Yin and Buyuktahtakin (2021a) introduces an MSP
model for a fair allocation of resources, such as beds and treatment centers,
during the West African Ebola outbreak, balancing efficiency and equity in
epidemic responses.

Similarly to the literature, our proposed model allows decision makers to
evaluate the trade-offs between fairness and efficiency in allocating critical
healthcare resources. Unique to our approach is the incorporation of
uncertainty in VH, which introduces additional complexity to resource
allocation strategies. We analyze the ‘price of fairness’, the efficiency loss
incurred when striving for equitable resource allocation over time and across
different geographic regions. The MSP compartmental model adapts to the
evolving landscape of healthcare needs, ensuring that resource allocation
remains responsive and equitable even under unpredictable changes in public
sentiment towards vaccination.

##  3 Method

We propose a compartmental MSP model to evaluate how uncertainty in VH, and
disease progression over time and space impact decisions about allocating
healthcare resources. Section 3.1 focuses on the generation of scenarios used
in the compartmental MSP model. Detailed descriptions of the compartmental
model and the compartmental MSP model are presented in Sections 3.2 and 3.3,
respectively.

###  3.1 _Uncertainty representation_

Motivation: Our proposed model considers unccertainty in VH. This is motivated
by observations made from the COVID-19 Trends and Impact Survey (CTIS)
conducted by the Delphi Group at Carnegie Mellon University, in collaboration
with Facebook (Salomon et al., 2021), and from data about vaccination rates
collected by the Centers for Disease Control and Prevention CDC (2021).

Figure 1(a) presents VH estimates of every county in the states of Missouri
(MO) and Rhode Island (RI), providing a snapshot of the willingness of the
public to vaccinate. The VH estimate is derived from responses to the
following question in the CTIS survey: “If a vaccine to prevent COVID-19 was
offered to you today, would you choose to get vaccinated?”. Participants were
asked to select one of the following answers “Yes, definitely,” “Yes,
probably,” “No, probably not,” or “No, definitely not.” The last three
responses are considered indicative of VH. Each thin line represents VH
progression of an individual county. The thicker line depicts the average VH
trend for the respective state.

Figure 1(b) shows the vaccination rates in counties of Missouri and Rhode
Island based on data from CDC (2021). The data show fully vaccinated
individuals who have received the complete dose - either two doses of a two-
dose vaccine or one dose of a single-dose vaccine. These data highlight a
general decline in VH over time and underscore the varied vaccination
sentiments in counties. Figure 1(c) illustrates a negative correlation between
average vaccination rates and average VH estimates with correlation
coefficients of $\rho=-0.85$ for Rhode Island and $\rho=-0.93$ for Missouri.
The high correlation coefficients in both states indicate that higher
vaccination rates are associated with lower VH.

Figure 1: (a) CTIS-derived vaccine hesitancy, (b) vaccination uptake trends,
and (c) correlation plot between VH and vaccination rate of counties in Rhode
Island (RI) and Missouri (MO).

Within a state, varying attitudes towards vaccination and differing
vaccination rates significantly impact disease transmission and the resulting
need for healthcare resources. We capture variations in VH over time and space
via a scenario tree that we describe next.

Scenario tree: Scenario trees provide a systematic approach to represent
uncertainty that is suitable for quantitative models. They provide likely
scenarios of the future with associated probabilities.

Each scenario $\omega$ ($\omega\in\Omega$) represents the changes of VH for a
specific population group, over time. We derive temporal changes in VH from
the historical data provided by the CTIS. We fit a distribution to the rate of
change of VH. Numerical analysis shows that the (monthly) rate of change of VH
is normally distributed. The corresponding mean ($\mu_{t}$) and standard
deviation ($\sigma_{t}$) change from one period (month) to the next (see
Figure 12). To construct a decision tree of manageable size, we approximate
the normal distribution with a discrete triangular distribution that assumes
the values $\mu_{t},\mu_{t}-\sigma_{t},\mu_{t}+\sigma_{t}$ with probabilities
of 0.158, 0.684, and 0.158, respectively.

Figure 2 presents the scenario tree for a problem with $T=3$ stages. The
resulting scenario tree has $3^{2}=9$ unique scenarios. Let $p_{\omega}$
denote the probability associated with each of these scenarios. The tree
showcases a baseline (green path), an optimistic (blue path), and a
pessimistic (orange path) scenario. In the optimistic scenario, the VH changes
by $\mu_{t}-\sigma_{t}$ in each stage. In the pessimistic scenario, the VH
changes by $\mu_{t}+\sigma_{t}$ at each stage. In the baseline scenario, the
VH changes by $\mu_{t}$ at each stage. Consequently, the probability
associated with the pessimistic and optimistic scenarios is $0.158^{2}=0.024$.
The probability associated with the baseline scenarios is $0.684^{2}=0.468$.

Figure 2: A scenario tree illustrating three stages, with each stage branching
into three paths, resulting in nine unique scenarios. Branch annotations
denote the rate of change and corresponding probabilities.

###  3.2 _A compartmental model of COVID-19_

The SEIR model divides the population into four groups (compartments):
Susceptible ($\boldsymbol{S}$), Exposed ($\boldsymbol{E}$), Infectious
($\boldsymbol{I}$), and Recovered ($\boldsymbol{R}$). The model captures the
transition of individuals through these compartments over time. It also
accounts for the latency period before an individual becomes infectious, which
is important in modeling the spread of diseases such as COVID-19. In a
previous study, our team extended this model to address the specific
complexities of COVID-19, such as vaccine availability, VH, and the scarcity
of critical healthcare resources in the spread of the disease (Bui et al.,
2024). For the convenience of the reader, we summarize the model in the
following paragraphs.

A schematic representation of the proposed susceptible-vaccinated-exposed-
infectious-hospitalized-recovered (SVEIHR) model is shown in Figure 3. In this
model, susceptible individuals can be exposed or opt to vaccinate, moving
forward to the vaccinated compartment (V). The size of compartment V is
dictated by $\rho$, representing the maximum vaccination rate achievable in
the absence of VH. The willingness of people to get vaccinated is represented
by $(1-\tilde{h})$, where, $\tilde{h}$ is a stochastic parameter that
represents the proportion of the population reluctant to vaccinate.
Furthermore, the model considers that vaccination does not confer absolute
immunity; therefore, vaccinated individuals can progress to an Exposed but
Vaccinated compartment ($\boldsymbol{EV}$), with a probability determined by
vaccine efficacy, $\epsilon$.

Figure 3: SVEIHR: COVID-19 compartmental model that incorporates vaccination
dynamics, VH, and healthcare resource constraints.

The model captures different postexposure health trajectories based on
incubation rates $\alpha$ and the severity of the disease. Individuals
progress to Mild ($\boldsymbol{I^{\text{m}}}$), Severe
($\boldsymbol{I^{\text{s}}}$), or Critical Condition
($\boldsymbol{H^{\text{c}}}$) compartments with probabilities
$p^{\text{m}},p^{\text{s}},$ and $p^{\text{c}}$, respectively. Vaccinated
individuals are less likely to experience severe/critical health outcomes;
thus, several infected individuals are likely to recover and move to the
compartment ($\boldsymbol{R}$).

The model factors the availability of critical healthcare resources, such as
hospital beds (Bed) and ventilators (Vent), in determining patients’ health
outcomes. Refer to the example presented in Appendix A for a demonstration of
the impact of the number of beds and ventilators on the model outcomes.
Recovery ($\gamma$) and mortality ($\mu$) rates vary depending on the severity
of the disease and the availability of healthcare resources. The model
accounts for the transition of critically ill patients to the compartment
($\boldsymbol{D}$) when ventilators are unavailable, emphasizing the dire
consequences of resource limitations during a pandemic.

###  3.3 _A compartmental multi-stage stochastic optimization model_

This section outlines the formulation of the proposed multi-stage stochastic
program (MSP) model. The model determines an optimal sequence of decisions
related to ventilator allocation across various regions over time. The MSP
minimizes the total number of deaths by dynamically adjusting the distribution
of healthcare resources in response to changing VH rates and state transitions
within the SVEIHR model. Resource allocation decisions take into account not
just immediate health outcomes, but also the need for equitable access to
care. By integrating the compartmental model with the MSP model, we now better
understand how VH and disease dynamics affects healthcare systems, and use
this information to optimize the allocation of ventilators. We refer to this
proposed model as the compartmental MSP model. For a detailed description of
the model notations used, please see Appendix C.

Model assumptions: Our model operates on three primary sets of assumptions
that are critical to its structure and function:

_Assumption 1: Uncertainty in vaccine hesitancy (VH):_ The model considers
three potential outcomes of VH in each decision stage. Let ${h}_{t-1}$
represent the proportion of the population reluctant to vaccinate in period
$t-1$. Then, in period $t$, the proportion of the population reluctant to
vaccinate for the optimistic, expected, and pessimistic scenarios are
$h_{t-1}(1+\mu_{t}-\sigma_{t});h_{t-1}(1+\mu_{t})$; and
$h_{t-1}(1+\mu_{t}+\sigma_{t})$. Each of these outcomes is associated with a
fixed probability that does not change over the course of the study. These
probabilities remain constant and unaffected by changing decisions. We note
that decision makers can customize these outcomes (number of outcomes and
values) based on available data.

_Assumption 2: Disease dynamics_ : Since SVEIHR is a compartmental
epidemiological model, the following assumptions of SEIR model apply. ($i$)
Every individual in a population is equally likely to interact with others. As
a result, disease progression and transitions between compartments are
consistent throughout the population. However, it should be noted that in a
real-world setting, elements such as social networks and age demographics can
significantly influence disease transmission patterns. ($ii$) The size of the
total population does not change due to natural births, non-COVID-19-related
deaths, or migration. Since the time frame considered in our study spans
several months, recovered individuals are assumed to have full and lasting
immunity to the disease. ($iii$) Infection, recovery, and exposure rates are
assumed to be deterministic. However, we recognize that in practice, they can
fluctuate for many reasons, from policy interventions and demographic
differences to viral mutations.

_Assumption 3: Healthcare system:_ The available healthcare resources (i.e.,
beds, ventilators) are only used to treat COVID-19 patients. Furthermore, our
proposed MSP model does not account for potential reallocations of healthcare
resources.

Mathematical model formulation: Classical SIR are continuous-time models.
However, several researchers propose discrete-time approximations of these
models that are easier to solve (Wacker and Schlüter, 2020). These studies
show that many of the desired properties of the time-continuous SIR models
remain valid in the time-discrete approximations. We use a discrete-time model
to mimic the dynamic changes in the size of compartments of the SVEIHR model
over time. Let $\mathcal{T}$ represent the set of time periods in the planning
horizon. Let $\mathcal{J}$ represent the set of decision periods (stages).
Note that, between every two consecutive decision periods, there are several
time periods, thus $\mathcal{J}\subset\mathcal{T}$. We denote by $\mathcal{R}$
the set of regions studied. We distinguish between a decision stage $j$ and a
time period $t$ to keep the discrete approximation of the SEIR model
effective, ensuring short time intervals between transitions, while avoiding
too many decision stages for the MSP model.

In a decision period (stage), we determine the allocation of ventilators. Our
decision variables are $x_{\omega,j,r}$, which represent the number of
ventilators assigned to the region $r\in\mathcal{R}$ at the beginning of
period $j\in\mathcal{J}$ in scenario $\omega\in\Omega$. Let
$\mathcal{X}_{\omega,t,r}$ represent the cumulative number of ventilators
allocated to region $r$ by time $t$.

The decisions made affect the state variables, which are the number of
susceptible individuals who are unvaccinated ($S_{\omega,t,r}$) and vaccinated
($V_{\omega,t,r}$); the number of exposed individuals who are unvaccinated
($E_{\omega,t,r}$) and vaccinated ($EV_{\omega,t,r}$); the number of
infectious individuals with mild ($I^{\text{m}}_{\omega,t,r}$) and and severe
($I^{\text{s}}_{\omega,t,r}$) symptoms; the number of hospitalized individuals
in severe ($H^{\text{s}}_{\omega,t,r}$) and in critical
($H^{\text{c}}_{\omega,t,r}$) condition; the number of individuals admitted to
a hospital and in severe ($\mathcal{A}^{\text{s}}_{\omega,t,r}$) or critical
($\mathcal{A}^{\text{c}}_{\omega,t,r}$ ) condition; the number of individuals
in severe ( $\mathcal{K}^{\text{s}}_{\omega,t,r}$) and critical
($\mathcal{K}^{\text{c}}_{\omega,t,r}$) condition not admitted to the hospital
due to capacity constraints; the number of recovered individuals
($R_{\omega,t,r}$); and the number of deaths ($D_{\omega,t,r}$). Let
$\mathcal{L}=\\{S,V,E,EV,I^{\text{m}},I^{\text{s}},H^{\text{s}},H^{\text{c}},R,D\\}$
represent the set of our state variables.

Next, we provide detailed descriptions of the objective function and
constraints of the proposed compartmental MSP.

_Objective function:_ The objective is to minimize the total expected number
of deaths across regions during the planning horizon.

$\mathcal{Z}=\min\sum_{r\in\mathcal{R}}\sum_{\omega\in\Omega}p_{\omega}D_{\omega,|\mathcal{T}|,r}$
(1)

_Constraints:_ We organize our constraints into five groups (G1, G2,…, G5),
each tailored to fulfill a specific function. This structure enhances the
model’s clarity and highlights the essential contribution of each constraint
category to the model’s objective.

(G1:) _SVEIHR compartmental model:_ This group of constraints captures the
dynamic changes of state variables over time. We begin by initializing the
state variables $\ell\in\mathcal{L}$ at $t=0$. Let $\boldsymbol{\pi}_{\ell,r}$
represent the corresponding initial values.

$\ell_{\omega,0,r}=\boldsymbol{\pi}_{\ell,r}\hskip
7.22743pt\forall\ell\in\mathcal{L},r\in\mathcal{R},\omega\in\Omega$ (2)

Constraints (3)-(12) describe the transition of the population between
different compartments of the epidemiological model over time (as illustrated
in Figure 3). Constraint (3) updates the size of the susceptible population in
compartment S. To determine the size of S at the end of period $t+1$, we
update the size of S at the end of period $t$ by subtracting the number of
individuals exposed, vaccinated and subtracting (adding) those who immigrated
(emigrated) to the region from (to) other regions considered in the analysis.
In Appendix F, we provide descriptions of the effects of migration in our
model. Constraint (4) updates the size of the vaccinated population,
compartment V. To determine the size of V at the end of the period $t+1$, we
update the size of V at the end of period $t$ by subtracting the number of
individuals exposed and adding the number of individuals vaccinated. The
efficacy rate of the vaccine, $\epsilon$, affects the size of this
compartment. Constraint (5) updates the size of the unvaccinated population
who were exposed, compartment E. To determine the size of E at the end of
period $t+1$, we update the size of E at the end of period $t$ subtracting the
number of individuals exposed who became infected, and by adding the number of
individuals exposed during period $t+1$. Constraint (6) updates the size of
the vaccinated population who are exposed, compartment EV. This constraint
updates the size of EV by taking a similar approach to (5). The difference
between constraints (5) and (6) is that vaccinated individuals can recover
after exposure and do not become critically ill. Constraints (7) and (8)
update the size of compartments $I^{\text{m}}$ and $I^{\text{c}}$,
respectively. To determine the size of compartment $I^{\text{m}}$ at the end
of $t+1$, constraint (7) considers new infections, and recoveries. To
determine the size of compartment $I^{\text{s}}$ at the end of $t+1$,
constraint (8) considers new infections, recoveries, deaths and transitions to
the hospitalization compartment. Constraints (9) and (10) update the size of
compartments $H^{\text{s}}$ and $H^{\text{c}}$. To determine the number of
severely hospitalized individuals (the size of compartment $H^{\text{s}}$) at
the end of period $t+1$, constraint (9) considers new admissions, recoveries,
and deaths. To determine the number of critically hospitalized individuals
(the size of compartment $H^{\text{c}}$) at the end of period $t+1$,
constraint (10) also considers new admissions, recoveries, and deaths.
Constraint (11) determines the number of recovered individuals (the size of
the compartment R) at the end of the period $t+1$ by aggregating the number of
those recovered from various compartments (of infected and hospitalized
individuals), each with specific recovery rates. Lastly, the constraint (12)
determines the size of the compartment D. It compiles the death toll within
the region, including deaths from the inability to hospitalize individuals in
severe and critical conditions, and deaths among individuals in critical
conditions, adjusting to their respective mortality rates.

Constraints (3)-(12) are for all
$t\in\\{1,\ldots,|\mathcal{T}|-1\\},r\in\mathcal{R},\omega\in\Omega$:

$\displaystyle S_{\omega,t+1,r}$
$\displaystyle=S_{\omega,t,r}-\theta_{\omega,t,r}S_{\omega,t,r}-\rho_{r}(1-h_{\omega,t,r})S_{\omega,t,r}+\sum_{r^{\prime}\in\mathcal{R}}\nu_{r^{\prime}\rightarrow
r}S_{\omega,t,r^{\prime}}-\sum_{r^{\prime}\in\mathcal{R}}\nu_{r\rightarrow
r^{\prime}}S_{\omega,t,r}$ (3) $\displaystyle V_{\omega,t+1,r}$
$\displaystyle=V_{\omega,t,r}+\rho_{r}(1-h_{\omega,t,r})S_{\omega,t,r}-(1-\epsilon)\theta_{\omega,t,r}V_{\omega,t,r}$
(4) $\displaystyle E_{\omega,t+1,r}$
$\displaystyle=E_{\omega,t,r}+\theta_{\omega,t,r}S_{\omega,t,r}-\alpha_{r}E_{\omega,t,r}$
(5) $\displaystyle EV_{\omega,t+1,r}$
$\displaystyle=EV_{\omega,t,r}+(1-\epsilon)\theta_{\omega,t,r}V_{\omega,t,r}-\left(\alpha_{r}(p^{\text{m},v}_{r}+p^{\text{s},v}_{r})+\gamma^{v}_{r}p^{\text{r},v}_{r}\right)EV_{\omega,t,r}$
(6) $\displaystyle I^{\text{m}}_{\omega,t+1,r}$
$\displaystyle=I^{\text{m}}_{\omega,t,r}+\alpha_{r}p^{\text{m}}_{r}E_{\omega,t,r}+\alpha_{r}p^{\text{m},v}_{r}EV_{\omega,t,r}-\gamma^{\text{m}}_{r}I^{\text{m}}_{\omega,t,r}$
(7) $\displaystyle I^{\text{s}}_{\omega,t+1,r}$
$\displaystyle=I^{\text{s}}_{\omega,t,r}+\alpha_{r}p^{\text{s}}_{r}E_{\omega,t,r}+\alpha_{r}p^{\text{s},v}_{r}EV_{\omega,t,r}-\sigma_{r}I^{\text{s}}_{\omega,t,r}$
(8) $\displaystyle H^{\text{s}}_{\omega,t+1,r}$
$\displaystyle=H^{\text{s}}_{\omega,t,r}+\mathcal{A}^{\text{s}}_{\omega,t,r}-\gamma^{\text{s}}_{r}H^{\text{s}}_{\omega,t,r}$
(9) $\displaystyle H^{\text{c}}_{\omega,t+1,r}$
$\displaystyle=H^{\text{c}}_{\omega,t,r}+\mathcal{A}^{\text{c}}_{\omega,t,r}-(1-\varsigma^{\text{c}}_{r})\mu^{\text{c}}_{r}H^{\text{c}}_{\omega,t,r}-\varsigma^{\text{c}}_{r}\gamma^{\text{c}}_{r}H^{\text{c}}_{\omega,t,r}$
(10) $\displaystyle R_{\omega,t+1,r}$
$\displaystyle=R_{\omega,t,r}+\gamma^{v}_{r}p^{\text{r},v}_{r}EV_{\omega,t,r}+\gamma^{\text{m}}_{r}I^{\text{m}}_{\omega,t,r}+\gamma^{\text{s}}_{r}H^{\text{s}}_{\omega,t,r}+\varsigma^{\text{ks}}_{r}\gamma^{\text{ks}}_{r}\mathcal{K}^{\text{s}}_{\omega,t,r}+\varsigma^{c}_{r}\gamma^{\text{c}}_{r}H^{\text{c}}_{\omega,t,r}$
(11) $\displaystyle D_{\omega,t+1,r}$
$\displaystyle=D_{\omega,t,r}+(1-\varsigma^{\text{ks}}_{r})\mu^{\text{ks}}_{r}\mathcal{K}^{\text{s}}_{\omega,t,r}+(1-\varsigma^{\text{c}}_{r})\mu^{\text{c}}_{r}H^{\text{c}}_{\omega,t,r}+\mathcal{K}^{\text{c}}_{\omega,t,r}$
(12)

(G2) _Resource capacity constraints:_ Constraints (13) to (15) determine the
values of the decision variables. The initial number of ventilators in
hospitals is limited. Consider that $\Delta_{j}$ is the number of new
ventilators available for allocation in period (stage) $j\in\mathcal{J}$.
Constraints (13) determine the allocation of these ventilators to hospitals in
different regions. Constraints (14) determine the total number of ventilators
available in region $r$ at period $j$. Since the number of time periods
between every two consecutive decision periods is greater than one, we ensure
(via constraints (15)) that the total number of available ventilators does not
change until the next decision period.

$\displaystyle\sum_{r\in\mathcal{R}}x_{\omega,j,r}\leq\Delta_{j}$
$\displaystyle\forall j\in\mathcal{J},\omega\in\Omega,$ (13)
$\displaystyle\mathcal{X}_{\omega,j,r}=\mathcal{X}_{\omega,j-1,r}+x_{\omega,j,r}$
$\displaystyle\forall j\in\mathcal{J},r\in\mathcal{R},\omega\in\Omega,$ (14)
$\displaystyle\mathcal{X}_{\omega,t,r}=\mathcal{X}_{\omega,t-1,r}$
$\displaystyle\forall
t\in\mathcal{T}\setminus\mathcal{J},r\in\mathcal{R},\omega\in\Omega.$ (15)

Constraint (16) connects the decision variables with the state variables. This
constraint determines the number of critically ill patients admitted to a
hospital, $\mathcal{A}^{\text{c}}_{\omega,t,r}$. This number depends on the
demand for critical care, the number of ventilators available, and the number
of beds available.

Constraint (17) determines the number of severely ill patients admitted to a
hospital, $\mathcal{A}^{\text{s}}_{\omega,t,r}$. This number depends on the
demand for care and beds available after the admission of patients in critical
condition. In Appendix D we detail the linearization of “min" functions in
constraints (16) and (17). The restrictions (18) and (19) determine the number
of patients with severe and critical conditions who were not admitted to a
hospital due to limited resources. Constraints (16) to (19) are for all
$t\in\mathcal{T},r\in\mathcal{R},\omega\in\Omega$.

$\displaystyle\mathcal{A}^{\text{c}}_{\omega,t,r}$
$\displaystyle=\min\Bigl{\\{}\alpha_{r}p^{\text{c}}_{r}E_{\omega,t,r},\>\mathcal{X}_{\omega,t,r}-H^{\text{c}}_{\omega,t,r},b_{r}-H^{\text{c}}_{\omega,t,r}-H^{\text{s}}_{\omega,t,r}\Bigr{\\}}$
(16) $\displaystyle\mathcal{A}^{\text{s}}_{\omega,t,r}$
$\displaystyle=\min\Biggl{\\{}\sigma
I^{\text{s}}_{\omega,t,r},\>b_{r}-H^{\text{c}}_{\omega,t,r}-H^{\text{s}}_{\omega,t,r}-\mathcal{A}^{\text{c}}_{\omega,t,r}\Biggr{\\}}$
(17) $\displaystyle\mathcal{K}^{\text{c}}_{\omega,t,r}$
$\displaystyle=\alpha_{r}p^{\text{c}}_{r}E_{\omega,t,r}-\mathcal{A}^{\text{c}}_{\omega,t,r}$
(18) $\displaystyle\mathcal{K}^{\text{s}}_{\omega,t,r}$ $\displaystyle=\sigma
I^{\text{s}}_{\omega,t,r}-\mathcal{A}^{\text{s}}_{\omega,t,r}$ (19)

(_G3_ :) _Resource allocation constraints:_ Equity in healthcare resource
distribution, especially during emergencies, is vital for public health
policy. Defining and measuring equity in public health decisions is a complex
task with no widely accepted standard. We present several approaches used by
the literature to ensure fair allocation of healthcare resources (see our
discussion of (20)-(22)). Each constraint, when applied individually, offers a
distinct approach to resource allocation. Constraint (20), adapted from (Yin
and Buyuktahtakin, 2021a), uses information about the proportion of the
expected critically ill individuals in a region (the first term in this
constraint) and the corresponding population proportion of the region (the
second term in this constraint) to allocate ventilators. It uses an equity
threshold parameter, denoted by $k\in(0,1]$, to gauge need-based equity in
resource allocation. As $k$ approaches zero, equity in resource allocation is
enforced. As the value of $k$ approaches one, the equity requirements become
less stringent. In Appendix D, you will find our approach to linearizing this
constraint.

$\displaystyle\left\lvert\frac{\sum_{t\in\mathcal{T}}\sum_{\omega\in\Omega}p_{\omega}H^{\text{c}}_{\omega,t,r}}{\sum_{t\in\mathcal{T}}\sum_{r^{\prime}\in\mathcal{R}}\sum_{\omega\in\Omega}p_{\omega}H^{\text{c}}_{\omega,t,r^{\prime}}}-\frac{n_{r}}{\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}}\right\rvert\leq
k$ (20)

The restrictions (21) present a variation in the population-proportional
method for the allocation of resources (Dangerfield et al., 2019; World Health
Organization, 2020). We include an equity threshold parameter $\zeta$. This
constraint ensures that each region receives ventilators in proportion to its
share of the total population, with $\zeta$ providing adjustable scaling to
fine-tune how strictly the allocation of resources follows population size.
For $\zeta=1$, the additional ventilators available for use during period $j$,
denoted as $\Delta_{j}$, are distributed according to the population size of
each region. As $\zeta$ decreases, this constraint relaxes, and the equity
requirements become less stringent.

$\displaystyle\sum_{\omega\in\Omega}p_{\omega}x_{\omega,j,r}\geq\left(\frac{n_{r}}{\sum_{r^{\prime}\in\mathcal{R}}N_{r^{\prime}}}\right)\zeta\Delta_{j},\hskip
7.22743pt\forall j\in\mathcal{J},r\in\mathcal{R}.$ (21)

Both constraints (20) and (21), focus on ensuring an equitable distribution of
critical resources. One of the constraints uses a need-based approach, and the
other population-based approach to resource allocation.

Constraints (22) ensure that ventilators are distributed equally, ignoring
differences in population size or the specific needs of each region. As a
result, the same number of ventilators is allocated to every region.

$\displaystyle\sum_{\omega\in\Omega}p_{\omega}x_{\omega,j,r}\geq\left\lfloor\dfrac{\Delta_{j}}{|\mathcal{R}|}\right\rfloor,\hskip
7.22743pt\forall j\in\mathcal{J},r\in\mathcal{R}.$ (22)

Notice that we solve the compartmental MSP model with either constraints (20),
or (21), or (22) to enforce equity or equality in resource allocation. We also
solve the compartmental MSP model without constraints (20), (21), and (22) to
demonstrate the effects of a utilitarian approach that aims to minimize the
overall cost of resource allocation.

(G4) _Non-anticipativity constraints (NACs):_ These constraints ensure that
decisions made at a given time (decision stage) are based solely on the
information available up to that point. This prevents the use of future, yet
unknown, information in decision making. Specifically, constraints (23) ensure
that decisions that share an identical scenario path up to time $j$ must be
consistent and identical across shared paths. In (23),
$\mathcal{S}_{\omega,j,r}$ represents the set of scenarios that are
indistinguishable from $\omega$ at time $t$. In Appendix E, we provide an
example that demonstrates the implementation of these constraints.

$\displaystyle x_{\omega,j,r}-x_{\omega^{\prime},j,r}$ $\displaystyle=0,$
$\displaystyle\forall\omega\in\Omega,\omega^{\prime}\in\mathcal{S}_{\omega,j,r},j\in\mathcal{J},r\in\mathcal{R}.$
(23)

(G5) _Other constraints:_ Constraints (24) restrict the decision variables to
be non-negative integers. Constraints (25) and (26) restrict the state
variables to be nonnegative. The following constraints are also for all
$\omega\in\Omega,r\in\mathcal{R}$.

$\displaystyle x_{\omega,j,r}\in\mathbb{Z}_{\geq 0}$ $\displaystyle\forall
j\in\mathcal{J},$ (24) $\displaystyle\ell_{\omega,t,r}\geq 0$
$\displaystyle\forall\ell\in\mathcal{L},t\in\mathcal{T},$ (25)
$\displaystyle\mathcal{A}^{\text{s}}_{\omega,t,r},\mathcal{A}^{\text{c}}_{\omega,t,r},\mathcal{K}^{\text{s}}_{\omega,t,r},\mathcal{K}^{\text{c}}_{\omega,t,r}\geq
0$ $\displaystyle\forall t\in\mathcal{T}.$ (26)

##  4 Data Collection and Model Verification and Validation

To validate the compartmental MSP model, we developed a case study using
county-level data from Arkansas during January to May 2021. During this period
of the COVID-19 pandemic, vaccines were available to the public. Next, we
discuss our data collection, data processing, and our approach to model
validation and verification. The SVEIHR and MSP compartmental models presented
in Section 3 are modeled using Python 3.10. The optimization model in Section
3.3 is solved using the DOcplex API from CPLEX version 22.1.1. We performed
our computational analysis on a system equipped with a 32-core CPU and 192 GB
of RAM.

###  4.1 _Data collection and processing_

Data collection: We collect data about county-level vaccination uptake from
January to May 2021 from the CDC (2021), and data about VH from the CTIS
dataset (Salomon et al., 2021). Data on healthcare resources come from the
CovidCareMap (2020) project. Data are organized by facility, county, referral
region, and state. Table 1 summarizes the data for Arkansas.

Table 1: Healthcare capacity by region in Arkansas. Region | # Counties | Population (Ratio) | Licensed Beds | ICU Beds
---|---|---|---|---
R1 | 18 | 394,446 (14.45%) | 1,405 | 77
R2 | 14 | 1,551,512 (56.84%) | 9,357 | 566
R3 | 13 | 171,946 (6.30%) | 490 | 18
R4 | 30 | 611,686 (22.41%) | 2,032 | 105
Total | 75 | 2,729,590 | 13,284 | 766

Figure 4(a) presents the county-level distribution of hospitals and licensed
beds in Arkansas at the beginning of the COVID-19 pandemic. Figure 4(b)
illustrates the relationship between population size and the number of
licensed beds in each county. We cluster counties of Arkansas into four
regions based on their VH. Below, we provide details of our clustering
approach. Figure 4(c) summarizes the total number of licensed beds and
intensive care unit (ICU) beds in each region. We assume that each ICU bed has
a ventilator. The results of these figures emphasize the variations in
healthcare capacity across the state.

Figure 4: (a) County-level distribution of licensed hospital beds in Arkansas.
(b) Relationship between population size and number of licensed hospital beds
by county. (c) The total number of licensed and ICU beds, aggregated by region
of study.

Clustering of the data: In an effort to reduce the size of the decision tree
of the MSP model, we cluster counties of Arkansas into four regions. The
counties are clustered according to their VH similarity. We measure similarity
using the Euclidean distance between county-level VH time series data. We use
the SciPy Python dendrogram tool to visualize the clustering of counties that
belong to the same region (Figure 5(a)) (Virtanen et al., 2020). We use as a
threshold a distance of 0.4 and identify four distinct regions. Figures
5(b)-(e) illustrate VH trends within each region. Each thin line represents
the VH of a county, and the dotted line represents the average VH of that
region. In particular, Region 3, marked red, starts with the highest VH and,
despite the declining trend, remains the highest throughout the study period.
In contrast, Region 2, marked green, begins with a VH of 0.39, which then
drops dramatically, ending as the region with the lowest VH. Figure 5(f) maps
the regions geographically, highlighting that Arkansas’s major cities,
Fayetteville, Fort Smith, and Little Rock, are in Region 2, which has the
lowest VH.

Figure 5: (a) Hierarchical clustering of CTIS VH data for counties in
Arkansas. (b)-(e) Region-specific VH trends over five months (January to May
2021). (f) Spatial distribution of the four regions.

Other data processing: Our model employs the parameter $\Delta_{j}$ to
indicate the number of ventilators available for distribution among different
regions at the start of each decision period, $j$. The COVID-19 pandemic
significantly increased the demand for ventilators, a critical medical device
that saved many lives. In response, under Defense Production Act (DPA),
various manufacturers, including auto manufacturers, quickly began to increase
the production of ventilators or repurpose their facilities for this task
(Albergotti and Siddiqui, 2020). For example, General Motors (GM) was
contracted to supply 30,000 medical ventilators for the national stockpile
(Massey and Cole, 2020). In April 2020, AdvaMed reported that its member
companies were ramping up production from an average of 2,000 to 3,000
ventilators per week to an expected 5,000 to 7,000 per week. This was a
notable increase from the 700 ventilators they produced weekly for the U.S.
market in 2019 AdvaMed (2020). Due to the lack of real-world data related to
the additional ventilators allocated to Arkansas after the enactment of DPA,
we use the following approach to generate these data.

We use equations (27) to (29) to determine $\Delta_{j}$. Let $\underline{p}$
represent the initial stockpile of ventilators in Arkansas prior to enactment
of DPA, and let $\overline{p}$ represent an upper bound on the total number of
available ventilators. Let $p_{j}$ represent the additional ventilators made
available in period $j$.

$\displaystyle\Delta_{1}$ $\displaystyle=\underline{p}$ (27)
$\displaystyle\Delta_{j}$ $\displaystyle=\Delta_{j-1}+p_{j}$
$\displaystyle\forall j\in\mathcal{J}\setminus\\{1\\}$ (28)
$\displaystyle\Delta_{j}$ $\displaystyle\leq\overline{p}$
$\displaystyle\forall j\in\mathcal{J}$ (29)

We begin with an initial stockpile ($\underline{p}$) of 100 ventilators. In
each subsequent time period, we increase the number of ventilators ($p_{j}$)
by 50. As a result, the number of ventilators available for distribution in
Arkansas follows this pattern: $\Delta_{1}=100$, $\Delta_{5}=150$,
$\Delta_{9}=200$, $\Delta_{13}=250$, and $\Delta_{17}=300$.

Please refer to tables in Appendix G for a complete list of the parameter
values used in our MSP compartmental model. The table also presents the tuned
parameters, which are outcomes of the model calibration process detailed in
Section 4.2.

###  4.2 _Calibration of the compartmental model_

The calibration of SVEIHR focuses on fine-tuning the model parameters to
ensure that it accurately represents the dynamics of disease progression in
each region of study. We collected data from COVID-19 Open Data (Wahltinez et
al., 2020). This dataset provides county-level data. For the purpose of this
study, we aggregate county-level data for Arkansas into the four
aforementioned regions.

We use data from two critical time periods in the U.S. The first spans 20
weeks, from October 25, 2020 to March 7, 2021. This period includes the peak
of the wave of coronavirus infections during the Fall of 2020. We use these
data to evaluate SVEIHR’s power in predicting new infections. The second spans
20 weeks, from January 10 to May 23, 2021, and corresponds to the initial
stages of vaccine rollout. We use these data to evaluate SVEIHR’s power in
predicting vaccinations. These study periods do not coincide with the start of
the pandemic. Therefore, we adjust the size of the initial susceptible
population to account for individuals who were infected, recovered, or died
prior to the start of the study period.

Figure 6 summarizes the results of our model calibration of the baseline
scenario. Each plot in this figure presents the actual data and SVEIHR model
predictions of the weekly number of infections and vaccination intakes for
each region of study. New infections at time $t$ are calculated using this
equation:
$\alpha(p^{\text{m}}+p^{\text{s}})E_{t,r}+\alpha(p^{\text{m},v}+p^{\text{s},v})EV_{t,r}$.
The cumulative number of vaccinations up to week $t$ is calculated using
$\sum_{t}(1-h_{t,r})\rho_{r}S_{t,r}$. We employ the Optuna framework to
optimize parameter tuning (Akiba et al., 2019). Table 2 summarizes the
different performance metrics that we use to evaluate the performance of
SVEIHR model of the baseline scenario, including root-mean-square error
(RMSE), mean absolute error (MAE), and mean absolute percentage error (MAPE).
The low value of MAPE highlights the accuracy of SVEIHR in tracking new
infections and vaccination trends in each region.

Figure 6: Comparison of predictions from the calibrated SVEIHR model and actual data across four regions for the baseline scenario, illustrating new infections and cumulative vaccinations over time. Table 2: Statistical analysis of prediction error from the calibrated SVEIHR model for the baseline scenario. | Infections: Fig. 6(a) - (d) | Vaccinations: Fig. 6 (e)-(h)
---|---|---
Region | RMSE | MAE | MAPE | RMSE | MAE | MAPE
R1 | 705.43 | 545.28 | 0.10% | 17,094.55 | 15,883.84 | 2.00%
R2 | 2,288.45 | 1,877.74 | 0.16% | 62,287.10 | 57,541.79 | 1.50%
R3 | 394.94 | 292.13 | 0.10% | 4,531.79 | 4,177.05 | 4.06%
R4 | 765.10 | 620.95 | 0.06% | 26,119.66 | 24,296.23 | 2.11%

###  4.3 _Value of stochastic solution_

The value of stochastic solution (VSS) quantifies the benefit of incorporating
uncertainty into decision-making processes by comparing the results of the
compartmental MSP with those of its deterministic counterpart. To compute the
VSS for each decision stage $j\in\mathcal{J}$, we adopt the methodology
proposed by Escudero et al. (2007). In this method, the expected value (EV)
model, which operates on the expected values of uncertain parameters, is used
as the deterministic benchmark. Let $\mathit{EEV}_{j}$ denote the optimal
value of the proposed compartmental MSP model when the solution up to the
decision period immediately preceding $j$ is fixed to the solution from the
deterministic model. Let $\mathcal{Z}^{*}$ denote the optimal objective
function value of the compartmental MSP model. Therefore, VSS at stage $j$,
represented as $\mathit{VSS}_{j}=\mathit{EEV}_{j}-\mathcal{Z}^{*}$, captures
the value derived from accommodating uncertainty in the model up to decision
stage $j$.

The evolution of the VSS across the decision stages is illustrated in Figure
7(a). For the model used in this illustration, there are 5 decision stages
that correspond to weeks 1, 5, 9, 13, and 17. The results of this figure
highlight the benefits derived from the proposed compartmental MSP model. The
VSS starts at zero, which is expected since
$\mathit{EEV}_{1}=\mathcal{Z}^{*}$. As the model progresses through the
stages, the VSS trends upward, reflecting the benefits of incorporating
uncertainty into strategic decision making. We notice a decrease in the rate
of change in VSS from week 13 to 17. This is mainly because the mean change in
VH is lower during these weeks. The results in Figure 12 show that the mean
change in VH from week 13 to week 17 ranged from 0.03 to 0.05, which is lower
than 0.06-0.09 seen from weeks 9 to 13.

Note that $VSS_{j}=87$ in the illustrative example in Figure 7(a) appears to
be a rather small number compared to the optimal solution
$\mathcal{Z}^{*}=43,900.$ In this example, a critically ill patient uses a
ventilator for three weeks, the recovery rate is 51% and $\overline{p}=300$.
This leads to a total of 1,000 ventilators available for allocation and use
over five months, and a maximum of 729 critically ill patients to be saved
using these ventilators. Thus, a $VSS_{j}=87$ represents 12.21% of this total,
which is pretty high and highlights the value of using the compartmental MSP.

Figure 7: (a) Value of stochastic solution of the proposed model. (b) Expected
ventilator allocation to different regions.

Note. Compartmental MSP model include constraint (21) with $\zeta$ is set to
0.8. The total number of ventilators available in each decision period:
$\Delta_{1}=100,\Delta_{5}=150,\Delta_{9}=200,\Delta_{13}=250,\Delta_{17}=300$.

##  5 Discussion of the Results

RQ1: _Can we develop an effective resource allocation model of critical
resources needed during disease outbreaks?_

To respond to this research question, we develop the compartmental MSP model
that facilitates resource allocation decisions under uncertainty. We make
three observations: ($i$) _The compartmental MSP model facilitates capacity
allocation decisions given the progression and spread of a disease._ One could
also simulate the SVEIHR model for different values of the number of
ventilators available during the planning period and for different assignments
of these ventilators to regions of study. However, such an approach would be
time-consuming. ($ii$) _The MSP model, rather than its deterministic
counterpart, leads to better decision making under uncertainty._ We
demonstrate this via our discussion of the VSS in Section 4.3. ($iii$) _The
MSP model facilitates incorporating fairness in capacity allocation
decisions._

Let us demonstrate these observations using the example in Figure 7(b). The
figure shows the allocation of ventilators among our four study regions over a
20-week planning period. The decision periods are weeks 1, 5, 9, 13 and 17. In
this experiment, we use an equity threshold $\zeta=0.8$, which enforces that
80% of the ventilators are allocated based on the population size of each
region. The remaining 20% of the ventilators are allocated in a way that
minimizes the objective function. To help visualize the results, we use a
different color for each region. The hatched portions of the bars represent
the allocation of the remaining 20% of the ventilators. These results indicate
that in the early stages (weeks 5 and 9), the model allocates the additional
ventilators to region 2 since this region’s population size and density are
the highest. Its population corresponds to 56.84% of the total. Allocating
additional ventilators to this region has the greatest impact on reducing the
number of deaths. Later in the planning period, additional ventilators are
allocated to regions 1 and 4. This is due to an increase in COVID-19 cases and
increased vulnerability within these regions. These dynamic adjustments of
resource allocation emphasize the model’s ability to respond to changing needs
throughout the pandemic, providing an efficient and fair distribution of
resources.

These observations highlight the value that the MSP model brings to decision
makers during health crises, enabling them to balance between efficiently and
equitably allocating scarce critical resources.

RQ2: _What should be effective resource allocation strategies at the start of
an outbreak?_

To answer this research question, we designed the following experiment. We
create three sets of problems. In each set, the time at which additional
ventilators are deployed is different. In set 1, additional ventilators are
deployed beginning in period 1. In set 2, that happens in period 5, and in set
3, it happens in period 9 of a 20-week planning period. Each set consists of 9
problems, one for different combinations of the initial stockpile
($\underline{p}$) and the number by which we increase the number of
ventilators deployed in each period ($p_{j}$). These experiments help us
understand the impact of the timing and scope of decisions related to
increasing the availability of critical resources (such as the timing and
scope of the DPA) on the total expected number of deaths. A summary of the
results from these experiments is presented in Figure 8. Detailed solution
results for a problem from set 1 with $\underline{p}=100,p_{j}=50$, are
presented in Figure 7(b).

We make the following observations: ($i$) _The timing of decisions related to
the deployment of critical resources greatly impacts the expected number of
deaths._ Based on our experiments, beginning the deployment of ventilators in
the 1${}^{\text{st}}$ rather than in the $9^{\text{th}}$ week of the planning
period could lead to an 1.30% average of decrease in the expected number of
deaths, which corresponds to 571 individuals. ($ii$) _The size of the initial
stockpile and the number by which we increase the availability of our critical
resource every period greatly impact the expected number of deaths._ We
develop a simple linear regression model to explore the relationships between
the time of deployment of additional ventilators (${s}$), initial stockpile
($\underline{p}$), the number of additional ventilators deployed ($p_{j}$),
and expected deaths ($\overline{D}$). Since in our experiments, we keep
$p_{j}$ constant over time, below, we use $p$ instead. The model is described
by:

$\overline{D}=44527.66+285.41\;s-5.44\;\underline{p}-4.81\;p$ (30)

Figure 8: Impact of varying timing of the deployment and the number of new
ventilators on the expected death toll.

Note. ∗ The ventilator allocation is shown in Fig. 7.

This regression model has an adjusted $R^{2}=0.9477$, which shows that the
timing of ventilator deployment, initial stockpile, and the amount by which we
increase resource availability explain approximately 94.77% of the variation
in the expected number of deaths. Based on our experiments, delaying the start
of the deployment of additional ventilators by one month could lead to an
average increase in the expected number of deaths by 285.41/month,
highlighting the importance of prompt action. Each additional ventilator in
the initial stockpile and each additional ventilator deployed lead to a
decrease in the expected number of deaths by 5.44 (equivalent to 1.09/month)
and 4.81 (0.962/month), respectively, highlighting the importance of
maintaining a large stockpile and scalable production responses.

RQ3: _What are the trade-offs between equity, equality, and effectiveness of
resource allocation at the onset of a pandemic?_

To answer this research question, we designed the following experiment. We
create four problems, P1 to P4. The compartmental MSP model without
constraints (20), (21) and (22) determines a resource allocation strategy that
minimizes the expected number of deaths. This is our first problem (P1), which
maximizes efficiency. The compartmental MSP model with constraints (20)
focuses on the equitable distribution of resources. As we change the value of
the tolerance level $k\in(0,1]$ from 0 to 1, the equity requirements become
less stringent. This is our second problem (P2). The compartmental MSP model
with constraints (21) is our third problem (P3), which also focuses on the
equitable distribution of resources. As we change the value of the tolerance
level $\zeta\in(0,1)$ from 0 to 1, the equity requirements become rigid. The
compartmental MSP model with constraints (22) is our fourth problem (P4),
which focuses on the equal distribution of resources. The decision periods and
the availability of resources are similar to the experiment described in RQ1.

We make the following observations: ($i$) _Strategies that focus on equality
underperform those that focus on efficiency and equity._ Figure 9(e) presents
the total expected number of deaths as the value of $k$ changes from $10^{-3}$
to 1 (P2), and $\zeta$ changes from 0.97 to 0 (P3). Straight lines at 44,594
and 43,900 represent the expected deaths when ventilators are distributed
equally (P4) and efficiently (P1), respectively. Based on these results, equal
distribution of resources leads to the highest expected number of deaths.

In our experiments, as the value of $k$ increases, the equity restrictions
become less stringent, and the expected number of deaths decreases from 44,453
to 43,900, representing a 1.24% decrease. Decreasing $\zeta$ makes equity
restrictions less stringent, reducing the expected number of deaths from
44,208 to 43,900, which represents a 0.70% decrease. The expected number of
deaths is lowest when $k=1$, or $\zeta=0$, which represents a model that
focuses on efficiency (P1).

($ii$) _The price of ensuring equitable allocation of resources changes over
time and space. This price is highest during the peak of a disease outbreak
and in highly populated regions._ Figures 9(a) - (d) present the expected
number of deaths in each region as the value of $k$ increases from $10^{-3}$
to 1, and the value of $\zeta$ decreases from 0.97 to 0. As $k$ increases and
as $\zeta$ decreases, the equity restrictions become less stringent. These
changes affect the regions of study differently. Region 3, which has the
smallest population size and the lowest density, benefits the most when equity
is enforced (the value of $k$ is close to $10^{-3}$). Region 2, which has the
largest population size and the highest density, benefits the most when equity
is not enforced (the value of $k$ is close to 1). Table 3 summarizes the
expected number of deaths for each problem. The table also summarizes the
change in the number of deaths for all problems compared to the results of P1.

Figure 9: (a)-(d) Expected deaths at different equity levels for each region. (e) Total expected deaths for different resource allocation strategies. Table 3: Total expected deaths when applying different allocation strategies. Strategy | Expected Death
---|---
R1 | R2 | R3 | R4 | Total | Diff from (P1)
(P1) Max efficiency | 8,427 | 15,524 | 3,585 | 16,364 | 43,900 | _
(P2) Equity ($k$) | 0.001 | 8,477 | 15,988 | 3,569 | 16,419 | 44,453 | 553
0.05 | 8,475 | 15,855 | 3,567 | 16,420 | 44,317 | 417
1 | 8,425 | 15,527 | 3,585 | 16,364 | 43,900 | 0
(P3) Equity ($\zeta$) | 0.97 | 8,375 | 16,004 | 3,556 | 16,272 | 44,208 | 308
0.8 | 8,384 | 15,912 | 3,561 | 16,286 | 44,143 | 243
0.6 | 8,395 | 15,808 | 3,567 | 16,303 | 44,073 | 173
(P4) Equal allocation | 8,349 | 16,480 | 3,505 | 16,260 | 44,594 | 694

Figure 10: Evaluating the timing of resource allocation decisions on the
expected number of deaths.

Note. Summary of results from solving P2 at different time $t^{\prime}$.

To evaluate the impact that the timing of resource allocations has on the
expected number of deaths, we performed the following experiment. We solve P2
several times by changing the time period over which we enforce the equity
constraints (20). Thus, instead of adding for ${t\in\mathcal{T}}$, we add over
${t\in\\{\mathcal{T}|t\geq t^{\prime}\\}}$, which means that the equity in
resource allocation is enforced after $t^{\prime}$. Prior to $t^{\prime}$,
resource allocation is focused on maximizing efficiency. We solve P2 for
different values of $t^{\prime}$ and for different values of $k.$ The results
of these experiments are summarized in Figure 10. These results show a
decreasing trend in the expected number of deaths as $k$ approaches 1 and
$t^{\prime}$ approaches $|\mathcal{T}|.$

This analysis is important for decision makers considering the appropriate
level of equity tolerance that harmonizes both public health objectives and
ethical standards. Our analysis provides a data-driven framework for such
decisions, emphasizing the necessity of a comprehensive approach that aims to
balance the imperative of saving lives and the principles of fairness and
justice in resource allocation.

##  6 Conclusion and future research directions

Summary of the Proposed Research: This study proposes a modeling framework
that integrates a compartmental model within a multi-stage stochastic program
(MSP) to tackle the complex challenge of allocating scarce healthcare
resources during infectious disease outbreaks. The compartmental model
simulates the dynamic changes of vaccine hesitancy (VH) over time and space,
improving the ability of MSP to adapt resource allocation strategies to the
dynamics of disease spread and changes in population behavior toward vaccines.
This compartmental MSP model also allowed us to analyze the price of fairness
in resource allocation. Through our numerical analysis, we gain a deeper
understanding of the impacts that the deployment time and size of the resource
allocation, the effectiveness in managing the spread of the disease, and the
price of fairness in resource allocation have on the expected number of
deaths. The results of our analysis highlight the importance of advanced
planning and preparedness for future outbreaks, effective management of
critical resources, and the balance of equity and effectiveness in resource
allocation.

Research Contributions: The main contribution of this work is the development
of a compartmental MSP model that captures the impact of VH in the
effectiveness of allocating critical healthcare resources during a disease
outbreak. The results of our numerical analysis show how public behaviors
towards the COVID-19 vaccine impacted disease dynamics and the need for
healthcare resources.

Our review of the literature showed that there is very limited research on
compartmental MSP models. Our proposed data-driven MSP compartmental model
fills this gap in the literature and provides a robust framework for planning
responses to disease outbreaks. The model supports dynamic decision-making
over time and space, taking into account the unpredictable nature of disease
spread and resource availability. The model shows how dynamic, evidence-based
strategies can enhance public health responses. The model addresses the price
of fairness in the allocation of scarce healthcare resources. The model
supports effective vaccination strategies and resource allocations that ensure
greater immunity and mitigate the impact of outbreaks.

Research Findings: Key observations from our study include: ($i$) The proposed
compartmental MSP model offers a more effective framework for resource
allocation than an MSP model, as it captures the impact of disease progression
(accounting for varying levels of VH) on the need for critical healthcare
resources. ($ii$) The proposed model supports resource allocation in different
regions over time by considering both population size and immediate needs for
medical resources, and balancing fairness and efficiency in decision making
processes. ($iii$) The size of the initial stockpile, the timing of deployment
of additional resources, and the number of additional resources have a
significant impact on the expected number of deaths. Based on our numerical
analysis, each unit increase in the initial stockpile and each unit increase
of available ventilators could lead to a decrease in the expected number of
deaths by 1.09/month and 0.962/month, respectively. Based on our numerical
analysis, delaying the deployment of critical resources by a month could
increase the expected number of deaths by 285.41/month. ($iv$) An equal
resource allocation strategy results in the highest number of expected deaths
(that is, 44,594 over 5 months, assuming that 1,000 become available during
this period), while a utilitarian resource allocation results in the smallest
(i.e., 43,900 over 5 months). ($v$) The timing of the implementation of equity
in resource allocation greatly impacts health outcomes. Delaying its
implementation reduces deaths to levels comparable to those achieved via a
utilitarian approach.

Future Research Directions: The model presented in this study can be extended
in multiple ways to address other relevant problems related to the allocation
of critical healthcare resources. The compartmental MSP models could be
extended to capture the interactions between various health interventions
(i.e., masking mandates, travel restrictions, school closures, etc.) and the
allocation of critical resources under uncertainty. These interventions affect
the distribution of VH and the infection rate over time and space. Assume that
the VH and infection rate are the uncertain problem parameters considered. In
this case, the distribution of random variables is affected “endogenously" by
decisions made by the model. MSPs with endogenous uncertainty are more
difficult to solve compared to models with exogenous uncertainty.

Our proposed model can be extended to support other comparative studies across
different regions or demographics to understand the effectiveness of
strategies for allocating critical healthcare resources. In addition to the
resources considered in our study, the proposed model could be extended to
evaluate the impact of allocation strategies for other critical health
resources (such as personal protective equipment, testing kits, and hospital
beds) on the spread of the disease and the expected number of deaths.

These model extensions provide opportunities to develop novel solution
approaches and address relevant healthcare-related problems.

##  Appendix A Impact of Healthcare Resources on Disease Outcomes

The availability of healthcare resources such as beds and ventilators plays a
critical role in managing the outcomes of infectious diseases. We performed
two simulations of the SVEIHR model to explore these relationships. The first
simulation assumes unlimited healthcare resources, and the second introduces
bed and ventilator availability constraints to reflect a real-world scenario
with limited resource availability. The outcomes of the first simulation,
depicted in Figure 11(a), show a low number of hospitalizations and deaths due
to the assumption of infinite resources. However, the results of the second
simulation, presented in Figure 11(b), reveal a higher number of deaths due to
limited healthcare resources. The increase in the cumulative number of deaths
becomes apparent when the number of individuals in critical condition exceeds
the number of beds available. The contrast is notable. When the number of beds
is unlimited, the expected number of deaths is 5,239, whereas when the number
of beds is limited, the model results in 28,507 deaths.

Figure 11: Comparative impact of healthcare resource availability on the
SVEIHR model outcomes, illustrating differences in hospitalization and death
rates under unlimited and limited resources scenarios.

##  Appendix B Modeling Uncertainty in Vaccine Hesitancy

Figure 12 shows that the rate of change in vaccine hesitancy fits the Normal
Distribution. The corresponding mean and standard deviation changes over time.

Figure 12: (a)-(d) Distribution of the monthly rate of change of vaccine
hesitancy for different clusters in the study, with each row representing
transitions from the previous month.

##  Appendix C Model Notation

Tables 4, 5, and 6 represent the model notations that are used throughout the
paper.

Table 4: Sets and indices in the proposed MSP. Set | Index | Description
---|---|---
$\mathcal{T}$ | $t$ | Set of time periods.
$\mathcal{J}$ | $j$ | Set of decision time periods, $\mathcal{J}\subset\mathcal{T}$.
$\mathcal{R}$ | $r$ | Set of regions.
$\Omega$ | $\omega$ | Set of scenarios.
$\mathcal{L}$ | $\ell$ | Set of state variables after fitting the SVEIHR model with actual data.
$\mathcal{S}_{\omega,j,r}$ |  | Set of scenarios that are indistinguishable with $\omega$ at decision stage $j$.
Table 5: Variables in the proposed MSP.

Variable |  Description (for $r\in\mathcal{R}$, at time $t\in\mathcal{T}$, $j\in\mathcal{J},$ and under scenario $\omega\in\Omega$)
---|---
${x_{\omega,j,r}}^{(1)}$ | Ventilators allocated at the start of decision period $j$.
$\mathcal{X}_{\omega,t,r}$ | Cumulative numbers of ventilators at the end of time $t$.
$S_{\omega,t,r},\;V_{\omega,t,r}$ | Number of susceptible individuals who are unvaccinated and vaccinated, respectively.
$E_{\omega,t,r},\;EV_{\omega,t,r}$ | Number of exposed individuals who are unvaccinated and vaccinated, respectively.
$I^{\text{m}}_{\omega,t,r},\;I^{\text{s}}_{\omega,t,r}$ | Number of infectious individuals with mild and severe symptoms, respectively.
$H^{\text{s}}_{\omega,t,r},\;H^{\text{c}}_{\omega,t,r}$ | Number of hospitalized individuals with severe and critical conditions, respectively.
$\mathcal{A}^{\text{s}}_{\omega,t,r},\;\mathcal{A}^{\text{c}}_{\omega,t,r}$ | Number of hospital admissions for severe and critical cases, respectively.
$\mathcal{K}^{\text{s}}_{\omega,t,r},\;\mathcal{K}^{\text{c}}_{\omega,t,r}$ | Number of severe and critical cases not admitted to the hospital due to capacity constraints.
$R_{\omega,t,r}$ | Cumulative number of recovered individuals up to period $t$.
$D_{\omega,t,r}$ | Cumulative number of deaths up to period $t$.

  * (1)

decision variable

Table 6: Parameters in the proposed MSP. Parameter | Description
---|---
${h_{\omega,t,r}}$ | Vaccine hesitancy rate at region $r$ and time $t$ under scenario $\omega$.
$\beta_{r}$ | COVID-19 contact rate in region $r$. $\theta_{\omega,t,r}=\beta_{r}S_{\omega,t,r}\left(I^{\text{m}}_{\omega,t,r}+I^{\text{s}}_{\omega,t,r}\right)/n_{r}$.
$\rho_{r}$ | Maximum vaccination rate in region $r$.
$\epsilon$ | Vaccine efficacy.
$\alpha$ | Incubation rate.
$\gamma^{\text{m}}_{r},\;\gamma^{\text{s}}_{r}$ | Recovery rates for mild and severe cases in region $r$.
$\gamma^{\text{ks}}_{r},\;\gamma^{\text{c}}_{r}$ | Recovery rates for non-hospitalized severe and critical cases in region $r$.
$\sigma_{r}$ | Rate from symptoms to hospitalization in region $r$.
$\mu^{\text{ks}}_{r},\;\mu^{\text{c}}_{r}$ | Mortality rates for non-hospitalized severe and hospitalized critical cases in region $r$.
$\varsigma^{\text{ks}}_{r},\;\varsigma^{\text{c}}_{r}$ | Recovery rates for non-hospitalized severe and hospitalized critical cases in the region $r$.
$k$ | Equity threshold in constraint (20).
$\zeta$ | Equity threshold in constraint (21).
$p^{\text{m}},\;p^{\text{m},v},\;p^{\text{s}},\;p^{\text{s},v}$ | Transition percentages to mild and severe cases for unvaccinated and vaccinated individuals.
$p^{\text{r},v}$ | Recovery rate for vaccinated individuals.
$\boldsymbol{\pi}_{l,r}$ | Initial count of individuals in state $\ell$ (compartment in SVEIHR model) at time $t=0$ in region $r$.
$n_{r}$ | Population size of region $r$.
$b_{r}$ | Hospital bed count in region $r$.
$\Delta_{t}$ | Additional ventilators available to allocate at start of time $t$.

##  Appendix D Constraint Linearization

The proposed compartmental MSP model presented in Section 3.3 contains three
non-linear constraints. They can be linearized as follows:

We linearize the ‘min’ functions in constraints (16) and (17) using the big-M
method. To illustrate, let’s consider the linearization of constraint (17).
For ease of demonstration, let us declare a few auxiliary variables:
$a_{1}=\sigma I^{\text{s}}_{\omega,t,r}$ representing the number of severe
case admissions and
$a_{2}=b_{r}-H^{\text{c}}_{\omega,t,r}-H^{\text{s}}_{\omega,t,r}-\mathcal{A}^{\text{c}}_{\omega,t,r}$
representing the remaining bed capacity after accounting for current
hospitalizations.

The original constraint can be written as
$\mathcal{A}^{\text{s}}_{\omega,t,r}=\min\\{a_{1},a_{2}\\}$, indicating that
the number of severe cases admitted is the lesser of the demand or the
remaining capacity. We use a binary variable ($y_{\omega,t,r}$) and big-M
constant ($M$) to enforce $A\geq\min\\{a_{1},a_{2}\\}$ . The original ‘min’
constraint is then replaced by the following system of linear inequalities
(17a)-(17d):

$\displaystyle\mathcal{A}^{\text{s}}_{\omega,t,r}$ $\displaystyle\leq a_{1}$
(17a) $\displaystyle\mathcal{A}^{\text{s}}_{\omega,t,r}$ $\displaystyle\leq
a_{2}$ (17b) $\displaystyle\mathcal{A}^{\text{s}}_{\omega,t,r}$
$\displaystyle\geq a_{1}-M(1-y_{\omega,t,r})$ (17c)
$\displaystyle\mathcal{A}^{\text{s}}_{\omega,t,r}$ $\displaystyle\geq
a_{2}-My_{\omega,t,r}$ (17d)

Next, we linearize the absolute function in constraint (20) using the
following procedure. Again, for demonstration purposes, let us introduce the
following auxiliary variables:
$a_{1}=\sum_{t\in\mathcal{T}}\sum_{\omega\in\Omega}p_{\omega}H^{\text{c}}_{\omega,t,r}$,
and
$a_{2}=\sum_{t\in\mathcal{T}}\sum_{r^{\prime}\in\mathcal{R}}\sum_{\omega\in\Omega}p_{\omega}H^{\text{c}}_{\omega,t,r^{\prime}}$.
Since $a_{2}>0,\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}>0,$ the original
constraint can be written as
$\left\lvert\frac{a_{1}}{a_{2}}-\frac{n_{r}}{\sum_{r^{\prime}\in\mathcal{R}}N_{r^{\prime}}}\right\rvert\leq
k$ and is equivalent to:

$\displaystyle\begin{dcases*}\frac{a_{1}}{a_{2}}-\frac{n_{r}}{\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}}\leq
k\\\
\frac{a_{1}}{a_{2}}-\frac{n_{r}}{\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}}\geq-k\end{dcases*}$
$\displaystyle\Leftrightarrow$
$\displaystyle\begin{cases}a_{1}\left(\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}\right)-a_{2}n_{r}\leq
ka_{2}\left(\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}\right)\\\
a_{1}\left(\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}\right)-a_{2}n_{r}\geq-
ka_{2}\left(\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}\right)\end{cases}$

Constraint (20) can be replaced by the following linear inequalities
(20a)-(20b):

$\displaystyle
a_{1}\left(\sum_{r^{\prime}\in\mathcal{R}}N_{r^{\prime}}\right)-a_{2}n_{r}-ka_{2}\left(\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}\right)$
$\displaystyle\leq 0$ (20a) $\displaystyle
a_{1}\left(\sum_{r^{\prime}\in\mathcal{R}}N_{r^{\prime}}\right)-a_{2}n_{r}+ka_{2}\left(\sum_{r^{\prime}\in\mathcal{R}}n_{r^{\prime}}\right)$
$\displaystyle\geq 0$ (20b)

##  Appendix E Illustration of non-anticipativity constraints (NACs)

Figure 13(a) shows a scenario tree with three stages and nine possible
scenarios. Figures 13(a) and (b) help the reader understand the implications
of the nonanticipativity constraints (NAC). non-anticipativity constraints
ensure that decisions made at a certain stage are the same across all
scenarios that share the same history up to that stage. For example, the root
node (stage 1) in Figure 13(a) is colored cyan. This node is common in all
scenarios. As such, decisions made in stage 1 affect decisions made along all
scenario paths. This is why in Figure 13(b), the nodes in the top layer are
all colored cyan. Therefore, the set of scenarios that are indistinguishable
from $\omega=1$ at time $t=1$ is $\mathcal{S}_{1,1,r}=\\{2,3,\cdots,9\\}$. We
use the following eight NACs to ensure this:
$x_{1,1,r}=x_{2,1,r};\;x_{1,1,r}=x_{3,1,r};\dots;\;x_{1,1,r}=x_{9,1,r}$.

Figure 13: (a) Illustration of a scenario tree with three stages and nine
scenarios, (b) an alternative form to visualize the NACs, which are shown in
red.

##  Appendix F Population Flow

We utilized data related to county-to-county migration flows from 2016 to
2020, sourced from (U.S. Census Bureau, 2020). We estimate the annual movement
of individuals between counties in Arkansas. The results are shown in Figure
14(a). The black dots in this figure represent county centroids and the shaded
areas represent the four different regions of study. Major cities like
Fayetteville, Fort Smith, and Little Rock are labeled, signifying possible
hubs of migration. The simplified flow between regions is depicted in Fig.
14(b). These data are integrated into the compartmental MSP model to update
the susceptible individuals at each time step $t$ and region $r$ in constraint
(3). The weekly movement rate are shown in Table 7. The expression
$\sum_{r^{\prime}\in\mathcal{R}}\nu_{r^{\prime}\rightarrow
r}S_{\omega,t,r^{\prime}}-\sum_{r^{\prime}\in\mathcal{R}}\nu_{r\rightarrow
r^{\prime}}S_{\omega,t,r}$ represents the net movement of susceptible
individuals between regions. Since most migration occurs within regions, the
flow between regions has a negligible impact on the model.

Figure 14: County and regional migration patterns in Arkansas: Analyzing annual inter-county flows and aggregate regional migration trends. Table 7: Weekly population movement rate between studied regions. From \To | R1 | R2 | R3 | R4
---|---|---|---|---
R1 | - | 5.34e-04 | 4.50e-05 | 2.48e-04
R2 | 1.28e-04 | - | 1.92e-05 | 1.45e-04
R3 | 1.58e-04 | 3.34e-04 | - | 3.41e-04
R4 | 1.45e-04 | 4.19e-04 | 7.76e-05 | -

##  Appendix G Model Data

Table 8 presents the parameter values used in the proposed compartmental MSP
model. Table 9 details the parameter values for the four different regions
included in the study. The upper portion of the table includes the tuned
parameters derived from the calibration process, which is described in Section
4.2.

Table 8: Parameter values used in the proposed model. Parameter | Value | Source | Parameter | Value | Source
---|---|---|---|---|---
$|\mathcal{T}|$ | 20 |  | $\gamma^{\text{s}},\gamma^{\text{ks}},\gamma^{\text{c}}$ | 0.637, 0.233, 0.333 | (Zhou et al., 2020b)
$|\mathcal{R}|$ | 4 |  | $\mu^{\text{ks}}$ | 1.000 |
$|\Omega|$ | 81 |  | $p^{\text{m}},p^{\text{s}},p^{\text{c}}$ | 0.800, 0.150, 0.050 | (Wu and McGoogan, 2020)
$\mathcal{J}$ | $\\{1,5,9,13,17\\}$ |  | $p^{\text{r},v},p^{\text{m},v},p^{\text{s},v}$ | 0.850, 0.100, 0.050 | (Barnes and Others, 2023)
$\epsilon$ | 0.950 | (Katella, 2023) | $\varsigma^{\text{ks}}$ | 0.400 |
$\gamma^{v}$ | 1.000 |  | $\varsigma^{\text{c}}$ | 0.510 | (Wu and McGoogan, 2020)
Table 9: Parameter values for each region used in the proposed model.

Parameter | Region ($\boldsymbol{r}$)
---|---
R1 | R2 | R3 | R4
$\beta_{r}$ (1) | 1.873 | 2.162 | 1.805 | 1.781
$\rho_{r}$ (1) | 0.051 | 0.051 | 0.043 | 0.050
$\gamma^{\text{m}}_{r}$ (1) | 0.704 | 0.704 | 0.704 | 0.699
$\sigma_{r}$ (1) | 0.943 | 0.877 | 0.952 | 0.901
$S_{0,r}$ (1) | 299,000 | 1,107,000 | 129,000 | 502,000
$E_{0,r}$ (1) | 3,100 | 3,200 | 100 | 3,500
$I^{\text{m}}_{0,r}$ (1) | 200 | 5,900 | 1,300 | 3,600
$I^{\text{s}}_{0,r}$ (1) | 2,500 | 1,470 | 80 | 60
$b_{r}$ | 1,405 | 9,357 | 490 | 2,032
$\mathcal{X}_{0,r}$ | 77 | 566 | 18 | 105
${n_{r}}/{\sum n_{r}}$ | 0.145 | 0.568 | 0.063 | 0.224

  * (1)

Tuned parameters using Optuna framework, an open-source optimization library
that automates hyperparameter tuning.

## References

  * Sy et al. [2021] Karla Therese L. Sy, Laura F. White, and Brooke E. Nichols.  Population density and basic reproductive number of COVID-19 across united states counties.  _PLOS ONE_ , 16:e0249271, 4 2021.  ISSN 1932-6203.  doi:10.1371/journal.pone.0249271.
  * Moghadas et al. [2021] Seyed M Moghadas, Thomas N Vilches, Kevin Zhang, Chad R Wells, Affan Shoukat, Burton H Singer, Lauren Ancel Meyers, Kathleen M Neuzil, Joanne M Langley, Meagan C Fitzpatrick, and Alison P Galvani.  The impact of vaccination on Coronavirus disease 2019 (COVID-19) outbreaks in the United States.  _Clinical Infectious Diseases_ , 73:2257–2264, 12 2021.  ISSN 1058-4838.  doi:10.1093/cid/ciab079.
  * Brooks and Butler [2021] John T. Brooks and Jay C. Butler.  Effectiveness of mask wearing to control community spread of SARS-CoV-2.  _JAMA_ , 325:998, 3 2021.  ISSN 0098-7484.  doi:10.1001/jama.2021.1505.
  * Chu et al. [2020] Derek K Chu, Elie A Akl, Stephanie Duda, Karla Solo, Sally Yaacoub, Holger J Schünemann, Derek K Chu, Elie A Akl, Amena El-harakeh, Antonio Bognanni, Tamara Lotfi, Mark Loeb, Anisa Hajizadeh, Anna Bak, Ariel Izcovich, Carlos A Cuello-Garcia, Chen Chen, David J Harris, Ewa Borowiack, Fatimah Chamseddine, Finn Schünemann, Gian Paolo Morgano, Giovanna E U Muti Schünemann, Guang Chen, Hong Zhao, Ignacio Neumann, Jeffrey Chan, Joanne Khabsa, Layal Hneiny, Leila Harrison, Maureen Smith, Nesrine Rizk, Paolo Giorgi Rossi, Pierre AbiHanna, Rayane El-khoury, Rosa Stalteri, Tejan Baldeh, Thomas Piggott, Yuan Zhang, Zahra Saad, Assem Khamis, Marge Reinap, Stephanie Duda, Karla Solo, Sally Yaacoub, and Holger J Schünemann.  Physical distancing, face masks, and eye protection to prevent person-to-person transmission of SARS-CoV-2 and COVID-19: a systematic review and meta-analysis.  _The Lancet_ , 395:1973–1987, 6 2020.  ISSN 01406736.  doi:10.1016/S0140-6736(20)31142-9.
  * Blasioli et al. [2023] Emanuele Blasioli, Bahareh Mansouri, Srinivas Subramanya Tamvada, and Elkafi Hassini.  Vaccine allocation and distribution: A review with a focus on quantitative methodologies and application to equity, hesitancy, and COVID-19 pandemic.  _Operations Research Forum_ , 4:27, 3 2023.  ISSN 2662-2556.  doi:10.1007/s43069-023-00194-8.
  * Kermack et al. [1927] William Ogilvy Kermack, A. G. McKendrick, and Gilbert Thomas Walker.  A contribution to the mathematical theory of epidemics.  _Proceedings of the Royal Society of London. Series A, Containing Papers of a Mathematical and Physical Character_ , 115(772):700–721, 1927.  doi:10.1098/rspa.1927.0118.
  * Kong et al. [2022] Lingcai Kong, Mengwei Duan, Jin Shi, Jie Hong, Zhaorui Chang, and Zhijie Zhang.  Compartmental structures used in modeling COVID-19: a scoping review.  _Infectious Diseases of Poverty_ , 11:72, 6 2022.  ISSN 2049-9957.  doi:10.1186/s40249-022-01001-y.
  * Xiang et al. [2021] Yue Xiang, Yonghong Jia, Linlin Chen, Lei Guo, Bizhen Shu, and Enshen Long.  COVID-19 epidemic prediction and the impact of public health interventions: A review of COVID-19 epidemic models.  _Infectious Disease Modelling_ , 6:324–342, 2021.  ISSN 24680427.  doi:10.1016/j.idm.2021.01.001.
  * Buonomo et al. [2022] Bruno Buonomo, Rossella Della Marca, Alberto d’Onofrio, and Maria Groppi.  A behavioural modelling approach to assess the impact of COVID-19 vaccine hesitancy.  _Journal of Theoretical Biology_ , 534:110973, 2 2022.  ISSN 00225193.  doi:10.1016/j.jtbi.2021.110973.
  * Mesa et al. [2023] Daniela Olivera Mesa, Peter Winskill, Azra C Ghani, and Katharina Hauck.  The societal cost of vaccine refusal: A modelling study using measles vaccination as a case study.  _Vaccine_ , 41:4129–4137, 6 2023.  ISSN 0264410X.  doi:10.1016/j.vaccine.2023.05.039.
  * Mehrotra et al. [2020] Sanjay Mehrotra, Hamed Rahimian, Masoud Barah, Fengqiao Luo, and Karolina Schantz.  A model of supply-chain decisions for resource sharing with an application to ventilator allocation to combat COVID-19.  _Naval Research Logistics (NRL)_ , 67:303–320, 8 2020.  ISSN 0894-069X.  doi:10.1002/nav.21905.
  * Bhavani et al. [2021] Sivasubramanium V Bhavani, Yuan Luo, William D Miller, Lazaro N Sanchez-Pinto, Xuan Han, Chengsheng Mao, Burhaneddin Sandıkçı, Monica E Peek, Craig M Coopersmith, Kelly N Michelson, and William F Parker.  Simulation of ventilator allocation in critically ill patients with COVID-19.  _American journal of respiratory and critical care medicine_ , 204:1224–1227, 11 2021.  ISSN 1535-4970.  doi:10.1164/rccm.202106-1453LE.
  * Koonin et al. [2020] Lisa M. Koonin, Satish Pillai, Emily B. Kahn, Danielle Moulia, and Anita Patel.  Strategies to inform allocation of stockpiled ventilators to healthcare facilities during a pandemic.  _Health Security_ , 18:69–74, 4 2020.  ISSN 2326-5094.  doi:10.1089/hs.2020.0028.
  * Bertsimas et al. [2021] Dimitris Bertsimas, Leonard Boussioux, Ryan Cory-Wright, Arthur Delarue, Vassilis Digalakis, Alexandre Jacquillat, Driss Lahlou Kitane, Galit Lukin, Michael Li, Luca Mingardi, Omid Nohadani, Agni Orfanoudaki, Theodore Papalexopoulos, Ivan Paskov, Jean Pauphilet, Omar Skali Lami, Bartolomeo Stellato, Hamza Tazi Bouardi, Kimberly Villalobos Carballo, Holly Wiberg, and Cynthia Zeng.  From predictions to prescriptions: A data-driven response to COVID-19.  _Health Care Management Science_ , 24:253–272, 6 2021.  ISSN 1386-9620.  doi:10.1007/s10729-020-09542-0.
  * Birge and Louveaux [2011] John R Birge and Francois Louveaux.  _Introduction to stochastic programming_.  Springer Science & Business Media, 2011.
  * Yin and Buyuktahtakin [2021a] Xuecheng Yin and I. Esra Buyuktahtakin.  A multi-stage stochastic programming approach to epidemic resource allocation with equity considerations.  _Health Care Management Science_ , 24:597–622, 9 2021a.  ISSN 1386-9620.  doi:10.1007/s10729-021-09559-z.
  * Yin et al. [2023] Xuecheng Yin, I. Esra Buyuktahtakin, and Bhumi P. Patel.  COVID-19: Data-driven optimal allocation of ventilator supply under uncertainty and risk.  _European Journal of Operational Research_ , 304:255–275, 1 2023.  ISSN 03772217.  doi:10.1016/j.ejor.2021.11.052.
  * Albergotti and Siddiqui [2020] Reed Albergotti and Faiz Siddiqui.  Ford and GM are racing to build coronavirus ventilators, but their efforts may be too late, 04 2020.  URL https://www.washingtonpost.com/business/2020/04/04/ventilators-coronavirus-ford-gm/.  (Accessed on 2023-10-23).
  * Mishra et al. [2021] Vaibhav Mishra, Golnoush Seyedzenouzi, Ahmad Almohtadi, Tasnim Chowdhury, Arwa Khashkhusha, Ariana Axiaq, Wing Yan Elizabeth Wong, and Amer Harky.  Health inequalities during COVID-19 and their effects on morbidity and mortality.  _Journal of Healthcare Leadership_ , Volume 13:19–26, 1 2021.  ISSN 1179-3201.  doi:10.2147/JHL.S270175.
  * Wu et al. [2020] Joseph T Wu, Kathy Leung, and Gabriel M Leung.  Nowcasting and forecasting the potential domestic and international spread of the 2019-nCoV outbreak originating in Wuhan, China: a modelling study.  _The Lancet_ , 395:689–697, 2 2020.  ISSN 01406736.  doi:10.1016/S0140-6736(20)30260-9.
  * Adiga et al. [2020] Aniruddha Adiga, Devdatt Dubhashi, Bryan Lewis, Madhav Marathe, Srinivasan Venkatramanan, and Anil Vullikanti.  Mathematical models for COVID-19 pandemic: A comparative analysis.  _Journal of the Indian Institute of Science_ , 100:793–807, 10 2020.  ISSN 0970-4140.  doi:10.1007/s41745-020-00200-6.
  * Rodríguez et al. [2021] Alexander Rodríguez, Anika Tabassum, Jiaming Cui, Jiajia Xie, Javen Ho, Pulak Agarwal, Bijaya Adhikari, and B. Aditya Prakash.  DeepCOVID: An operational deep learning-driven framework for explainable real-time COVID-19 forecasting.  _Proceedings of the AAAI Conference on Artificial Intelligence_ , 35:15393–15400, 5 2021.  ISSN 2374-3468.  doi:10.1609/aaai.v35i17.17808.
  * Zhou et al. [2020a] Tao Zhou, Quanhui Liu, Zimo Yang, Jingyi Liao, Kexin Yang, Wei Bai, Xin Lu, and Wei Zhang.  Preliminary prediction of the basic reproduction number of the Wuhan novel coronavirus 2019-nCoV.  _Journal of Evidence-Based Medicine_ , 13:3–7, 2 2020a.  ISSN 1756-5383.  doi:10.1111/jebm.12376.
  * Chang et al. [2021a] Meng-Chun Chang, Rebecca Kahn, Yu-An Li, Cheng-Sheng Lee, Caroline O. Buckee, and Hsiao-Han Chang.  Variation in human mobility and its impact on the risk of future COVID-19 outbreaks in Taiwan.  _BMC Public Health_ , 21:226, 12 2021a.  ISSN 1471-2458.  doi:10.1186/s12889-021-10260-7.
  * Kerr et al. [2021] Cliff C. Kerr, Robyn M. Stuart, Dina Mistry, Romesh G. Abeysuriya, Katherine Rosenfeld, Gregory R. Hart, Rafael C. Núñez, Jamie A. Cohen, Prashanth Selvaraj, Brittany Hagedorn, Lauren George, Michał Jastrzębski, Amanda S. Izzo, Greer Fowler, Anna Palmer, Dominic Delport, Nick Scott, Sherrie L. Kelly, Caroline S. Bennette, Bradley G. Wagner, Stewart T. Chang, Assaf P. Oron, Edward A. Wenger, Jasmina Panovska-Griffiths, Michael Famulare, and Daniel J. Klein.  Covasim: An agent-based model of COVID-19 dynamics and interventions.  _PLOS Computational Biology_ , 17:e1009149, 7 2021.  ISSN 1553-7358.  doi:10.1371/journal.pcbi.1009149.
  * Chang et al. [2021b] Serina Chang, Emma Pierson, Pang Wei Koh, Jaline Gerardin, Beth Redbird, David Grusky, and Jure Leskovec.  Mobility network models of COVID-19 explain inequities and inform reopening.  _Nature_ , 589:82–87, 1 2021b.  ISSN 0028-0836.  doi:10.1038/s41586-020-2923-3.
  * Kimathi et al. [2021] Mark Kimathi, Samuel Mwalili, Viona Ojiambo, and Duncan Kioi Gathungu.  Age-structured model for COVID-19: Effectiveness of social distancing and contact reduction in kenya.  _Infectious Disease Modelling_ , 6:15–23, 2021.  ISSN 24680427.  doi:10.1016/j.idm.2020.10.012.
  * Tuite et al. [2020] Ashleigh R. Tuite, David N. Fisman, and Amy L. Greer.  Mathematical modelling of COVID-19 transmission and mitigation strategies in the population of ontario, canada.  _Canadian Medical Association Journal_ , 192:E497–E505, 5 2020.  ISSN 0820-3946.  doi:10.1503/cmaj.200476.
  * Foy et al. [2021] Brody H. Foy, Brian Wahl, Kayur Mehta, Anita Shet, Gautam I. Menon, and Carl Britto.  Comparing COVID-19 vaccine allocation strategies in india: A mathematical modelling study.  _International Journal of Infectious Diseases_ , 103:431–438, 2 2021.  ISSN 12019712.  doi:10.1016/j.ijid.2020.12.075.
  * Lekone and Finkenstädt [2006] Phenyo E. Lekone and Bärbel F. Finkenstädt.  Statistical inference in a stochastic epidemic SEIR model with control intervention: Ebola as a case study.  _Biometrics_ , 62:1170–1177, 12 2006.  ISSN 0006-341X.  doi:10.1111/j.1541-0420.2006.00609.x.
  * Gray et al. [2011] A. Gray, D. Greenhalgh, L. Hu, X. Mao, and J. Pan.  A stochastic differential equation SIS epidemic model.  _SIAM Journal on Applied Mathematics_ , 71:876–902, 1 2011.  ISSN 0036-1399.  doi:10.1137/10081856X.
  * Faranda and Alberti [2020] Davide Faranda and Tommaso Alberti.  Modeling the second wave of COVID-19 infections in France and Italy via a stochastic SEIR model.  _Chaos: An Interdisciplinary Journal of Nonlinear Science_ , 30, 11 2020.  ISSN 1054-1500.  doi:10.1063/5.0015943.
  * Ganyani et al. [2021] Tapiwa Ganyani, Christel Faes, and Niel Hens.  Simulation and analysis methods for stochastic compartmental epidemic models.  _Annual Review of Statistics and Its Application_ , 8:69–88, 3 2021.  ISSN 2326-8298.  doi:10.1146/annurev-statistics-061120-034438.
  * Maltsev and Stern [2021] Alexander V. Maltsev and Michael D. Stern.  Social heterogeneity drives complex patterns of the COVID-19 pandemic: Insights from a novel stochastic heterogeneous epidemic model (SHEM).  _Frontiers in Physics_ , 8, 1 2021.  ISSN 2296-424X.  doi:10.3389/fphy.2020.609224.
  * Gatto et al. [2020] Marino Gatto, Enrico Bertuzzo, Lorenzo Mari, Stefano Miccoli, Luca Carraro, Renato Casagrandi, and Andrea Rinaldo.  Spread and dynamics of the COVID-19 epidemic in Italy: Effects of emergency containment measures.  _Proceedings of the National Academy of Sciences_ , 117:10484–10491, 5 2020.  ISSN 0027-8424.  doi:10.1073/pnas.2004978117.
  * Giordano et al. [2020] Giulia Giordano, Franco Blanchini, Raffaele Bruno, Patrizio Colaneri, Alessandro Di Filippo, Angela Di Matteo, and Marta Colaneri.  Modelling the COVID-19 epidemic and implementation of population-wide interventions in Italy.  _Nature Medicine_ , 26:855–860, 6 2020.  ISSN 1078-8956.  doi:10.1038/s41591-020-0883-7.
  * Kretzschmar et al. [2020] Mirjam E Kretzschmar, Ganna Rozhnova, Martin C J Bootsma, Michiel van Boven, Janneke H H M van de Wijgert, and Marc J M Bonten.  Impact of delays on effectiveness of contact tracing strategies for COVID-19: a modelling study.  _The Lancet Public Health_ , 5:e452–e459, 8 2020.  ISSN 24682667.  doi:10.1016/S2468-2667(20)30157-2.
  * Kucharski et al. [2020] Adam J Kucharski, Timothy W Russell, Charlie Diamond, Yang Liu, John Edmunds, Sebastian Funk, Rosalind M Eggo, Fiona Sun, Mark Jit, James D Munday, Nicholas Davies, Amy Gimma, Kevin van Zandvoort, Hamish Gibbs, Joel Hellewell, Christopher I Jarvis, Sam Clifford, Billy J Quilty, Nikos I Bosse, Sam Abbott, Petra Klepac, and Stefan Flasche.  Early dynamics of transmission and control of COVID-19: a mathematical modelling study.  _The Lancet Infectious Diseases_ , 20:553–558, 5 2020.  ISSN 14733099.  doi:10.1016/S1473-3099(20)30144-4.
  * WHO [2019] WHO.  Ten threats to global health in 2019, 2019.  https://www.who.int/news-room/spotlight/ten-threats-to-global-health-in-2019 (Accessed on 2121-06-14).
  * Choi and Shim [2020] Wongyeong Choi and Eunha Shim.  Optimal strategies for vaccination and social distancing in a game-theoretic epidemiologic model.  _Journal of Theoretical Biology_ , 505:110422, 11 2020.  ISSN 00225193.  doi:10.1016/j.jtbi.2020.110422.
  * Awad et al. [2021] Susanne F. Awad, Godfrey Musuka, Zindoga Mukandavire, Dillon Froass, Neil J. MacKinnon, and Diego F. Cuadros.  Implementation of a vaccination program based on epidemic geospatial attributes: COVID-19 pandemic in Ohio as a case study and proof of concept.  _Vaccines_ , 9:1242, 10 2021.  ISSN 2076-393X.  doi:10.3390/vaccines9111242.
  * Jing et al. [2023] Shuanglin Jing, Russell Milne, Hao Wang, and Ling Xue.  Vaccine hesitancy promotes emergence of new SARS-CoV-2 variants.  _Journal of Theoretical Biology_ , 570:111522, 8 2023.  ISSN 00225193.  doi:10.1016/j.jtbi.2023.111522.
  * Bui et al. [2024] Hieu Bui, Sandra Eksioglu, and Ruben Proano.  Evaluating the impact of vaccine hesitancy on the allocation of vital resources during COVID-19 pandemic, 03 2024.  URL arXiv:2403.07971.  (Accessed on 3/15/2024).
  * Abdin et al. [2023] Adam F. Abdin, Yi-Ping Fang, Aakil Caunhye, Douglas Alem, Anne Barros, and Enrico Zio.  An optimization model for planning testing and control strategies to limit the spread of a pandemic – the case of COVID-19.  _European Journal of Operational Research_ , 304:308–324, 1 2023.  ISSN 03772217.  doi:10.1016/j.ejor.2021.10.062.
  * Pereira and Pinto [1991] M. V. F. Pereira and L. M. V. G. Pinto.  Multi-stage stochastic optimization applied to energy planning.  _Mathematical Programming_ , 52:359–375, 5 1991.  ISSN 0025-5610.  doi:10.1007/BF01582895.
  * Goh et al. [2007] Mark Goh, Joseph Y.S. Lim, and Fanwen Meng.  A stochastic model for risk management in global supply chain networks.  _European Journal of Operational Research_ , 182:164–173, 10 2007.  ISSN 03772217.  doi:10.1016/j.ejor.2006.08.028.
  * Zhang and Cardin [2017] Sizhe Zhang and Michel-Alexandre Cardin.  Flexibility and real options analysis in emergency medical services systems using decision rules and multi-stage stochastic programming.  _Transportation Research Part E: Logistics and Transportation Review_ , 107:120–140, 11 2017.  ISSN 13665545.  doi:10.1016/j.tre.2017.09.003.
  * Li and Grossmann [2021] Can Li and Ignacio E. Grossmann.  A review of stochastic programming methods for optimization of process systems under uncertainty.  _Frontiers in Chemical Engineering_ , 2, 1 2021.  ISSN 2673-2718.  doi:10.3389/fceng.2020.622241.
  * Laporte and Louveaux [1993] Gilbert Laporte and François V. Louveaux.  The integer L-shaped method for stochastic integer programs with complete recourse.  _Operations Research Letters_ , 13:133–142, 4 1993.  ISSN 01676377.  doi:10.1016/0167-6377(93)90002-X.
  * Guignard [2003] Monique Guignard.  Lagrangean relaxation.  _Top_ , 11:151–200, 12 2003.  ISSN 1134-5764.  doi:10.1007/BF02579036.
  * Christian and Cremaschi [2015] Brianna Christian and Selen Cremaschi.  Heuristic solution approaches to the pharmaceutical R&D pipeline management problem.  _Computers & Chemical Engineering_, 74:34–47, 3 2015.  ISSN 00981354.  doi:10.1016/j.compchemeng.2014.12.014.
  * Yin and Buyuktahtakin [2021b] Xuecheng Yin and I. Esra Buyuktahtakin.  Risk-averse multi-stage stochastic programming to optimizing vaccine allocation and treatment logistics for effective epidemic response.  _medRxiv_ , 2021b.  doi:10.1101/2021.05.28.21258003.
  * Hosseini-Motlagh et al. [2023] Seyyed-Mahdi Hosseini-Motlagh, Mohammad Reza Ghatreh Samani, and Shamim Homaei.  Design of control strategies to help prevent the spread of COVID-19 pandemic.  _European Journal of Operational Research_ , 304:219–238, 1 2023.  ISSN 03772217.  doi:10.1016/j.ejor.2021.11.016.
  * Tanner et al. [2008] Matthew W. Tanner, Lisa Sattenspiel, and Lewis Ntaimo.  Finding optimal vaccination strategies under parameter uncertainty using stochastic programming.  _Mathematical Biosciences_ , 215:144–151, 10 2008.  ISSN 00255564.  doi:10.1016/j.mbs.2008.07.006.
  * Yarmand et al. [2014] Hamed Yarmand, Julie S. Ivy, Brian Denton, and Alun L. Lloyd.  Optimal two-phase vaccine allocation to geographically different regions under uncertainty.  _European Journal of Operational Research_ , 233:208–219, 2 2014.  ISSN 03772217.  doi:10.1016/j.ejor.2013.08.027.
  * WHO [2004] WHO.  Guidance on ethics and equitable access to HIV treatment and care.  https://www.who.int/publications/i/item/9241592583, June 2004.  (Accessed on 03/15/2024).
  * Lane et al. [2017] Haylee Lane, Mitchell Sarkies, Jennifer Martin, and Terry Haines.  Equity in healthcare resource allocation decision making: A systematic review.  _Social Science & Medicine_, 175:11–27, 2 2017.  ISSN 02779536.  doi:10.1016/j.socscimed.2016.12.012.
  * Sun et al. [2014] Li Sun, Gail W. DePuy, and Gerald W. Evans.  Multi-objective optimization models for patient allocation during a pandemic influenza outbreak.  _Computers & Operations Research_, 51:350–359, 11 2014.  ISSN 03050548.  doi:10.1016/j.cor.2013.12.001.
  * Liu et al. [2015] Ming Liu, Zhe Zhang, and Ding Zhang.  A dynamic allocation model for medical resources in the control of influenza diffusion.  _Journal of Systems Science and Systems Engineering_ , 24:276–292, 9 2015.  ISSN 1004-3756.  doi:10.1007/s11518-015-5276-y.
  * Liu and Zhang [2016] Ming Liu and Ding Zhang.  A dynamic logistics model for medical resources allocation in an epidemic control with demand forecast updating.  _Journal of the Operational Research Society_ , 67:841–852, 6 2016.  ISSN 0160-5682.  doi:10.1057/jors.2015.105.
  * Zaric and Brandeau [2001] Gregory S. Zaric and Margaret L. Brandeau.  Resource allocation for epidemic control over short time horizons.  _Mathematical Biosciences_ , 171:33–58, 5 2001.  ISSN 00255564.  doi:10.1016/S0025-5564(01)00050-5.
  * Day et al. [2020] R. Thomas Day, Bradley S. Guidry, Brian C. Drolet, and Ellen W. Clayton.  From ventilators to vaccines: Reframing the ethics of resource allocation.  _The American Journal of Bioethics_ , 20:W15–W16, 7 2020.  ISSN 1526-5161.  doi:10.1080/15265161.2020.1782530.
  * Orgut et al. [2016] Irem Sengul Orgut, Julie Ivy, Reha Uzsoy, and James R. Wilson.  Modeling for the equitable and effective distribution of donated food under capacity constraints.  _IIE Transactions_ , 48:252–266, 3 2016.  ISSN 0740-817X.  doi:10.1080/0740817X.2015.1063792.
  * Davis et al. [2015] Ashley E. Davis, Sanjay Mehrotra, John J. Friedewald, Mark S. Daskin, Anton I. Skaro, Michael M. Abecassis, and Daniela P. Ladner.  Improving geographic equity in kidney transplantation using alternative kidney sharing and optimization modeling.  _Medical Decision Making_ , 35:797–807, 8 2015.  ISSN 0272-989X.  doi:10.1177/0272989X14557696.
  * Dönmez et al. [2022] Zehranaz Dönmez, Serkan Turhan, Özlem Karsu, Bahar Y. Kara, and Oya Karaşan.  Fair allocation of personal protective equipment to health centers during early phases of a pandemic.  _Computers & Operations Research_, 141:105690, 5 2022.  ISSN 03050548.  doi:10.1016/j.cor.2021.105690.
  * Salomon et al. [2021] Joshua A. Salomon, Alex Reinhart, Alyssa Bilinski, Eu Jing Chua, Wichada La Motte-Kerr, Minttu M. Rönn, Marissa B. Reitsma, Katherine A. Morris, Sarah LaRocca, Tamer H. Farag, Frauke Kreuter, Roni Rosenfeld, and Ryan J. Tibshirani.  The US COVID-19 Trends and Impact Survey: Continuous real-time measurement of COVID-19 symptoms, risks, protective behaviors, testing, and vaccination.  _Proceedings of the National Academy of Sciences_ , 118(51):e2111454118, 2021.
  * CDC [2021] CDC.  COVID-19 vaccinations in the United States, county, 2021.  URL https://data.cdc.gov/Vaccinations/COVID-19-Vaccinations-in-the-United-States-County/8xkx-amqh.  (Accessed on 2021-12-25).
  * Wacker and Schlüter [2020] Benjamin Wacker and Jan Schlüter.  Time-continuous and time-discrete SIR models revisited: theory and applications.  _Advances in Difference Equations_ , 2020(1):556, 2020.
  * Dangerfield et al. [2019] Ciara E. Dangerfield, Martin Vyska, and Christopher A. Gilligan.  Resource allocation for epidemic control across multiple sub-populations.  _Bulletin of Mathematical Biology_ , 81:1731–1759, 6 2019.  ISSN 0092-8240.  doi:10.1007/s11538-019-00584-2.
  * World Health Organization [2020] World Health Organization.  A global framework to ensure equitable and fair allocation of COVID-19 products and potential implications for COVID-19 vaccination.  Technical report, World Health Organization, 2020.  URL https://apps.who.int/gb/COVID-19/pdf_files/02_07/Global_Allocation_Framework.pdf.  Accessed: 2024-02-24.
  * CovidCareMap [2020] CovidCareMap.  COVID-19 health system capacity.  https://github.com/covidcaremap/covid19-healthsystemcapacity/tree/master/data/published, 2020.  URL https://github.com/covidcaremap/covid19-healthsystemcapacity/tree/master/data/published.  Accessed: 2024-01-12.
  * Virtanen et al. [2020] Pauli Virtanen, Ralf Gommers, Travis E. Oliphant, Matt Haberland, Tyler Reddy, David Cournapeau, Evgeni Burovski, Pearu Peterson, Warren Weckesser, Jonathan Bright, Stéfan J. van der Walt, Matthew Brett, Joshua Wilson, K. Jarrod Millman, Nikolay Mayorov, Andrew R. J. Nelson, Eric Jones, Robert Kern, Eric Larson, C J Carey, İlhan Polat, Yu Feng, Eric W. Moore, Jake VanderPlas, Denis Laxalde, Josef Perktold, Robert Cimrman, Ian Henriksen, E. A. Quintero, Charles R. Harris, Anne M. Archibald, Antônio H. Ribeiro, Fabian Pedregosa, Paul van Mulbregt, and SciPy 1.0 Contributors.  SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python.  _Nature Methods_ , 17:261–272, 2020.  doi:10.1038/s41592-019-0686-2.
  * Massey and Cole [2020] Delano Massey and Devan Cole.  HHS to work with GM under Defense Production Act to produce 30,000 ventilators for national stockpile | cnn politics, 04 2020.  URL https://www.cnn.com/2020/04/08/politics/general-motors-ventilators-defense-production-act-coronavirus/index.html.  (Accessed on 10/06/2023).
  * AdvaMed [2020] AdvaMed.  Tens of thousands of life-saving ventilators on deck to combat COVID-19 - advamed, 04 2020.  URL https://www.advamed.org/industry-updates/news/tens-of-thousands-of-life-saving-ventilators-on-deck-to-combat-covid-19/.  (Accessed on 11/06/2023).
  * Wahltinez et al. [2020] O. Wahltinez et al.  COVID-19 Open-Data: curating a fine-grained, global-scale data repository for SARS-CoV-2.  2020\.  URL https://goo.gle/covid-19-open-data.
  * Akiba et al. [2019] Takuya Akiba, Shotaro Sano, Toshihiko Yanase, Takeru Ohta, and Masanori Koyama.  Optuna: A next-generation hyperparameter optimization framework.  In _Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining_ , 2019.
  * Escudero et al. [2007] Laureano Escudero, Araceli Garín, María Merino, and Gloria Pérez.  The value of the stochastic solution in multistage problems.  _TOP: An Official Journal of the Spanish Society of Statistics and Operations Research_ , 15:48–64, 02 2007.
  * U.S. Census Bureau [2020] U.S. Census Bureau.  County-to-county migration flows: 2016-2020 acs.  https://www.census.gov/data/tables/2020/demo/geographic-mobility/county-to-county-migration-2016-2020.html, 2020.  Accessed: 2024-01-12.
  * Zhou et al. [2020b] Fei Zhou, Ting Yu, Ronghui Du, Guohui Fan, Ying Liu, Zhibo Liu, Jie Xiang, Yeming Wang, Bin Song, Xiaoying Gu, Lulu Guan, Yuan Wei, Hui Li, Xudong Wu, Jiuyang Xu, Shengjin Tu, Yi Zhang, Hua Chen, and Bin Cao.  Clinical course and risk factors for mortality of adult inpatients with COVID-19 in Wuhan, China: a retrospective cohort study.  _The Lancet_ , 395:1054–1062, 3 2020b.  ISSN 01406736.  doi:10.1016/S0140-6736(20)30566-3.
  * Wu and McGoogan [2020] Zunyou Wu and Jennifer M. McGoogan.  Characteristics of and important lessons from the Coronavirus disease 2019 (COVID-19) outbreak in China.  _JAMA_ , 323:1239, 4 2020.  ISSN 0098-7484.  doi:10.1001/jama.2020.2648.
  * Barnes and Others [2023] Eleanor Barnes and Others.  SARS-CoV-2-specific immune responses and clinical outcomes after COVID-19 vaccination in patients with immune-suppressive disease.  _Nature Medicine_ , 29:1760–1774, 7 2023.  ISSN 1078-8956.  doi:10.1038/s41591-023-02414-4.
  * Katella [2023] Kathy Katella.  Comparing the COVID-19 vaccines: How are they different?  https://www.yalemedicine.org/news/covid-19-vaccine-comparison, 10 2023.  (Accessed on 01/13/2024).

OR

    operations research
SIR

    susceptible-infectious-recovered
SEIR

    susceptible-exposed-infectious-recovered
SVEIHR

    susceptible-vaccinated-exposed-infectious-hospitalized-recovered
SP

    stochastic programming
PPE

    personal protective equipment
ICU

    intensive care unit
WHO

    World Health Organization
MSP

    multi-stage stochastic program
VH

    vaccine hesitancy
CTIS

    COVID-19 Trends and Impact Survey
ODE

    ordinary differential equation
DTW

    dynamic time warping
CDC

    Centers for Disease Control and Prevention
NAC

    non-anticipativity constraint
RMSE

    root-mean-square error
MAPE

    mean absolute percentage error
MAE

    mean absolute error
GM

    General Motors
DPA

    Defense Production Act
HHS

    U.S. Department of Health and Human Services
VSS

    value of stochastic solution
EV

    expected value
  *[VH]: vaccine hesitancy
  *[MSP]: multi-stage stochastic program
  *[SIR]: susceptible-infectious-recovered
  *[SEIR]: susceptible-exposed-infectious-recovered
  *[SP]: stochastic programming
  *[DPA]: Defense Production Act
  *[WHO]: World Health Organization
  *[PPE]: personal protective equipment
  *[CTIS]: COVID-19 Trends and Impact Survey
  *[SVEIHR]: susceptible-vaccinated-exposed-infectious-hospitalized-recovered
  *[ICU]: intensive care unit
  *[GM]: General Motors
  *[RMSE]: root-mean-square error
  *[MAE]: mean absolute error
  *[MAPE]: mean absolute percentage error
  *[VSS]: value of stochastic solution
  *[EV]: expected value
  *[NACs]: non-anticipativity constraint