3.1 Metric Temporal Logic (MTL)

In this subsection, we briefly review metric temporal logic (MTL) [44] interpreted over discrete-time trajectories. The state x belongs to the domain $$ \mathcal{X}\subset {\mathbb{R}}^n$$. The time set is $$ \mathbb{T}={\mathbb{R}}_{\ge 0}$$. The domain $$ \mathbb{B}=\{\text{True},\text{False}\}$$ is the Boolean domain, and the time index set is $$ \mathbb{I}=\{0,1,\dots \}$$. We use $$ t\left[k\right]\in \mathbb{T}$$ to denote the time instant at time index $$ k\in \mathbb{I}$$ and $$ x\left[k\right]\stackrel{▵}{=}x\left(t\right[k\left]\right)$$ to denote the value of x at time instant t[k]. We use ξ to denote a trajectory as a function from $$ \mathbb{T}$$ to $$ \mathcal{X}$$. A set AP is a set of atomic propositions, each mapping $$ \mathcal{X}$$ to $$ \mathbb{B}$$. For example, the state x can be in the form of [x₁, x₂], where x₁ and x₂ represent the number of deaths from COVID-19 and the number of recovered patients from COVID-19 in a certain geographic region, respectively. Then, an atomic proposition could be in the form of (x₁ ≤ 0.01), which means “the number of deaths from COVID-19 in the region does not exceed 0.01 million”, or in the form of (x₂ ≥ 6), which means “the number of recovered patients from COVID-19 in the region is at least 6 million” (in this paper we assume that the unit in the state is million, unless otherwise indicated). In such context, if x₁ = 0.002 and x₂ = 4, then (x₁ ≤ 0.01) is True, and (x₂ ≥ 6) is False.

The syntax of MTL is defined recursively as follows:

We define the set of states that satisfy the atomic proposition π as $$ \mathcal{O}\left(\pi \right)\subset \mathcal{X}$$. We denote ⟨⟨φ⟩⟩(ξ, k) = ⊤ if the trajectory ξ satisfies the formula φ at discrete-time instant t[k] ($$ k\in \mathbb{I}$$). Then the Boolean semantics of MTL are defined recursively as follows [45]:

3.2 COVID-19 models with control strategies

In this subsection, we study three models for COVID-19 epidemic [6, 7, 14] and introduce the corresponding models with vaccination control, shield immunity control and quarantine control.

COVID-19 SEIR model with vaccination control. The susceptible, exposed, infectious, recovered (SEIR) model has been frequently used in epidemic analyses. As shown in Fig 1, the total population is divided into five subgroups:

The susceptible population S: everyone is susceptible to the disease by birth since immunity is not hereditary;

The exposed population E: the individuals who have been exposed to the disease, but are still not infectious;

The infectious population I: the individuals who are infectious;

The immune (recovered) population R: the individuals who are vaccinated or recovered from the disease, i.e., the population who are immune to the disease;

The dead population D: the dead individuals from the disease.

Fig 1

Block diagram of the COVID-19 SEIR model with vaccination control.

We consider a COVID-19 SEIR model [14, 15] with vaccination control [19] as follows.

Remark 1. Note that one difference between this model and the vaccination control model in [19] is that we control V as the number of vaccinated individuals per day (constrained to be less than the susceptible population S), while in [19] the control input is the ratio of the vaccinated individuals per day to the average born population per day. We found it more convenient this way for computational convenience in the control synthesis in Section 3.3.

COVID-19 SEIR model with shield immunity control. Shield immunity is a strategy recently proposed in [6] to limit the transmission of SARS-CoV-2. The basic idea of this strategy is to increase the proportion of interactions with recovered individuals as opposed to the other individuals in the population. The effectiveness of this strategy is based on the assumption that recovered individuals (virus-negative and antibody-positive) can safely interact with both susceptible and infectious individuals without getting infected with the disease.

As the model used in [6] is modified from an SIR model, we consider a corresponding SEIR model with shield immunity control as follows (see Fig 2 as an illustration).

Fig 2

Block diagram of the COVID-19 SEIR model with shield immunity control.

COVID-19 SUQC model with quarantine control. The susceptible, un-quarantined infected, quarantined infected, confirmed infected (SUQC) model was recently proposed in [7] based on the COVID-19 data in Wuhan, China. As shown in Fig 3, we consider four subgroups in the population:

The susceptible population S: everyone is susceptible to the disease by birth since immunity is not hereditary;

The un-quarantined infected population U: the individuals who are infected and un-quarantined, and they can be either asymptomatic or symptomatic;

The quarantined infected population Q: the individuals who are infectious and quarantined (the un-quarantined infected become quarantined infected by isolation or hospitalization, and the quarantined infected lose the ability of infecting the susceptible individuals);

The confirmed infected population C: the individuals who are confirmed to be infected with the disease (i.e., the positive cases).

Fig 3

Block diagram of the COVID-19 SUQC model with quarantine control.

We consider the SUQC model with quarantine control as follows.

3.3 Control synthesis of COVID-19 epidemic with metric temporal logic specifications

In this subsection, we present the control synthesis methods for the three COVID-19 epidemic models in Section 3.2 with vaccination control, shield immunity control and quarantine control, respectively.

Vaccination control. For the COVID-19 SEIR model with vaccination control, we discretize the model in (3) as follows.

Following the notations in Section 3.1, we use x_V = [I, E, S, R, D] to denote the state of (6) and $$ {\xi }_{·;x_{\mathrm{V}}^{init},V}^{\mathrm{V}}$$ to denote the trajectory of (6) starting from $$ x_{\mathrm{V}}^{init}=[I\left[0\right],E\left[0\right],S\left[0\right],R\left[0\right],D\left[0\right]]$$ and vaccination control input V[⋅].

Problem 1 (Vaccination control). Given the SEIR model with vaccination control in (6) and an MTL specification φ_V, compute the control input V[⋅] that minimizes the vaccination control efforts ‖V[⋅]‖ while satisfying $$ ⟨⟨{\phi }_{\mathrm{V}}⟩⟩({\xi }_{·;x_{\mathrm{V}}^{init},V}^{\mathrm{V}},0)=\top $$, i.e., the trajectory $$ {\xi }_{·;x_{\mathrm{V}}^{init},V}^{\mathrm{V}}$$ satisfies the MTL specification φ_V.

The vaccination control synthesis problem can be formulated as a constrained optimization problem as follows.

The above optimization problem is generally a mixed-integer non-linear programming problem. We refer the readers to [38] for a detailed description of how the constraint $$ ⟨⟨{\phi }_{\mathrm{V}}⟩⟩({\xi }_{·;x_{\mathrm{V}}^{init},V}^{\mathrm{V}},0)=\top $$ is encoded to satisfy an MTL specification φ_V. The integer variables are introduced when a big-M formulation [46] is needed to satisfy MTL specifications such as ♢_[0,10] φ (φ should hold true for at least one day during the first 10 days) or φ₁ ∨ φ₂ (at least one of the MTL formulas φ₁, φ₂ should hold true). As the change of total population is relatively small compared to the multiplication of the susceptible population and the infectious population, we approximate the term T_s βS[k]I[k]/N[k] with T_s βS[k]I[k]/N₀. With such an approximation, the optimization problem becomes a mixed-integer bi-linear programming problem, which can be more efficiently solved using techniques such as McCormick’s relaxation [47, 48]. Furthermore, if the MTL specification φ consists of only conjunctions (∧) and the always operator (□), the integers in the optimization problem can be eliminated [38] and the problem becomes a bi-linear programming problem.

Shield immunity control. For the COVID-19 SEIR model with shield immunity control, we discretize the model in (4) as follows.

Following the notations in Section 3.1, we use x_S = [I, E, S, R, D] to denote the state of (8) and $$ {\xi }_{·;x_{\mathrm{S}}^{init},\chi }^{\mathrm{S}}$$ to denote the trajectory of (8) starting from $$ x_{\mathrm{S}}^{init}=[I\left[0\right],E\left[0\right],S\left[0\right],R\left[0\right],D\left[0\right]]$$ and shield immunity control input χ[⋅].

Problem 2 (Shield immunity control). Given the SEIR model with shield immunity control in (8) and an MTL specification φ_S, compute the control input χ[⋅] that minimizes the shield immunity control efforts ‖χ[⋅]‖ while satisfying $$ ⟨⟨{\phi }_{\mathrm{S}}⟩⟩({\xi }_{·;x_{\mathrm{S}}^{init},\chi }^{\mathrm{S}},0)=\top $$, i.e., the trajectory $$ {\xi }_{·;x_{\mathrm{S}}^{init},\chi }^{\mathrm{S}}$$ satisfies the MTL specification φ_S.

The shield immunity control synthesis problem can be formulated as a constrained optimization problem as follows.

The above optimization problem is generally a mixed-integer fractional constrained programming problem. If the MTL specification φ consists of only conjunctions (∧) and the always operator (□), the integers in the optimization problem can be eliminated [38] and the problem becomes a fractional constrained programming problem.

Quarantine control. For the COVID-19 SUQC model with quarantine control, we discretize the model in (5) as follows.

Following the notations in Section 3.1, we use x_Q = [S, U, Q, C] to denote the state of (10) and $$ {\xi }_{·;x_{\mathrm{Q}}^{init},q}^{\mathrm{Q}}$$ to denote the trajectory of (10) starting from $$ x_{\mathrm{Q}}^{init}=[S\left[0\right],U\left[0\right],Q\left[0\right],C\left[0\right]]$$ and quarantine control input q[⋅].

Problem 3 (Quarantine control). Given the SUQC model with quarantine control in (10) and an MTL specification φ_Q, compute the control input q[⋅] that minimizes the quarantine control efforts ‖q[⋅]‖ while satisfying $$ ⟨⟨{\phi }_{\mathrm{Q}}⟩⟩({\xi }_{·;x_{\mathrm{Q}}^{init},q}^{\mathrm{Q}},0)=\top $$, i.e., the trajectory $$ {\xi }_{·;x_{\mathrm{Q}}^{init},q}^{\mathrm{Q}}$$ satisfies the MTL specification φ_Q.

The quarantine control synthesis problem can be formulated as a constrained optimization problem as follows.

The above optimization problem is generally a mixed-integer non-linear programming problem. As the change of total population is relatively small compared to the multiplication of the susceptible population and the un-quarantined infectious population, we approximate the term T_s β₀ U[k]S[k]/N[k] with $$ T_{\mathrm{s}}{\beta }_{0}U\left[k\right]S\left[k\right]/{\widehat{N}}_0$$ (we use $$ {\widehat{N}}_0$$ to denote the initial population in the region in the scenario with quarantine control). With such an approximation, the optimization problem becomes a mixed-integer bi-linear programming problem, which can be more efficiently solved using techniques such as McCormick’s relaxation [47, 48]. Furthermore, if the MTL specification φ consists of only conjunctions (∧) and the always operator (□), the integers in the optimization problem can be eliminated [38] and the problem becomes a bi-linear programming problem.

4.1 COVID-19 SEIR model with vaccination control

The parameters of the COVID-19 SEIR model are shown in Table 1. They were estimated using the data from February 23 to March 16, 2020 in Lombardy, Italy with no isolation measures (see Section 4.1 in [14] for details). Specifically, in Table 1, the values of λ and μ are from “μ⁻¹ ≈ 83 years” (as 83 years = 30295 days) in Section 4 in [14] and our assumption that λ = μ, the values of α, β, ϵ and γ are from the initial values listed in Table 1 in [14], and the value of N₀ is from Lombardy’s population listed in Wikipedia (https://en.wikipedia.org/wiki/Lombardy). The start time for the simulations in this subsection is February 23, 2020. We consider three MTL specifications as shown in Table 2. For example, $$ {\phi }_{\mathrm{V}}^1={\square }_{[0,100]}(\Delta D\le 0.001)\wedge {\square }_{[0,100]}(D\le 0.05)\wedge {♢}_{[40,60]}(R\ge 6)$$, which means “the deaths from COVID-19 should never exceed 0.001 million (i.e., one thousand) per day and 0.05 million (i.e., 50 thousand) in total within the next 100 days, and the immune population should eventually exceed 6 million within the next 40 to 60 days”. Following Section 4.1 in [14, we choose the initial values of the states as I[0] = 1000 (people), E[0] = 0.02 million, S[0] = 9.979 million, R[0] = 0 and D[0] = 0, with S[0] + E[0] + I[0] + R[0] + D[0] = N₀ = 10 million. Fig 4 shows the simulation results without any vaccination. It can be seen that the three MTL specifications $$ {\phi }_{\mathrm{V}}^1$$, $$ {\phi }_{\mathrm{V}}^2$$ and $$ {\phi }_{\mathrm{V}}^3$$ are all violated in such a situation. Note that as isolation measures (i.e., home isolation, social distancing and partial national lockdown) were taken since March 16 in Lombardy, Italy, the real situation was better than those shown in Fig 4. Now we investigate the hypothetical scenario where the isolation measures are replaced by vaccination.

Fig 4

Simulation results for COVID-19 SEIR model estimated from data from Lombardy, Italy with no isolation measures.

Table 1

Parameters of COVID-19 SEIR model estimated from data from Lombardy, Italy from February 23 to March 16 (2020) with no isolation measures [14].

parameter	value	parameter	value
λ	1/30295	ϵ	0.2/day
μ	1/30295	γ	0.2/day
α	0.006/day	N0	10 million
β	0.75/day	Ts	1 day

Table 2

MTL specifications and simulation results for vaccination control (Section 4.1).

MTL specification	control effort	computation time
$$	1.28	1.365 s
$$	1.927	2.276 s
$$	6.934	3.289 s

We use the solver GEKKO [11] to solve the optimization problems formulated in Section 3.3. Fig 5 and Table 2 show the simulation results for vaccination control of COVID-19 SEIR model with MTL specifications $$ {\phi }_{\mathrm{V}}^1$$, $$ {\phi }_{\mathrm{V}}^2$$ and $$ {\phi }_{\mathrm{V}}^3$$, respectively. The results show that the MTL specifications $$ {\phi }_{\mathrm{V}}^1$$, $$ {\phi }_{\mathrm{V}}^2$$ and $$ {\phi }_{\mathrm{V}}^3$$ are satisfied with the synthesized vaccination control inputs respectively. It can be seen that vaccination within the first 40 days after the outbreak can mitigate the spread of SARS-CoV-2 in the most efficient manner. The results also show that the control effort for satisfying $$ {\phi }_{\mathrm{V}}^1$$ is less than that for satisfying $$ {\phi }_{\mathrm{V}}^2$$, which is still less than that for satisfying $$ {\phi }_{\mathrm{V}}^3$$. This is consistent with the fact that $$ {\phi }_{\mathrm{V}}^2$$ implies $$ {\phi }_{\mathrm{V}}^1$$, and $$ {\phi }_{\mathrm{V}}^3$$ implies both $$ {\phi }_{\mathrm{V}}^1$$ and $$ {\phi }_{\mathrm{V}}^2$$. For all three specifications, the computations are completed within 4 seconds on a MacBook laptop with 1.40-GHz Core i5 CPU and 16-GB RAM.

Fig 5

Simulation results for the COVID-19 SEIR model with vaccination control and MTL specifications $$ {\phi }_{\mathrm{V}}^1$$ (first row), $$ {\phi }_{\mathrm{V}}^2$$ (second row) and $$ {\phi }_{\mathrm{V}}^3$$ (third row).

The red dotted lines indicate the thresholds in the atomic propositions of the MTL specifications $$ {\phi }_{\mathrm{V}}^1$$, $$ {\phi }_{\mathrm{V}}^2$$ and $$ {\phi }_{\mathrm{V}}^3$$.

4.2 COVID-19 SEIR model with shield immunity control

We use the same parameters of the COVID-19 SEIR model as shown in Table 1 (in the same setting as in Section 4.1 for Lombardy, Italy, with no isolation measures). Following Section 4.1 in [14], we choose the initial values of the states as I[0] = 1000 (people), E[0] = 0.02 million, S[0] = 9.979 million, R[0] = 0 and D[0] = 0, with S[0] + E[0] + I[0] + R[0] + D[0] = N₀ = 10 million. We set χ_max = 100. The start time for the simulations in this subsection is February 23, 2020. We set the three MTL specifications $$ {\phi }_{\mathrm{S}}^1$$, $$ {\phi }_{\mathrm{S}}^2$$ and $$ {\phi }_{\mathrm{S}}^3$$ (as shown in Table 3) to be less stringent than the MTL specifications with the vaccination control, as shield immunity is generally less effective than vaccination. It can be shown that without any control strategies the three MTL specifications $$ {\phi }_{\mathrm{S}}^1$$, $$ {\phi }_{\mathrm{S}}^2$$ and $$ {\phi }_{\mathrm{S}}^3$$ are all violated. Now we investigate the hypothetical scenario where the isolation measures are replaced by shield immunity control.

Table 3

MTL specifications and simulation results for shield immunity control (Section 4.2).

MTL specification	control effort	computation time
$$	16879.53	2.112 s
$$	45595.10	2.881 s
$$	67786.88	5.323 s

Fig 6 and Table 3 show the simulation results for shield immunity control of the COVID-19 SEIR model with MTL specifications $$ {\phi }_{\mathrm{S}}^1$$, $$ {\phi }_{\mathrm{S}}^2$$ and $$ {\phi }_{\mathrm{S}}^3$$, respectively. The results show that the MTL specifications $$ {\phi }_{\mathrm{S}}^1$$, $$ {\phi }_{\mathrm{S}}^2$$ and $$ {\phi }_{\mathrm{S}}^3$$ are satisfied with the synthesized shield immunity control inputs respectively. We observe that with the three MTL specifications, the synthesized shield immunity control inputs all increase to a peak after approximately 20 to 40 days and then gradually decrease. These observations indicate that shield immunity at early days of COVID-19 is more efficient than shield immunity at later days. The results also show that the control effort for satisfying $$ {\phi }_{\mathrm{S}}^1$$ is less than that for satisfying $$ {\phi }_{\mathrm{S}}^2$$, which is still less than that for satisfying $$ {\phi }_{\mathrm{S}}^3$$. This is consistent with the fact that $$ {\phi }_{\mathrm{S}}^2$$ implies $$ {\phi }_{\mathrm{S}}^1$$, and $$ {\phi }_{\mathrm{S}}^3$$ implies both $$ {\phi }_{\mathrm{S}}^1$$ and $$ {\phi }_{\mathrm{S}}^2$$.

Fig 6

Simulation results for the COVID-19 SEIR model with shield immunity control and MTL specifications $$ {\phi }_{\mathrm{S}}^1$$ (first row), $$ {\phi }_{\mathrm{S}}^2$$ (second row) and $$ {\phi }_{\mathrm{S}}^3$$ (third row).

The red dotted lines indicate the thresholds in the atomic propositions of the MTL specifications $$ {\phi }_{\mathrm{S}}^1$$, $$ {\phi }_{\mathrm{S}}^2$$ and $$ {\phi }_{\mathrm{S}}^3$$.

4.3 COVID-19 SUQC model with quarantine control

The parameters of the COVID-19 SUQC model are shown in Table 4. In Table 4, the value of the infection rate β₀ is from “Methods/Parameter inference” in [7] (the authors in [7] used α to denote the infection rate), the value of $$ {\widehat{N}}_0$$ is from Wuhan’s urban population listed in Wikipedia (https://en.wikipedia.org/wiki/Wuhan), the value of the confirmation rate γ₂ is from the confirmation rate listed under Stage I (January 20 to January 30, 2020) of Wuhan in Table 1 in [7], and the value of σ is from “Methods/SUQC model” in [7] when no other special approaches are used for additional tests. We choose the initial values of the states as S[0] = 8.9 million, U[0] = 0.001 million, Q[0] = 0 and C[0] = 0. We set q_max = 1. We consider three MTL specifications as shown in Table 5. For example, $$ {\phi }_{\mathrm{Q}}^1={\square }_{[0,200]}(\Delta C\le 0.001)\wedge {\square }_{[0,200]}(C\le 0.1)$$ means “the confirmed infected population should never exceed 0.001 million (i.e., one thousand) per day and 0.1 million (i.e., 100 thousand) in total within the next 200 days”. The start time for the simulations in this subsection is January 20, 2020. Fig 7 shows the simulation results for the COVID-19 SUQC model with parameters in Table 4 (estimated from data in Stage I of Wuhan, China). It can be seen that the three MTL specifications $$ {\phi }_{\mathrm{Q}}^1$$, $$ {\phi }_{\mathrm{Q}}^2$$ and $$ {\phi }_{\mathrm{Q}}^3$$ are all violated in such a situation (with quarantine rate being always 0.063). Now we investigate the scenario where the quarantine rate can be controlled to satisfy the MTL specifications.

Fig 7

Simulation results for the COVID-19 SUQC model estimated from data in Stage I of Wuhan, China.

Table 4

Parameters of the COVID-19 SUQC model estimated from data in Stage I (January 20 to January 30, 2020) of Wuhan, China [7].

parameter	value	parameter	value
β0	0.2967	γ2	0.05
$$	8.9 million	σ	0
Ts	1 day

Table 5

MTL specifications and simulation results for quarantine control (Section 4.3).

MTL specification	control effort	computation time
$$	15.146	2.296 s
$$	15.638	2.598 s
$$	15.894	4.578 s

Fig 8 and Table 5 show the simulation results for quarantine control of the COVID-19 SUQC model with MTL specifications $$ {\phi }_{\mathrm{Q}}^1$$, $$ {\phi }_{\mathrm{Q}}^2$$ and $$ {\phi }_{\mathrm{Q}}^3$$, respectively. The results show that the MTL specifications $$ {\phi }_{\mathrm{Q}}^1$$, $$ {\phi }_{\mathrm{Q}}^2$$ and $$ {\phi }_{\mathrm{Q}}^3$$ are satisfied with the synthesized quarantine control inputs respectively. The results also show that the control effort for satisfying $$ {\phi }_{\mathrm{Q}}^1$$ is less than that for satisfying $$ {\phi }_{\mathrm{Q}}^2$$, which is still less than that for satisfying $$ {\phi }_{\mathrm{Q}}^3$$. This is consistent with the fact that $$ {\phi }_{\mathrm{Q}}^2$$ implies $$ {\phi }_{\mathrm{Q}}^1$$, and $$ {\phi }_{\mathrm{Q}}^3$$ implies both $$ {\phi }_{\mathrm{Q}}^1$$ and $$ {\phi }_{\mathrm{Q}}^2$$. We observe that with $$ {\phi }_{\mathrm{Q}}^1$$, the synthesized quarantine control inputs first increase to a peak at approximately 90 days and then gradually decrease for most of the time; with $$ {\phi }_{\mathrm{Q}}^2$$, the synthesized quarantine control inputs first increase to a peak at approximately 50 days and then gradually decrease for most of the time; and with $$ {\phi }_{\mathrm{Q}}^3$$, the synthesized quarantine control inputs are at a peak from the beginning and gradually decrease for most of the time. These observations indicate that quarantine in the early days of COVID-19 can reduce the number of confirmed infected cases more efficiently than quarantine in the later days, and more stringent control specifications require stronger quarantine measures to be implemented.

Fig 8

Simulation results for the COVID-19 SUQC model with quarantine control and MTL specifications $$ {\phi }_{\mathrm{Q}}^1$$ (first row), $$ {\phi }_{\mathrm{Q}}^2$$ (second row) and $$ {\phi }_{\mathrm{Q}}^3$$ (third row).

The red dotted lines indicate the thresholds in the atomic propositions of the MTL specifications $$ {\phi }_{\mathrm{Q}}^1$$, $$ {\phi }_{\mathrm{Q}}^2$$ and $$ {\phi }_{\mathrm{Q}}^3$$.