# Handling plant-model mismatch in Koopman Lyapunov-based model predictive
control via offset-free control framework

Sang Hwan Son<EMAIL_ADDRESS>Abhinav Narasingam<EMAIL_ADDRESS>Joseph Sang-Il Kwon<EMAIL_ADDRESS>Artie McFerrin Department of Chemical
Engineering, Texas A&M University, College Station, TX 77845 USA  Texas A&M
Energy Institute, Texas A&M University, College Station, TX 77845 USA

###### Abstract

Koopman operator theory enables a global linear representation of a given
nonlinear dynamical system by transforming the nonlinear dynamics into a
higher dimensional observable function space where the evolution of observable
functions is governed by an infinite-dimensional linear operator. For
practical application of Koopman operator theory, various data-driven methods
have been developed to derive lifted state-space models via approximation to
the Koopman operator. Based on approximate models, several Koopman-based model
predictive control (KMPC) schemes have been proposed. However, since a finite-
dimensional approximation to the infinite-dimensional Koopman operator cannot
fully represent a nonlinear dynamical system, plant-model mismatch inherently
exists in these KMPC schemes and negatively influences the performance of
control systems. In this work, we present offset-free Koopman Lyapunov-based
model predictive control (KLMPC) framework that addresses the inherent plant-
model mismatch in KMPC schemes using an offset-free control framework based on
a disturbance estimator approach and ensures feasibility and stability of the
control system by applying Lyapunov constraints to the optimal control
problem. The zero steady-state offset condition of the developed framework is
mathematically examined. The effectiveness of the developed framework is also
demonstrated by comparing the closed-loop results of the proposed offset-free
KLMPC and the nominal KLMPC.

###### keywords:

Koopman operator, extended dynamic mode decomposition, Koopman Lyapunov-based
model predictive control, offset-free tracking

## Introduction

Koopman and Von Neumann showed that a nonlinear dynamical system can be
represented by an operator-theoretic perspective, which is to describe the
system in terms of functions of states (called observables). They also showed
that the temporal evolution of these observables can be described by an
infinite-dimensional linear operator called Koopman operator (Koopman, 1931;
Koopman and Neumann, 1932).

Recently, based on the availability of abundant time-series data and advances
in numerical techniques, several data-driven methods have been developed to
obtain a finite-dimensional approximation to the infinite-dimensional Koopman
operator for practical application of Koopman operator theory (Rowley, Mezić,
Bagheri, Schlatter, Henningson et al., 2009; Mauroy and Mezić, 2012; Peitz and
Klus, 2019). Dynamic mode decomposition (DMD) proposed in Williams, Kevrekidis
and Rowley (2015, 2016) produces an approximation of the Koopman operator
using a set of monomials. DMD has been applied in many fields by the benefit
of its ease of implementation. However, since DMD only utilizes a limited
space that consists of scalar observables, DMD cannot provide proper
prediction power for highly nonlinear systems. Extended dynamic mode
decomposition (EDMD) is developed in Williams, Kevrekidis and Rowley (2015,
2016) to handle this limitation of DMD. EDMD enriches the observable space by
using nonlinear functions and provides a proper numerical approximation of
Koopman operator.

A key feature of data-driven model identification methods based on Koopman
operator theory is the availability of linear predictors in lifted space for
nonlinear dynamical systems. These linear predictors allow for implementation
of well-established linear controller design methods to nonlinear systems such
as linear model predictive control (MPC) that has been applied in various
applications as a standard model-based control strategy (Christofides,
Scattolini, de la Peña and Liu, 2013; Son, Park, Oh, Kim and Lee, 2020b, a).
Based on this feature, various Koopman-based model predictive control (KMPC)
schemes have been developed and applied for control of several nonlinear
systems (Narasingam and Kwon, 2020a; Arbabi, Korda and Mezić, 2018; Sootla,
Mauroy and Ernst, 2018). Narasingam and Kwon (2019) proposed Koopman Lyapunov-
based model predictive control (KLMPC) that integrates KMPC with Lyapunov-
based MPC (LMPC) scheme presented in Mhaskar, El-Farra and Christofides
(2006); de la Peña and Christofides (2008) to guarantee feasibility and
stability of control systems. Additionally, Narasingam and Kwon (2020b)
presented a mathematical analysis for feasibility and stability of KLMPC
system in the Lyapunov sense.

Since a finite-dimensional approximation obtained from previously introduced
data-driven numerical methods (e.g., EMDM) cannot completely represent the
infinite-dimensional Koopman operator, plant-model mismatch is inevitable in
KLMPC. This inherent plant-model mismatch in KLMPC can significantly degrade
the closed-loop performance of the KLMPC systems, but to the best of our
knowledge, handling of this plant-model mismatch in KLMPC system has not been
rigorously studied yet. To this end, we propose an offset-free KLMPC framework
that integrates KLMPC with the offset-free MPC framework to address the
inherent plant-model mismatch in KLMPC while guaranteeing feasibility and
stability of a control system with Lyapunov constraints.

Offset-free MPC is one of the most representative plant-model mismatch
compensating schemes in the MPC field. There are two ways to achieve offset-
free tracking in the presence of plant-model mismatch. One method exploits the
integration of tracking error in a compensator block as in Johnson (1970);
Davison and Smith (1971). However, since the integrated error is independent
of the controller, this method can cause a windup problem where an overshoot
occurs due to an accumulated error even after the system reaches a desired
point. Second approach is to design a disturbance estimator to compensate for
plant-model mismatch by augmenting a disturbance model to the system model and
estimating disturbance-augmented state from measurement (Muske and Badgwell,
2002; Pannocchia and Rawlings, 2003). Since the disturbance estimator approach
has an anti-windup effect (Maeder and Morari, 2010), this approach is one of
the most popular approaches to accomplish offset-free tracking in MPC (Maeder,
Borrelli and Morari, 2009).

Therefore, in this study, we develop an offset-free KLMPC framework based on
the disturbance estimator approach. Specifically, a disturbance model is
introduced to consider the influence of plant-model mismatch, and is augmented
with a linear model in the lifted state space, identified by EDMD. Then, based
on the disturbance-augmented model, an estimator is designed to estimate the
lifted state and disturbance. A target problem is also designed to obtain
proper target values for the lifted states and inputs that will be applied to
the optimal control problem by considering the estimated disturbance.
Specifically, the optimal control problem is designed based on the
disturbance-augmented model and Lyapunov constraints to obtain an optimal
input to track the set-point while compensating for the influence of plant-
model mismatch on the controlled variables and ensuring the feasibility and
stability of control system. Unlike the nominal KLMPC scheme, since the target
lifted state and input values of the optimal control problem are continuously
updated through the target problem, the Lyapunov function and stabilizing
control law within the Lyapunov constraints are also continuously updated
during the closed-loop operation. A mathematical analysis for zero steady-
state offset condition of the proposed offset-free KLMPC framework is also
presented.

The rest of this paper is organized as follows. The next section introduces a
data-driven EDMD scheme to derive linear models in the lifted state-space for
nonlinear dynamical systems and presents a formulation of the nominal KLMPC
scheme that integrates EDMD and Lyapunov-based control. In Section 2, a
detailed formulation of the offset-free KLMPC framework and mathematical
analysis for zero stead-state offset condition are presented. In Section 3,
the performance of the proposed framework in handling plant-model mismatch
compared to the nominal KLMPC is shown with a closed-loop simulation of a CSTR
system. In the last section, we conclude with a few important remarks.

##  1 Preliminaries

###  1.1 Data-driven linear model identification of nonlinear systems

In this section, we introduce a data-driven model identification method to
obtain a linear model in a higher dimensional lifted space of nonlinear
systems, based on Koopman operator theory.

In Koopman (1931); Koopman and Neumann (1932), Koopman suggested an
alternative operator theoretic perspective to describe the dynamics of an
uncontrolled system in (1) in terms of the evolution of observable functions
$\psi:\mathbb{R}^{n_{x}}\rightarrow\mathbb{R}$.

$\displaystyle x(k+1)=f(x(k))$ (1)

where $x\in\mathbb{R}^{n_{x}}$ represents the state, and $n_{x}$ denotes the
dimension of state.

Koopman has also shown that an infinite-dimensional linear operator
$\mathcal{K}:\mathcal{F}\rightarrow\mathcal{F}$ exists which advances the
observables forward in time as follows:

$\displaystyle\mathcal{K}\psi(x)=\psi(f(x)).$ (2)

where $\mathcal{F}$ is a space of observables invariant under the action of
the Koopman operator (Budišić, Mohr and Mezić, 2012).

Since the Koopman operator theory was initially introduced based on
uncontrolled systems, several schemes have been suggested to generalize the
Koopman operator theory to controlled systems of the following form (Williams,
Hemati, Dawson, Kevrekidis and Rowley, 2016; Proctor, Brunton and Kutz, 2018;
Korda and Mezić, 2018):

$\displaystyle x(k+1)=f(x(k),u(k))$ (3)

where $u\in\mathcal{U}\subset\mathbb{R}^{n_{u}}$ represents the input, and
$n_{u}$ denotes the dimension of input.

In this study, we present a scheme in Korda and Mezić (2018) that generalizes
the Koopman operator theory for controlled systems by introducing an extended
state as

$\displaystyle\chi=\begin{bmatrix}x\\\ \mathbf{u}\end{bmatrix}$ (4)

where $\mathbf{u}:=\\{u_{i}\\}^{\infty}_{i=0}\in\ell(\mathcal{U})$
($u_{i}\in\mathcal{U}$) represents the input sequence and $\ell(\mathcal{U})$
denotes the space of all input sequences $\mathbf{u}$. The dynamics of the
extended state $\chi$ is described as

$\displaystyle f_{\chi}(\chi)=\begin{bmatrix}f(x,\mathbf{u}(0))\\\
\mathcal{S}\mathbf{u}\end{bmatrix}$ (5)

where $\mathbf{u}(i)$ denotes the $i^{th}$ element of $\mathbf{u}$ and
$\mathcal{S}$ represents the left shift operator:

$\displaystyle(\mathcal{S}\mathbf{u})(i)=\mathbf{u}(i+1).$ (6)

Then, the Koopman operator $\mathcal{K}:\mathcal{H}\rightarrow\mathcal{H}$
associated with the dynamics of the extended state in (5) can be defined based
on an extended observable
$\phi:\mathbb{R}^{n_{x}}\times\ell(\mathcal{U})\rightarrow\mathbb{R}$:

$\displaystyle\mathcal{K}\phi(\chi)=\phi(f_{\chi}(\chi))$ (7)

where $\mathcal{H}$ is an extended observable space.

Since the Koopman operator itself is not suitable for practical implementation
due to its infinite dimensional nature, several schemes have been developed
for finite-dimensional approximation to the Koopman operator. In this work, we
present a data-driven EDMD algorithm (Williams, Kevrekidis and Rowley, 2015,
2016).

In EDMD, a vector of observables $\phi$ is designed as

$\displaystyle\phi(x,\mathbf{u})=\begin{bmatrix}\psi(x)\\\
\mathbf{u}(0)\end{bmatrix}$ (8)
$\displaystyle\psi(x):=\begin{bmatrix}\psi_{1}(x),\cdots,\psi_{n_{z}}(x)\end{bmatrix}^{\top}$

where $n_{z}$ is the number of observables depending on $x$.

The scheme assumes that a collection of data
$((x_{j},\mathbf{u}_{j}),(x_{j}^{+},\mathbf{u}_{j}^{+}))$, $j=1,\cdots,N_{d}$
that satisfying (9) is available.

$\displaystyle\begin{bmatrix}x_{j}^{+}\\\
\mathbf{u}_{j}^{+}\end{bmatrix}=\begin{bmatrix}f(x_{j},\mathbf{u}_{j}(0))\\\
\mathcal{S}\mathbf{u}_{j}\end{bmatrix}$ (9)

where the superscript $+$ denotes the value at the next time step. Then, based
on the data, an approximation of the Koopman operator $\mathcal{A}$ is
obtained by minimizing

$\displaystyle
J(\mathcal{A})=\sum_{j=1}^{N_{d}}||\phi(x_{j}^{+},\mathbf{u}_{j}^{+})-\mathcal{A}\phi(x_{j},\mathbf{u}_{j})||^{2}$
(10)

where $||\cdot||$ denotes the Euclidean norm.

Since it is not essential to predict the future control sequence, we can
discard the last $n_{u}$ components of each
$\phi(x_{j}^{+},\mathbf{u}_{j}^{+})$ and the last $n_{u}$ rows of
$\mathcal{A}$. Additionally, let $\bar{\mathcal{A}}$ denote the remaining part
of $\mathcal{A}$ after discarding the part associated with the future input,
then we can decompose $\bar{\mathcal{A}}$ into $A\in\mathbb{R}^{n_{z}\times
n_{z}}$ and $B\in\mathbb{R}^{n_{z}\times n_{u}}$ as

$\displaystyle\bar{\mathcal{A}}=\begin{bmatrix}A&B\end{bmatrix}.$ (11)

Then, substituting (11) into the remaining part of (10) yields

$\displaystyle\bar{J}(A,B)=\sum_{j=1}^{N_{d}}||\psi(x_{j}^{+})-A\psi(x_{j})-B\textbf{u}_{j}(0)||^{2}.$
(12)

Finally, we can obtain a linear predictor in lifted space as in (13) by
minimizing (12) over $A$ and $B$.

$\displaystyle z(k+1)=Az(k)+Bu(k)$ (13)

where $z\in\mathbb{R}^{n_{z}}$ denotes the lifted state:

$\displaystyle z=\psi(x).$ (14)

Additionally, let the output be given as a function of the state as

$\displaystyle y(k)=g(x(k)).$ (15)

In this case, the output matrix $C$ is obtained as the best projection of $y$
onto the span of observable functions in a least square sense by minimizing

$\displaystyle J_{y}(C)=\sum_{j=1}^{N_{d}}||y_{j}-C\psi(x_{j})||^{2}.$ (16)

Then, a linear model for $y$ can be derived using $C$ as follows:

$\displaystyle y(k)=Cz(k).$ (17)

###  1.2 KLMPC design

Based on the data-driven model identification methodology presented in the
previous section, Narasingam and Kwon (2019) developed KLMPC framework that
integrates LMPC scheme in Mhaskar, El-Farra and Christofides (2006); de la
Peña and Christofides (2008) with EDMD to ensure feasibility and stability of
control systems.

In LMPC, following Lyapunov constraints are applied to the optimal control
problem:

$\displaystyle V(x_{i+1}-\bar{x}_{s})\leq r,\quad i=0,\dots,N-1$ (18)
$\displaystyle V(x_{1}-\bar{x}_{s})\leq V(x_{1}^{h}-\bar{x}_{s})$ (19)

where $V(x-\bar{x}_{s})$ represents the Lyapunov function, $\bar{x}_{s}$
denotes the target steady-state point of the control system,
$x_{1},\cdots,x_{N-1}$ denote future states, and $x_{1}^{h}$ denotes the state
evolved from $x_{0}$ with a stabilizing control law. The first Lyapunov
constraint in (18) guarantees that future states stay within the stability
region which can be described as a sublevel set of the Lyapunov function
$\Omega_{r}:=\\{x\in{\mathbb{R}^{n_{x}}}:V(x-\bar{x}_{s})\leq r\\}$. Then, the
second Lyapunov constraint in (19) ensures that the Lyapunov function value at
the next time step induced by the MPC control law is smaller than or equal to
that induced by a stabilizing control law. These Lyapunov constraints ensure
the closed-loop stability of LMPC systems.

However, since these Lyapunov constraints are nonlinear, the optimal control
problem becomes a non-convex problem. Narasingam and Kwon (2019, 2020b)
effectively address this issue by transforming the nonlinear Lyapunov
constraints in (18) and (19) to linear constraints with EDMD. Specifically,
the Lyapunov function is incorporated into the observable function library and
described as one of the lifted states.

For example, when the Lyapunov function is a quadratic form, it can be
incorporated into the observable function library as the $j_{v}^{th}$
observable:

$\displaystyle\psi_{j_{v}}(x)=(x-\bar{x}_{s})^{\top}Q_{v}(x-\bar{x}_{s}).$
(20)

Then, after obtaining the linear model in the lifted space in (13) through
EDMD, the Lyapunov function value can be simply derived from the lifted state
as follows:

$\displaystyle V(x-\bar{x}_{s})=D_{v}z$ (21)

where
$D_{v}:=\begin{bmatrix}\mathbf{0}_{1,j_{v}-1},1,\mathbf{0}_{1,n_{z}-j_{v}}\end{bmatrix}$.
The Lyapunov constraints can be reformulated as (22) and (23) by substituting
(21) into (18) and (19).

$\displaystyle D_{v}z_{i+1}\leq r$ (22) $\displaystyle D_{v}z_{1}\leq
V(D_{x}(z_{1}^{h}-\bar{z}_{s}))$ (23)

where

$\displaystyle
D_{x}:=\begin{bmatrix}D_{{x_{1}}}^{\top},\cdots,D_{{x_{n_{x}}}}^{\top}\end{bmatrix}^{\top}$
$\displaystyle
D_{x_{j}}:=\begin{bmatrix}\mathbf{0}_{1,j-1},1,\mathbf{0}_{1,n_{z}-j}\end{bmatrix}\quad\mathrm{for}\quad
j=1,\cdots,n_{x}$ $\displaystyle\bar{z}_{s}:=\psi(\bar{x}_{s}).$

Then, the finite-horizon optimal control problem of KLMPC can be constructed
as a standard convex quadratic problem in (24) by applying the Lyapunov
constraints in (22) and (23):

$\displaystyle\underset{\small{u_{0},\cdots,u_{N-1}}}{\mathrm{min}}$
$\displaystyle\;\sum_{i=0}^{N-1}{||z_{i+1}-\bar{z}_{s}||^{2}_{Q_{z}}+||u_{i}-\bar{u}_{s}||^{2}_{Q_{u}}}$
(24a) $\displaystyle\mathrm{s.t.}\quad$ $\displaystyle z_{0}=z$ (24b)
$\displaystyle z_{i+1}=Az_{i}+Bu_{i}$ (24c) $\displaystyle
z_{1}^{h}=Az_{0}+Bh_{0}(z_{0},\bar{z}_{s})$ (24d) $\displaystyle u_{min}\leq
u_{i}\leq u_{max}$ (24e) $\displaystyle y_{min}\leq Cz_{i+1}\leq y_{max}$
(24f) $\displaystyle D_{v}z_{i+1}\leq r,\quad i=0,\dots,N-1$ (24g)
$\displaystyle D_{v}z_{1}\leq V(D_{x}(z_{1}^{h}-\bar{z}_{s}))$ (24h)

where $\bar{u}_{s}$ is the steady-state target input.

##  2 Offset-free Koopman Lyapunov-based MPC

The KLMPC scheme in the previous section provides controller feasibility and
stability by effectively integrating EDMD and LMPC. However, the closed-loop
performance of a KLMPC system can be degraded due to the plant-model mismatch
between a real plant dynamics in (3) and a data-driven linear model, obtained
by EDMD, in (13). To this end, we present an offset-free KLMPC framework that
integrates an offset-free MPC framework in Maeder and Morari (2010); Borrelli,
Bemporad and Morari (2017) with the KLMPC scheme to compensate for this plant-
model mismatch.

###  2.1 Offset-free KLMPC framework

We suppose a data-driven linear model for lifted state $z$ in (13) and (17)
has been identified by EDMD in Section 1.1 and the controlled variable is
described as

$\displaystyle y_{c}(k)=Hy(k)$ (25)

where $y_{c}\in\mathbb{R}^{n_{y_{c}}}$ is the controlled variable and
$n_{y_{c}}$ denote the dimension of the controlled variable.

In the disturbance estimator approach, a disturbance model is augmented with
the lifted state-space model to consider the effect of plant-model mismatch.
Though various disturbance-augmented models have been proposed in Betti,
Farina and Scattolini (2013); Pannocchia, Gabiccini and Artoni (2015); Wang
(2009), in this study, we focus on the following form (Pannocchia and
Rawlings, 2003):

$\displaystyle\begin{cases}z(k+1)=Az(k)+Bu(k)+B_{d}d(k)\\\ d(k+1)=d(k)\\\
y(k)=Cz(k)+C_{d}d(k)\end{cases}$ (26)

where $d\in\mathbb{R}^{n_{d}}$ is the disturbance, and
$B_{d}\in\mathbb{R}^{n_{z}\times n_{d}}$ and $C_{d}\in\mathbb{R}^{n_{y}\times
n_{d}}$ denote matrices that represent the influence of the disturbance
variable on the evolution of the lifted state and the output, respectively.

$B_{d}$ and $C_{d}$ can be freely determined under the restriction that the
disturbance-augmented system in (26) is observable. The observability
condition of the augmented system associated with $B_{d}$ and $C_{d}$ design
is given in the following proposition.

###### Proposition 1

The disturbance-augmented system in (26) is observable if and only if the
original linear dynamical system in (13) is observable and the following
condition holds. The detailed derivation of the condition is provided in
Proposition 1 of Maeder, Borrelli and Morari (2009).

$\displaystyle{\mathrm{rank}}\begin{bmatrix}I-A&-B_{d}\\\
C&C_{d}\end{bmatrix}=n_{z}+n_{d}.$ (27)

###### Remark 1

For the condition in (27) to be satisfied, the dimension of the disturbance
should be smaller than or equal to the number of available measurements
(Rawlings, Mayne and Diehl, 2017; Morari and Maeder, 2012):

$\displaystyle n_{d}\leq n_{y}.$ (28)

The disturbance-augmented model with properly designed $B_{d}$ and $C_{d}$
allows a model predictive controller to effectively handle the negative impact
of plant-model mismatch on the lifted state and the output. The following
lifted state and disturbance estimator is designed based on the constructed
model in (26).

$\displaystyle\begin{bmatrix}\hat{z}(k+1)\\\ \hat{d}(k+1)\end{bmatrix}$
$\displaystyle=\begin{bmatrix}A&B_{d}\\\
0&I\end{bmatrix}\begin{bmatrix}\hat{z}(k)\\\
\hat{d}(k)\end{bmatrix}+\begin{bmatrix}B\\\ 0\end{bmatrix}u(k)$
$\displaystyle+\begin{bmatrix}L_{z}\\\
L_{d}\end{bmatrix}(-y_{p}(k)+C\hat{z}(k)+C_{d}\hat{d}(k))$ (29)

where $\hat{z}$ and $\hat{d}$ represent the estimated lifted state and
disturbance, respectively, $L_{z}$ and $L_{d}$ are the estimator gains for
lifted state and disturbance, respectively, that make the estimator stable,
and $y_{p}$ denotes the measured output from the plant.

Then, the equilibrium target lifted state $\bar{z}$ and target input $\bar{u}$
that satisfy the following equations are obtained:

$\displaystyle\bar{z}=A\bar{z}+B\bar{u}+B_{d}\hat{d}$ (30)
$\displaystyle\bar{y}_{c}=HC\bar{z}+HC_{d}\hat{d}.$ (31)

We can construct the following equation from (30) and (31):

$\displaystyle\begin{bmatrix}A-I&B\\\
HC&0\end{bmatrix}\begin{bmatrix}\bar{z}\\\
\bar{u}\end{bmatrix}=\begin{bmatrix}-B_{d}\hat{d}\\\
\bar{y}_{c}-HC_{d}\end{bmatrix}.$ (32)

A solution of (32) exists with proper $H$ and $\bar{y}_{c}$, but it may not be
unique. Therefore, a set of target lifted state and target input
$(\bar{z},\bar{u})$ is obtained by solving the following target problem based
on desired target values $(\bar{z}_{s},\bar{u}_{s})$ and process constraints:

$\displaystyle\underset{\small{\bar{z},\bar{u}}}{\mathrm{min}}\quad$
$\displaystyle||\bar{u}-\bar{u}_{s}||^{2}_{Q_{\bar{u}}}+||\bar{z}-\bar{z}_{s}||^{2}_{Q_{\bar{z}}}$
(33a) $\displaystyle\mathrm{s.t.}\quad$ $\displaystyle\begin{bmatrix}A-I&B\\\
HC&0\end{bmatrix}\begin{bmatrix}\bar{z}\\\
\bar{u}\end{bmatrix}=\begin{bmatrix}-B_{d}\hat{d}\\\
\bar{y}_{c}-HC_{d}\hat{d}\end{bmatrix}$ (33b) $\displaystyle
u_{min}\leq\bar{u}\leq u_{max}$ (33c) $\displaystyle y_{min}\leq\bar{y}\leq
y_{max}.$ (33d)

where $\bar{y}=C\bar{z}+C_{d}\hat{d}$.

###### Remark 2

In the nominal offset-free MPC framework, the objective function of the target
problem in (33a) is commonly designed to drive the target output value
$\bar{y}$ to a desired value $\bar{y}_{s}$ (Muske and Badgwell, 2002;
Pannocchia and Bemporad, 2007; Rawlings, Mayne and Diehl, 2017). For example,
the following objective function is commonly utilized in nominal offset-free
MPC.

$\displaystyle\bar{J}_{0}=||\bar{u}-\bar{u}_{s}||^{2}_{Q_{\bar{u}}}+||\bar{y}-\bar{y}_{s}||^{2}_{Q_{\bar{y}}}$
(34)

This objective function induces the state target $\bar{x}$ to stay near a
desired state $\bar{x}_{s}$ that satisfies

$\displaystyle C_{x}\bar{x}_{s}+C_{d}\hat{d}=\bar{y}_{s}.$ (35)

However, in the case of the proposed offset-free KLMPC framework, since the
system for the lifted state in (26) from EDMD is higher dimensional than the
original system in (3), the influence of the output correction on the lifted
state target $\bar{z}$ is insignificant. Therefore, we replaced the output
correction term $||\bar{y}-\bar{y}_{s}||^{2}_{Q_{\bar{u}}}$ in the objective
function in (34) with the lifted state correction term
$||\bar{z}-\bar{z}_{s}||^{2}_{Q_{\bar{z}}}$ as in (33a) to directly drive the
target lifted state to a desired point. The derivation of a proper target
lifted state $\bar{z}$ with this setting is crucial because the derived target
value will be applied to a finite-horizon optimal control problem.

As described in Section 1.2, the Lyapunov constraints in (18) and (19) are
easily handled in the KLMPC scheme by including the Lyapunov function with
fixed target $\bar{x}_{s}$ as a component of the observable function library.
However, in the proposed offset-free KLMPC framework, since the target lifted
state $\bar{z}$ of the optimal control problem is continuously updated over
time through the target problem in (33) as the estimated disturbance $\hat{d}$
is continuously updated from the estimator in (2.1), the following Lyapunov
constraints also should be continuously updated.

$\displaystyle V(x_{i+1}-\bar{x})\leq r,\quad i=0,\dots,N-1$ (36)
$\displaystyle V(x_{1}-\bar{x})\leq V(x_{1}^{h}-\bar{x})$ (37)
$\displaystyle\bar{x}:=D_{x}\bar{z}$ $\displaystyle x_{1}^{h}:=D_{x}z_{1}^{h}$

where $z_{1}^{h}$ is the lifted state evolved from $z_{0}$ by the stabilizing
control law $h(z_{0},\hat{d},\bar{z})$:

$\displaystyle z_{1}^{h}=Az_{0}+Bh(z_{0},\hat{d},\bar{z})+B_{d}\hat{d}.$ (38)

Since the equilibrium point $\bar{x}$ in the Lyapunov function $V(x-\bar{x})$
is continuously updated, we cannot simply replace $V(x-\bar{x})$ with the
observable $\psi_{j_{v}}=(x-\bar{x}_{s})^{\top}Q_{v}(x-\bar{x}_{s})$ in (20)
which is the Lyapunov function with the fixed equilibrium point $\bar{x}_{s}$
as in the nominal KLMPC. However, we can describe $V(x-\bar{x})$ using
$\psi_{j_{v}}$ with other terms associated with the update of $\bar{x}$. More
specifically, suppose the Lyapynov function is in the following quadratic
form:

$\displaystyle V(x-\bar{x})=(x-\bar{x})^{\top}Q_{v}(x-\bar{x}).$ (39)

Then, $V(x-\bar{x})$ can be reformulated with $\psi_{j_{v}}$ as follows:

$\displaystyle
V(x-\bar{x})=\psi_{j_{v}}+2(\bar{x}_{s}-\bar{x})^{\top}Q_{v}x+\bar{x}^{\top}Q_{v}\bar{x}-\bar{x}_{s}^{\top}Q_{v}\bar{x}_{s}.$
(40)

Substituting $x=D_{x}z$, $\bar{x}_{s}=D_{x}\bar{z}_{s}$, and
$\psi_{j_{v}}=D_{v}z$ into (40) yields

$\displaystyle
V(x-\bar{x})=F_{v}z+\bar{z}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}-\bar{z}_{s}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}_{s}$
(41) $\displaystyle
F_{v}:=D_{v}+2(\bar{z}_{s}-\bar{z})^{\top}D_{x}^{\top}Q_{v}D_{x}.$

With (41), we can effectively update the Lyapunov function as the target
lifted state $\bar{z}$ is updated. Additionally, by applying (41) to (36) and
(37), the Lyapunov constraints can be described as linear inequalities as
follows:

$\displaystyle F_{v}z_{i+1}\leq
r-\bar{z}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}+\bar{z}_{s}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}_{s}$
(42) $\displaystyle F_{v}z_{1}\leq
V(D_{x}(z_{1}^{h}-\bar{z}))-\bar{z}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}+\bar{z}_{s}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}_{s}.$
(43)

Not only the Lyapunov function, but also the stabilizing control law
$h(z,\hat{d},\bar{z})$ in (38) should be continuously updated over time as the
target lifted state $\bar{z}$ and estimated disturbance $\hat{d}$ are updated.
Therefore, $h(z,\hat{d},\bar{z})$ is updated based on the following state
feedback law (Preuss, 1982) as follows:

$\displaystyle h(z,d,\bar{z})=\bar{N}\bar{z}-K_{z}z-K_{d}d$ (44)

where $\bar{N}$ and $K_{d}$ are gains associated with $\bar{z}$ and $d$,
respectively, and $K_{z}$ is the feedback gain for the lifted state that makes
the controlled system stable.

###### Proposition 2

The stabilizing control law in (44) yields perfect set-point tracking at
steady state with the following $\bar{N}$ and $K_{d}$.

$\displaystyle\bar{N}=[(I-A+BK_{z})^{-1}B]^{\dagger}$ (45) $\displaystyle
K_{d}=\bar{N}(I-A+BK_{z})^{-1}B_{d}$ (46)

where $\dagger$ denotes the pseudo inverse.

###### Proof 1

Substituting $h(z,d,\bar{z})$ in (44) into the linear dynamical system in (26)
yields the following controlled system:

$\displaystyle z(k+1)=(A-BK_{z})z(k)+B\bar{N}\bar{z}(k)+(B_{d}-BK_{d})d(k).$
(47)

Let $z_{\infty}^{h}$, $\bar{z}_{\infty}$, and $d_{\infty}$ represent the
steady-state lifted state, target lifted state, and disturbance of the
controlled system in (47). At steady state, applying $z_{\infty}^{h}$,
$\bar{z}_{\infty}$, and $d_{\infty}$ to (47) and rearranging yields

$\displaystyle z_{\infty}^{h}=($ $\displaystyle
I-A+BK_{z})^{-1}B\bar{N}\bar{z}_{\infty}$ (48)
$\displaystyle+(I-A+BK_{z})^{-1}(B_{d}-BK_{d})d_{\infty}.$

Then, substituting $\bar{N}$ and $K_{d}$ in (45) and (46) into (48) yields

$\displaystyle z_{\infty}^{h}=\bar{z}_{\infty}.$ (49)

Therefore, we can see the controlled system with the stabilizing control law
yields perfect set-point tracking at steady state.

Given the targets ($\bar{z}$, $\bar{u}$) from the target problem in (33) and
estimates ($\hat{z}$, $\hat{d}$) from the estimator in (2.1), a finite-horizon
optimal control problem for the offset-free KLMPC can be brought into a lifted
state-regulation problem $P(\hat{z},\hat{d},\bar{z},\bar{u})$ as follows:

$\displaystyle\underset{\small{u_{0},\cdots,u_{N-1}}}{\mathrm{min}}$
$\displaystyle\;\sum_{i=0}^{N-1}{||z_{i+1}-\bar{z}||^{2}_{Q_{z}}+||u_{i}-\bar{u}||^{2}_{Q_{u}}}$
(50a) $\displaystyle\mathrm{s.t.}\quad$ $\displaystyle z_{0}=\hat{z}$ (50b)
$\displaystyle z_{i+1}=Az_{i}+Bu_{i}+B_{d}\hat{d}$ (50c) $\displaystyle
z_{1}^{h}=Az_{0}+Bh(z_{0},\hat{d},\bar{z})+B_{d}\hat{d}$ (50d) $\displaystyle
u_{min}\leq u_{i}\leq u_{max}$ (50e) $\displaystyle y_{min}\leq
Cz_{i+1}+C_{d}\hat{d}\leq y_{max}$ (50f) $\displaystyle F_{v}z_{i+1}\leq
r-\bar{z}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}$ (50g)
$\displaystyle\qquad\quad+\bar{z}_{s}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}_{s},\quad
i=0,\dots,N-1$ $\displaystyle F_{v}z_{1}\leq
V(D_{x}(z_{1}^{h}-\bar{z}))-\bar{z}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}$ (50h)
$\displaystyle\qquad\quad+\bar{z}_{s}^{\top}D_{x}^{\top}Q_{v}D_{x}\bar{z}_{s}.$

The Lyapunov constraints in (50g) and (50h) are formulated based on (42),
(43), and (44). By solving the optimal control problem, we obtain an optimal
input to move the original system in (3) to steady state where the zero
steady-state offset is satisfied in the presence of plant-model mismatch.

###  2.2 Zero steady-state offset condition

In this section, we derive the condition for zero steady-state offset of the
proposed offset-free KLMPC by examining the closed-loop behavior of the
control system at steady state ($\hat{z}(k)\rightarrow\hat{z}_{\infty}$,
$\hat{d}(k)\rightarrow\hat{d}_{\infty}$,
$\bar{z}(k)\rightarrow\bar{z}_{\infty}$,
$\bar{u}(k)\rightarrow\bar{u}_{\infty}$, $u(k)\rightarrow u_{\infty}$, and
$y_{p}(k)\rightarrow y_{p,\infty}$ as $k\rightarrow\infty$).

By considering the estimator in (2.1) at steady state, (51) and (52) are
obtained.

$\displaystyle\hat{z}_{\infty}=A\hat{z}_{\infty}+Bu_{\infty}+B_{d}\hat{d}_{\infty}-L_{z}e_{y,\infty}$
(51) $\displaystyle 0=L_{d}e_{y,\infty}$ (52)

where $u_{\infty}$ is an optimal control law derived from the optimal control
problem in (50), and $e_{y,\infty}$ is the error between predicted output and
measured output:

$\displaystyle e_{y,\infty}=y_{p,\infty}-\hat{y}_{\infty}$ (53)

where $\hat{y}_{\infty}:=C\hat{z}_{\infty}+C_{d}\hat{d}_{\infty}$. Similarly,
(54) and (55) are derived by considering (30) and (31) at steady state as
follows:

$\displaystyle\bar{z}_{\infty}=A\bar{z}_{\infty}+B\bar{u}_{\infty}+B_{d}\hat{d}_{\infty}$
(54) $\displaystyle\bar{y}_{c}=HC\bar{z}_{\infty}+HC_{d}\hat{d}_{\infty}.$
(55)

Now, let
$P(\hat{z}_{\infty},\hat{d}_{\infty},\bar{z}_{\infty},\bar{u}_{\infty})$
denote the optimal control problem in (50) at steady state with the converged
lifted state estimate $\hat{z}_{\infty}$, disturbance estimate
$\hat{d}_{\infty}$, lifted state target $\bar{z}_{\infty}$, and input target
$\bar{u}_{\infty}$. Then, we derive the relation between the optimal control
input $u_{\infty}$ from
$P(\hat{z}_{\infty},\hat{d}_{\infty},\bar{z}_{\infty},\bar{u}_{\infty})$ and
the lifted state estimate $\hat{z}_{\infty}$ in the following lemma.

###### Lemma 3

Let $P_{\hat{z}_{\infty}}$ represent the optimal control problem
$P(\hat{z}_{\infty},\hat{d}_{\infty},\bar{z}_{\infty},\bar{u}_{\infty})$.
Without loss of generality, we assume that the set-point is far from the upper
and lower bounds of process constraints, and thus, the process constraints in
(50e) and (50f) and the first Lyapunov constraint in (50g) are inactive at the
optimum of $P_{\hat{z}_{\infty}}$. Then, the optimal control law of
$P(\hat{z}_{\infty},\hat{d}_{\infty},\bar{z}_{\infty},\bar{u}_{\infty})$ can
be described as

$\displaystyle u_{\infty}=K_{mpc}\hat{z}_{\infty}+c_{mpc}$ (56)

with

$\displaystyle(K_{mpc},c_{mpc})=\begin{cases}(K_{un},c_{un})&\mathrm{when\;(\ref{eq50h})\;is\;inactive}\\\
(K_{\ell},c_{\ell})&\mathrm{when\;(\ref{eq50h})\;is\;active}\end{cases}$ (57)

where $(K_{un},c_{un})$ denotes a pair of the lifted state gain and the
constant term of the unconstrained optimal solution of $P_{\hat{z}_{\infty}}$,
$(K_{\ell},c_{\ell})$ denotes a pair of the lifted state gain and the constant
term of the optimal solution of $P_{\hat{z}_{\infty}}$ when the Lyapunov
constraint in (50h) is active at the optimum.

###### Proof 2

Since the optimal control problem in (50) is a quadratic program, the
optimizer function is described as a piecewise affine function of the lifted
state $\hat{z}_{\infty}$ (Bemporad, Morari, Dua and Pistikopoulos, 2002;
Alessio and Bemporad, 2009). Based on this fact, we derive the optimal
solution for each case as a function of $\hat{z}_{\infty}$.

Case 1: We can reformulate the objective function in (50a) of
$P_{\hat{z}_{\infty}}$ as a quadratic function of
$U=[u_{0}^{\top},\cdots,u_{N-1}^{\top}]^{\top}$:

$\displaystyle
J_{\hat{z}_{\infty}}(U)=\frac{1}{2}U^{\top}H_{\infty}U+f_{\infty}^{\top}U$
(58)

When all the constraints are inactive at the optimum of
$P_{\hat{z}_{\infty}}$, the following unconstrained optimal solution
$U^{*}_{un}$ derived from the first-order optimality condition of
$J_{\hat{z}_{\infty}}(U)$ is feasible and optimal:

$\displaystyle U^{*}_{un}=-H_{\infty}^{-1}f_{\infty}$ (59)

where

$\displaystyle H_{\infty}:=\Psi^{\top}Q_{\textbf{z}}\Psi+Q_{\textbf{u}}$
$\displaystyle
f_{\infty}:=\Psi^{\top}Q_{\textbf{z}}(\Phi\hat{z}_{\infty}+\Psi_{d}\hat{\mathbf{d}}_{\infty}-\bar{\mathbf{z}}_{\infty})-Q_{\textbf{u}}\bar{\mathbf{u}}_{\infty}$
$\displaystyle\Psi:=\begin{bmatrix}B&0&\cdots&0\\\ AB&B&\cdots&0\\\
\vdots&\vdots&\ddots&\vdots\\\
A^{\textit{N}-1}B&A^{\textit{N}-2}B&\cdots&B\end{bmatrix}$
$\displaystyle\Psi_{d}:=\begin{bmatrix}B_{d}&0&\cdots&0\\\
AB_{d}&B_{d}&\cdots&0\\\ \vdots&\vdots&\ddots&\vdots\\\
A^{\textit{N}-1}B_{d}&A^{\textit{N}-2}B_{d}&\cdots&B_{d}\end{bmatrix}$
$\displaystyle\Phi:=\begin{bmatrix}A\\\ A^{2}\\\ \vdots\\\
A^{N}\end{bmatrix},\;\hat{\mathbf{d}}_{\infty}:=\begin{bmatrix}\hat{d}_{\infty}\\\
\hat{d}_{\infty}\\\ \vdots\\\
\hat{d}_{\infty}\end{bmatrix},\;\bar{\mathbf{z}}_{\infty}:=\begin{bmatrix}\bar{z}_{\infty}\\\
\bar{z}_{\infty}\\\ \vdots\\\
\bar{z}_{\infty}\end{bmatrix},\;\bar{\mathbf{u}}_{\infty}:=\begin{bmatrix}\bar{u}_{\infty}\\\
\bar{u}_{\infty}\\\ \vdots\\\ \bar{u}_{\infty}\end{bmatrix}$ $\displaystyle
Q_{\textbf{z}}:=diag\\{Q_{z},Q_{z},\cdots,Q_{z}\\}$ $\displaystyle
Q_{\textbf{u}}:=diag\\{Q_{u},Q_{u},\cdots,Q_{u}\\}$

where $diag\\{Q_{1},Q_{2},\cdots,Q_{n}\\}$ represents a block diagonal matrix
that has $Q_{1},Q_{2},\cdots,Q_{n}$ as the main diagonal blocks. Then, the
control law $u_{\infty}$ can be described as the first component of
$U^{*}_{un}$:

$\displaystyle u_{\infty}=C_{1}U^{*}_{un}$ (60) $\displaystyle
C_{1}:=\begin{bmatrix}I_{n_{u}}&0&\cdots&0\end{bmatrix}.$

By substituting (59) into (60) and rearranging, we obtain

$\displaystyle u_{\infty}=K_{un}\hat{z}_{\infty}+c_{un}$ (61)

where

$\displaystyle
K_{un}:=-C_{1}(\Psi^{\top}Q_{\textbf{z}}\Psi+Q_{\textbf{u}})^{-1}\Psi^{\top}Q_{\textbf{z}}\Phi$
$\displaystyle
c_{un}:=C_{1}(\Psi^{\top}Q_{\textbf{z}}\Psi+Q_{\textbf{u}})^{-1}(\Psi^{\top}Q_{\textbf{z}}(\bar{\textbf{z}}_{\infty}-\Psi_{d}\hat{\textbf{d}}_{\infty})+Q_{\textbf{u}}\bar{\mathbf{u}}_{\infty})$

Case 2: When the Lyapunov constraint in (50h) is active at the optimum of
$P_{\hat{z}_{\infty}}$, the control law $u_{\infty}$ can be described as

$\displaystyle u_{\infty}=K_{\ell}\hat{z}_{\infty}+c_{\ell}$ (62)

where $K_{\ell}$ and $c_{\ell}$ are derived in Appendix A.

From (61) and (62), we can see that (56) holds.

In Lemma 3, we derived the relation between $\hat{z}_{\infty}$ and
$\hat{u}_{\infty}$ by analyzing $P_{\hat{z}_{\infty}}$. Now, we derive the
relation between the target lifted state $\bar{z}_{\infty}$ and target input
$\bar{u}_{\infty}$ at steady state in the following lemma by analyzing
$P(\bar{z}_{\infty},\hat{d}_{\infty},\bar{z}_{\infty},\bar{u}_{\infty})$ where
the estimated lifted state $\hat{z}_{\infty}$, which is set as the initial
lifted state $z_{0}$ in $P_{\hat{z}_{\infty}}$, is replaced with the target
lifted state $\bar{z}_{\infty}$.

###### Lemma 4

Let $P_{\bar{z}_{\infty}}$ represent the optimal control problem
$P(\bar{z}_{\infty},\hat{d}_{\infty},\bar{z}_{\infty},\bar{u}_{\infty})$. We
assume that the process constraints in (50e) and (50f) and the first Lyapunov
constraint in (50g) are inactive at the optimum of $P_{\bar{z}_{\infty}}$.
Then, $\bar{u}_{\infty}$ becomes the optimal control law of
$P_{\bar{z}_{\infty}}$, and the following equation holds:

$\displaystyle\bar{u}_{\infty}=K_{mpc}\bar{z}_{\infty}+c_{mpc}.$ (63)

###### Proof 3

Since the initial lifted state of $P_{\bar{z}_{\infty}}$ is set as
$z_{0}=\bar{z}_{\infty}$,
$\bar{\textbf{u}}_{\infty}:=[\bar{u}_{\infty}^{\top},\cdots,\bar{u}_{\infty}^{\top}]^{\top}$
makes the future lifted states remain at $\bar{z}_{\infty}$ (i.e.,
$z_{1},\cdots,z_{N}=\bar{z}_{\infty}$) from (54). Then, since the target
lifted state and target input of $P_{\bar{z}_{\infty}}$ are set as
$\bar{z}_{\infty}$ and $\bar{u}_{\infty}$, respectively, eventually the value
of the objective function in (50a) of $P_{\bar{z}_{\infty}}$ becomes 0 with
$\bar{\textbf{u}}_{\infty}$. Therefore, $\bar{u}_{\infty}$ is identical to the
unconstrained optimal control law of $P_{\bar{z}_{\infty}}$ that yields the
minimum objective value in the absence of the constraints.

Since all the parameters in $P_{\hat{z}_{\infty}}$ and $P_{\bar{z}_{\infty}}$
are identical except the initial lifted state values ($z_{0}=\hat{z}_{\infty}$
in $P_{\hat{z}_{\infty}}$ and $z_{0}=\bar{z}_{\infty}$ in
$P_{\bar{z}_{\infty}}$, respectively), the unconstrained control law of
$P_{\bar{z}_{\infty}}$ can be described with $\bar{z}_{\infty}$ using the same
$K_{un}$ and $c_{un}$ of $P_{\hat{z}_{\infty}}$ in (61):

$\displaystyle\bar{u}_{\infty}=K_{un}\bar{z}_{\infty}+c_{un}.$ (64)

Now, we analyze the behavior of the Lyapunov constraint in (50h) in
$P_{\bar{z}_{\infty}}$ when $U=\bar{\textbf{u}}_{\infty}$. For
$P_{\bar{z}_{\infty}}$, the stabilizing control law in (44) yields

$\displaystyle
h(\bar{z}_{\infty},\hat{d}_{\infty},\bar{z}_{\infty})=(\bar{N}-K_{z})\bar{z}_{\infty}-K_{d}\hat{d}_{\infty}.$
(65)

Substituting $z_{0}=\bar{z}_{\infty}$ and (65) into (50d) yields

$\displaystyle
z_{1}^{h}=(A-BK_{z}+B\bar{N})\bar{z}_{\infty}+(B_{d}-BK_{d})\hat{d}_{\infty}.$
(66)

We can reformulate (66) as

$\displaystyle
z_{1}^{h}-\bar{z}_{\infty}+(I-A+BK_{z})\bar{z}_{\infty}=B\bar{N}\bar{z}_{\infty}+(B_{d}-BK_{d})\hat{d}_{\infty}.$
(67)

Then, multiplying $(I-A+BK_{z})^{-1}$ on both the sides of (67) and
substituting $\bar{N}$ and $K_{d}$ in (45) and (46) into (67) yields

$\displaystyle(I-A+BK_{z})^{-1}(z_{1}^{h}-\bar{z}_{\infty})=0.$ (68)

Since $I-A+BK_{z}$ is nonsingular, (68) holds when (69) is true.

$\displaystyle z_{1}^{h}=\bar{z}_{\infty}.$ (69)

Since the lifted state evolved from $\bar{z}_{\infty}$ with $\bar{u}_{\infty}$
under the influence of $\hat{d}_{\infty}$ becomes $\bar{z}_{\infty}$ as in
(54), (70) holds when $z_{1}$ is derived from (50c) with
$u_{0}=\bar{u}_{\infty}$:

$\displaystyle z_{1}=\bar{z}_{\infty}.$ (70)

Now, by substituting $z_{1}^{h}$ and $z_{1}$ in (69) and (70) into (50h), both
sides of (50h) become 0 and the equality is satisfied. This means that the
Lyapunov constraint in (50h) is active when $U=\bar{\textbf{u}}_{\infty}$.
Therefore, $\bar{\textbf{u}}_{\infty}$ is the unconstrained optimal solution
of $P_{\bar{z}_{\infty}}$ as shown in (64) as well as the constrained optimal
solution of $P_{\bar{z}_{\infty}}$ when the Lyapunov constraint in (50h) is
active. Since $P_{\bar{z}_{\infty}}$ has the same formulation as
$P_{\hat{z}_{\infty}}$ with a different initial lifted state
$\bar{z}_{\infty}$ as described above, we can describe $\bar{u}_{\infty}$ with
the same $K_{\ell}$ and $c_{\ell}$ in (62):

$\displaystyle\bar{u}_{\infty}=K_{\ell}\bar{z}_{\infty}+c_{\ell}.$ (71)

From (64) and (71), we can see $\bar{u}_{\infty}$ is the optimal control law
of $P_{\bar{z}_{\infty}}$ in both cases $\mathcal{I}_{A}=\emptyset$ or
$\mathcal{I}_{A}=\\{i_{\ell}\\}$, and thus, (63) holds.

Based on Lemma 3 and Lemma 4, we can derive the relation between
$u_{\infty}-\bar{u}_{\infty}$ and $\hat{z}_{\infty}-\bar{z}_{\infty}$ by
subtracting (63) from (56):

$\displaystyle
u_{\infty}-\bar{u}_{\infty}=K_{mpc}(\hat{z}_{\infty}-\bar{z}_{\infty}).$ (72)

Now, the condition and proof for zero-offset at steady state of the proposed
offset-free KLMPC system are given in the following theorem.

###### Theorem 5

If $L_{z}$ and $L_{d}$ satisfy the null space condition in (73), then the
proposed offset-free KLMPC system yields zero-offset at steady state.

$\displaystyle\mathcal{N}(L_{d})\subseteq\mathcal{N}(H(I-C(I-A-
BK_{mpc})^{-1}L_{z})$ (73)

where $\mathcal{N}$ denotes null space.

###### Proof 4

First, we derive the relation between $\hat{z}_{\infty}$ and
$\bar{z}_{\infty}$. Subtracting (54) from (51) and rearranging yields

$\displaystyle(I-A)(\hat{z}_{\infty}-\bar{z}_{\infty})=B(u_{\infty}-\bar{u}_{\infty})-L_{z}e_{y,\infty}.$
(74)

Then, substituting (72) into (74) and rearranging yields

$\displaystyle\hat{z}_{\infty}-\bar{z}_{\infty}=-(I-A-
BK_{mpc})^{-1}L_{z}e_{y,\infty}.$ (75)

Now, we examine the offset of the measured controlled variable $Hy_{p,\infty}$
from the set-point $\bar{y}_{c}$ at steady state:

$\displaystyle e_{y_{c},\infty}=Hy_{p,\infty}-\bar{y}_{c}.$ (76)

By substituting (55) into (76), we obtain

$\displaystyle
e_{y_{c},\infty}=H(y_{p,\infty}-C\hat{z}_{\infty}-C_{d}\hat{d}_{\infty}+C(\hat{z}_{\infty}-\bar{z}_{\infty})).$
(77)

Then, substituting (53) and (75) into (77) yields

$\displaystyle e_{y_{c},\infty}=H[I-C(I-A-BK_{mpc})^{-1}L_{z}]e_{y,\infty}.$
(78)

Combining (52) and (78) and rearranging yields

$\displaystyle\begin{bmatrix}0\\\
I\end{bmatrix}e_{y_{c},\infty}=\begin{bmatrix}L_{d}\\\ H(I-C(I-A-
BK_{mpc})^{-1}L_{z})\end{bmatrix}e_{y,\infty}.$ (79)

Then, from the null space condition in (73), $e_{y,\infty}$ that satisfies
$L_{d}e_{y,\infty}=0$ also satisfies $H(I-C(I-A-
BK_{mpc})^{-1}L_{z})e_{y,\infty}=0$. Eventually, the second component of the
right hand side of (79) becomes 0, and thus the following equation holds:

$\displaystyle e_{y_{c},\infty}=0.$ (80)

From (80), we can see that the offset-free KLMPC system yields zero-offset at
steady state.

##  3 Numerical example

In this section, a numerical simulation is performed to demonstrate the
effectiveness of the proposed offset-free KLMPC framework in handling the
plant-model mismatch inherent in the nominal KLMPC. A continuous stirred-tank
reactor (CSTR) in Pannocchia and Rawlings (2003) where a first-order reaction
takes place is considered for the closed-loop simulation. The dynamics of the
CSTR is as follows:

$\displaystyle\frac{dc}{dt}=\dfrac{F_{0}(c_{0}-c)}{\pi
r^{2}h}-k_{0}e^{-\frac{E}{RT}}c$ (81a)
$\displaystyle\dfrac{dT}{dt}=\dfrac{F_{0}(T_{0}-T)}{\pi r^{2}h}-\dfrac{\Delta
H}{\rho C_{p}}k_{0}e^{-\frac{E}{RT}}c+\dfrac{2U}{r\rho C_{p}}(T_{c}-T)$ (81b)
$\displaystyle\dfrac{dh}{dt}=\dfrac{F_{0}-F}{\pi r^{2}}$ (81c)

where the parameters are given in Table 1.

Table 1: Parameters for the CSTR considered in the simulation. Parameter | Value | Unit
---|---|---
Inlet flowrate, $F_{0}$ | 0.1 | $\mathrm{m^{3}/min}$
Inlet temperature, $T_{0}$ | 350 | $\mathrm{K}$
Inlet concentration, $c_{0}$ | 1 | $\mathrm{kmol/m^{3}}$
Frequency factor, $k_{0}$ | 7.2$\times 10^{10}$ | $\mathrm{min^{-1}}$
Specific heat, $C_{p}$ | 0.239 | $\mathrm{kJ/kg\cdot K}$
Activation energy, $E$ | 7.275$\times 10^{4}$ | $\mathrm{kJ/kmol}$
Gas constant, $R$ | 8.314 | $\mathrm{kJ/kmol\cdot K}$
Heat coefficient, $U$ | 54.94 | $\mathrm{kJ/min\cdot m^{2}\cdot K}$
Radius of reactor, $r$ | 0.219 | $\mathrm{m}$
Density, $\rho$ | 1000 | $\mathrm{kg/m^{3}}$
Heat of reaction, $\Delta H$ | -5$\times 10^{4}$ | $\mathrm{kJ/kmol}$

The control objective of the control system is to track the set-point values
of the outlet concentration ($c$) and the reactor temperature ($T$) by
manipulating the jacket temperature ($T_{c}$) and the outlet flow rate ($F$)
of the reactor.

The closed-loop simulation is performed using the MATLAB ®R2019b with Intel
®Core TM i7-10510U CPU @ 1.80 GHz and 16GB RAM.

###  3.1 Data-driven model identification through EDMD

To demonstrate the application of the proposed offset-free KLMPC framework,
the nonlinear model of the CSTR in (81) is utilized as a virtual plant, but
the detailed dynamics of the CSTR is assumed to be not available.

To identify the linear state space model using EDMD in Section 1.1, we
collected input ($T_{c}$, $F$) and state ($c$, $T$, $h$) data from the virtual
CSTR. Specifically, 1000 operation trajectories are generated by solving the
ODEs in (81) for 500 min for each operation with a 1 min sampling time
interval. The initial condition for each operation is randomly assigned around
the following desired steady-state:

$\displaystyle x_{s}=[0.878,324.5,0.659].$ (82)

The input trajectory for each operation is also randomly generated from the
range determined by the following input constraints:

$\displaystyle 290\leq T_{c}\leq 315,\;0.04\leq F\leq 0.16.$

An observable function library $\psi(x)$ in (8) is constructed as follows:

$\displaystyle\psi(x)=[$ $\displaystyle
x_{1},x_{2},x_{3},x_{1}^{2},x_{2}^{2},x_{1}x_{2},$ $\displaystyle
x_{1}e^{-1/x_{2}},(x-x_{s})^{\top}P(x-x_{s})]^{\top}.$ (83)

Based on the discussion in Section 2.1, the state itself and the Lyapunov
function with the equilibrium point $x_{s}$ are included as observable
functions. Utilizing these observables, the Lyapunov constraints can be
transformed into linear forms as in (50g) and (50h).

Figure 1: Input trajectories applied for validation of the linear dynamical
model obtained via EDMD.

Then, a linear state-space model in the lifted space in (13) is identified
using the 500,000 sets of collected $((x_{j},u_{j}),(x_{j}^{+},u_{j}^{+}))$
data and the constructed observables. To validate the prediction power of the
lifted state-space model, we considered a series of step inputs in $T_{c}$ and
pulse inputs in $F$, with a duration of 15 min, as described in Fig. 1.

Figure 2: Output trajectories from the virtual CSTR and the lifted state-space
model induced by the input trajectories in Fig. 1.

Fig. 2 shows the output trajectories from the virtual CSTR in (81) and the
lifted state-space model in (13) induced by the input trajectories in Fig. 1.
While the lifted state-space model seems to predict the behavior of the
virtual CSTR for the most part utilizing the enriched observable space with
nonlinear functions as in (3.1), in Fig. 2, prediction error still exists over
operating time intervals [15 30] and [30 45] mins due to the inherent plant-
model mismatch caused by the finite-dimensional approximation of the infinite-
dimensional Koopman operator via EDMD.

To quantitatively examine the prediction performance of the lifted state-space
model, we derive a normalized root mean squared error (NRMSE) value for each
output variable as follows:

$\displaystyle NRMSE:=\frac{RMSE}{y_{p,max}-y_{p,min}}$ (84a) $\displaystyle
RMSE:=\left(\frac{\sum^{n_{D}}_{j=1}(\hat{y}_{m,j}-y_{p,j})^{2}}{n_{D}}\right)^{1/2}$
(84b)

where $y_{p}$ and $\hat{y}_{m}$ denote the output from the virtual CSTR and
the lifted state-space model, respectively. The resultant NRMSE value for each
output is obtained as

$\displaystyle NRMSE_{y}=[0.1319,0.0969,0.0142].$

The NRMSE value of $h$ is considerably smaller than those values of $C$ and
$T$, and we can also see the prediction accuracy of $h$ is much higher than
other variables in Fig. 2. This can be attributed to the relatively simple and
linear dynamics of $h$ in (81c) that allows EDMD to derive an accurate linear
predictor for $h$.

###  3.2 Closed-loop simulation results

Based on the lifted state-space model derived from data-driven EDMD, we design
an offset-free KLMPC system for reference tracking of the CSTR.

As previously discussed in Section 2, a disturbance model is augmented with
the Koopman-based lifted state space model with proper $B_{d}$ and $C_{d}$
that satisfy the condition in Proposition 1 to handle the plant-model mismatch
effect. Subsequently, a disturbance estimator in (2.1), a target problem in
(33), and an optimal control problem in (50) are designed based on the offset-
free KLMPC formulation in Section 2. The prediction horizon of the model
predictive controller is set as $N=10$. Additionally, the following process
constraints for output variables are applied to the optimal control problem:

$\displaystyle 0.81\leq c\leq 0.92,\;320\leq T\leq 330,\;0.4\leq h\leq 1.2.$

Figure 3: Reference tracking results from the proposed offset-free KLMPC and
the nominal KLMPC systems.

Fig. 3 presents the closed-loop trajectories of the output variables induced
by the proposed offset-free KLMPC and the nominal KLMPC schemes. The set-point
value for $c$ is changed for every 25 min over [0.85, 0.90] while the set-
point value for $T$ is fixed at 324.5 K. The measured $h$ value is applied to
the disturbance estimator, but it is not treated as a controlled variable. We
can see that $c$ and $T$ trajectories induced by the nominal KLMPC scheme show
considerable offsets from the set-point values due to the plant-model mismatch
between the linear model derived by EDMD in the previous section and the CSTR
process in (81). On the other hand, $c$ and $T$ trajectories induced by the
proposed offset-free KLMPC scheme show zero steady-state offset.

Figure 4: Estimated disturbance trajectories of the proposed offset-free KLMPC
scheme.

This is because the proposed offset-free KLMPC framework can effectively
compensate for the influence of plant-model mismatch by properly estimating
the disturbance signals as shown in Fig. 4 and applying them to the control
system. When the set-point changes, it takes some time to compute disturbance
values via the estimator as shown in Fig. 4, e.g., around 50, 75 and 150 mins.
Eventually, the trajectories of controlled variables $C$ and $T$ in Fig. 3
show considerable deviations from the set-point values near those sampling
instants. Nonetheless, as the control system reaches a steady state for each
set-point, the disturbance values derived from the properly designed stable
estimator converge to proper values as shown in Fig. 4. By applying these
disturbance values to the control system, the proposed offset-free KLMPC
scheme can derive an optimal input to drive the CSTR system to a steady state
where the zero steady-state offset is satisfied in the presence of plant-model
mismatch on the controlled variables.

## Conclusions

In this work, we proposed the offset-free KLMPC framework to address plant-
model mismatch in the KLMPC scheme which is inherently occurred in data-driven
model identification through EDMD caused by a finite-dimensional approximation
to the infinite-dimensional Koopman operator.

Specifically, a disturbance model is augmented with the lifted state-space
model derived from EDMD to consider the plant-model mismatch. Then, the
offset-free KLMPC framework that consists of an estimator, a target problem,
and an optimal control problem is developed. The estimator derives lifted
state and disturbance estimates from the measured prediction error, and the
target problem derives proper targets for lifted state and input considering
the effect of disturbance. Then, the optimal control problem with Lyapunov
constraints obtains an optimal input to track the set-point while compensating
for the plant-model mismatch, and ensuring feasibility and stability of the
control system.

Unlike the nominal KLMPC scheme, since the equilibrium point of the optimal
control problem (i.e., steady-state target) is continuously updated with the
disturbance estimate in the proposed framework, the formulations for updating
Lyapunov function and stabilizing control law in Lyapunov constraints are also
developed. Additionally, a mathematical analysis for the zero steady-state
offset condition of the proposed framework is performed.

The numerical simulation results showed the effectiveness of the proposed
offset-free KLMPC framework in addressing the plant-model mismatch compared to
the nominal KLMPC scheme. This is especially meaningful because the inaccuracy
in prediction due to the plant-model mismatch can be exacerbated when a system
state is lifted to the space of observables. In conclusion, the developed
offset-free KLMPC framework is expected to provide a valuable direction for
addressing the plant-model mismatch that is inherent in Koopman operator-based
predictive control methods.

## Acknowledgments

Financial support from the Artie McFerrin department of chemical engineering
and the Texas A&M Energy Institute are gratefully acknowledged.

##  Appendix A Derivation of (62)

We can reformulate $P_{\hat{z}_{\infty}}$ in a compact form by expressing the
predicted lifted states in the objective function and inequality constraints
explicitly with $U=[u_{0}^{\top},\cdots,u_{N-1}^{\top}]^{\top}$:

$\displaystyle P_{\hat{z}_{\infty}}:\;\underset{\small{U}}{\mathrm{min}}$
$\displaystyle\quad\frac{1}{2}U^{\top}H_{\infty}U+f_{\infty}^{\top}U$ (85a)
$\displaystyle\mathrm{s.t.}$ $\displaystyle\quad G_{\infty}U\leq
g_{\infty}+S_{\infty}\hat{z}_{\infty}$ (85b)

Then, we can derive the first-order Karush-Kuhn-Tucker (KKT) optimality
conditions of $P_{\hat{z}_{\infty}}$ as

$\displaystyle H_{\infty}U^{*}+f_{\infty}+G_{\infty}^{\top}v^{*}=0$ (86a)
$\displaystyle
v_{i}^{*}(G_{\infty,i}U^{*}-g_{\infty,i}-S_{\infty,i}\hat{z}_{\infty})=0$
(86b) $\displaystyle v_{i}^{*}\geq 0$ (86c) $\displaystyle
G_{\infty}U^{*}-g_{\infty}-S_{\infty}\hat{z}_{\infty}\leq 0$ (86d)

where $v$ denotes the dual variable, and $i$ denotes the index of constraints.

By rearranging (86a), we can obtain

$\displaystyle U^{*}=-H_{\infty}^{-1}(G_{\infty}^{\top}v^{*}+f_{\infty}).$
(87)

Then, substituting (87) into (86b) yields

$\displaystyle
v_{i}^{*}(-G_{\infty,i}H^{-1}_{\infty}(G_{\infty}^{\top}v^{*}+f_{\infty})-g_{\infty,i}-S_{\infty,i}\hat{z}_{\infty})=0.$
(88)

When the Lyapunov function constraint in (50h) is active at the optimum, (89)
holds by complementary slackness.

$\displaystyle-
G_{\infty,i_{\ell}}H^{-1}_{\infty}(G_{\infty,i_{\ell}}^{\top}v^{*}_{i_{\ell}}+f_{\infty})-g_{\infty,i_{\ell}}-S_{\infty,i_{\ell}}\hat{z}_{\infty}=0$
(89)

where

$\displaystyle G_{\infty,i_{\ell}}:=F_{v}BC_{1}$ $\displaystyle
g_{\infty,i_{\ell}}:=F_{v}B(\bar{N}\bar{z}_{\infty}-K_{d}\hat{d}_{\infty})$
$\displaystyle S_{\infty,i_{\ell}}:=-F_{v}BK_{z}$ $\displaystyle
F_{v}:=D_{v}+2(\bar{z}_{s}-\bar{z}_{\infty})^{\top}D_{x}^{\top}Q_{v}D_{x}$
$\displaystyle C_{1}:=\begin{bmatrix}I_{n_{u}},0,\cdots,0\end{bmatrix}.$

Rearranging (89) yields

$\displaystyle v^{*}_{i_{\ell}}=-$
$\displaystyle(G_{\infty,i_{\ell}}H^{-1}_{\infty}G_{\infty,i_{\ell}}^{\top})^{-1}$
(90)
$\displaystyle(G_{\infty,i_{\ell}}H^{-1}_{\infty}f_{\infty}+g_{\infty,i_{\ell}}+S_{\infty,i_{\ell}}\hat{z}_{\infty})$

Then, substituting (90) and $v^{*}_{i}=0\;\forall i\neq i_{\ell}$ into (87)
yields

$\displaystyle
U^{*}=-H_{\infty}^{-1}(G_{\infty,i_{\ell}}^{\top}v^{*}_{i_{\ell}}+f_{\infty}).$
(91)

Finally, we obtain $u_{\infty}$ from (91):

$\displaystyle
u_{\infty}=C_{1}U^{*}=K_{i_{\ell}}\hat{z}_{\infty}+c_{i_{\ell}}$ (92)

where

$\displaystyle K_{i_{\ell}}:=C_{1}(M_{1}\Psi^{\top}Q_{\mathbf{z}}\Phi-
M_{2}S_{\infty,i_{\ell}})$ $\displaystyle
c_{i_{\ell}}:=C_{1}(M_{1}(\Psi^{\top}Q_{\mathbf{z}}(\Psi_{d}\hat{\mathbf{d}}_{\infty}-\bar{\mathbf{z}}_{\infty})-Q_{\mathbf{u}}\bar{\mathbf{u}}_{\infty})-M_{2}g_{\infty,i_{\ell}})$
$\displaystyle
M_{1}:=H_{\infty}^{-1}G_{\infty,i_{\ell}}^{\top}(G_{\infty,i_{\ell}}H^{-1}_{\infty}G_{\infty,i_{\ell}}^{\top})^{-1}G_{\infty,i_{\ell}}H^{-1}_{\infty}-H_{\infty}^{-1}$
$\displaystyle
M_{2}:=-H_{\infty}^{-1}G_{\infty,i_{\ell}}^{\top}(G_{\infty,i_{\ell}}H^{-1}_{\infty}G_{\infty,i_{\ell}}^{\top})^{-1}.$

## References

  * Alessio and Bemporad (2009) Alessio, A., Bemporad, A., 2009\.  A survey on explicit model predictive control, in: Nonlinear model predictive control. Springer, pp. 345–369.
  * Arbabi et al. (2018) Arbabi, H., Korda, M., Mezić, I., 2018.  A data-driven Koopman model predictive control framework for nonlinear partial differential equations, in: 2018 IEEE Conference on Decision and Control (CDC), IEEE. pp. 6409–6414.
  * Bemporad et al. (2002) Bemporad, A., Morari, M., Dua, V., Pistikopoulos, E.N., 2002\.  The explicit linear quadratic regulator for constrained systems.  Automatica 38, 3–20.
  * Betti et al. (2013) Betti, G., Farina, M., Scattolini, R., 2013.  A robust MPC algorithm for offset-free tracking of constant reference signals.  IEEE Transactions on Automatic Control 58, 2394–2400.
  * Borrelli et al. (2017) Borrelli, F., Bemporad, A., Morari, M., 2017.  Predictive control for linear and hybrid systems.  Cambridge University Press.
  * Budišić et al. (2012) Budišić, M., Mohr, R., Mezić, I., 2012.  Applied Koopmanism.  Chaos: An Interdisciplinary Journal of Nonlinear Science 22, 047510.
  * Christofides et al. (2013) Christofides, P.D., Scattolini, R., de la Peña, D.M., Liu, J., 2013\.  Distributed model predictive control: A tutorial review and future research directions.  Computers & Chemical Engineering 51, 21–41.
  * Davison and Smith (1971) Davison, E., Smith, H., 1971\.  Pole assignment in linear time-invariant multivariable systems with constant disturbances.  Automatica 7, 489–498.
  * Johnson (1970) Johnson, C., 1970.  Further study of the linear regulator with disturbances–The case of vector disturbances satisfying a linear differential equation.  IEEE Transactions on Automatic Control 15, 222–228.
  * Koopman (1931) Koopman, B.O., 1931.  Hamiltonian systems and transformation in Hilbert space.  Proceedings of the National Academy of Sciences of the United States of America 17, 315–318.
  * Koopman and Neumann (1932) Koopman, B.O., Neumann, J.V., 1932\.  Dynamical systems of continuous spectra.  Proceedings of the National Academy of Sciences of the United States of America 18, 255–263.
  * Korda and Mezić (2018) Korda, M., Mezić, I., 2018\.  Linear predictors for nonlinear dynamical systems: Koopman operator meets model predictive control.  Automatica 93, 149–160.
  * Maeder et al. (2009) Maeder, U., Borrelli, F., Morari, M., 2009.  Linear offset-free model predictive control.  Automatica 45, 2214–2222.
  * Maeder and Morari (2010) Maeder, U., Morari, M., 2010\.  Offset-free reference tracking with model predictive control.  Automatica 46, 1469–1476.
  * Mauroy and Mezić (2012) Mauroy, A., Mezić, I., 2012\.  On the use of Fourier averages to compute the global isochrons of (quasi) periodic dynamics.  Chaos: An Interdisciplinary Journal of Nonlinear Science 22, 033112.
  * Mhaskar et al. (2006) Mhaskar, P., El-Farra, N.H., Christofides, P.D., 2006.  Stabilization of nonlinear systems with state and control constraints using Lyapunov-based predictive control.  Systems & Control Letters 55, 650–659.
  * Morari and Maeder (2012) Morari, M., Maeder, U., 2012\.  Nonlinear offset-free model predictive control.  Automatica 48, 2059–2067.
  * Muske and Badgwell (2002) Muske, K.R., Badgwell, T.A., 2002\.  Disturbance modeling for offset-free linear model predictive control.  Journal of Process Control 12, 617–632.
  * Narasingam and Kwon (2019) Narasingam, A., Kwon, J.S.I., 2019\.  Koopman Lyapunov-based model predictive control of nonlinear chemical process systems.  AIChE Journal 65, e16743.
  * Narasingam and Kwon (2020a) Narasingam, A., Kwon, J.S.I., 2020a.  Application of Koopman operator for model-based control of fracture propagation and proppant transport in hydraulic fracturing operation.  Journal of Process Control 91, 25–36.
  * Narasingam and Kwon (2020b) Narasingam, A., Kwon, J.S.I., 2020b.  Data-driven feedback stabilization of nonlinear systems: Koopman-based model predictive control.  arXiv preprint arXiv:2005.09741 .
  * Pannocchia and Bemporad (2007) Pannocchia, G., Bemporad, A., 2007\.  Combined design of disturbance model and observer for offset-free model predictive control.  IEEE Transactions on Automatic Control 52, 1048–1053.
  * Pannocchia et al. (2015) Pannocchia, G., Gabiccini, M., Artoni, A., 2015.  Offset-free MPC explained: Novelties, subtleties, and applications.  IFAC-PapersOnLine 48, 342–351.
  * Pannocchia and Rawlings (2003) Pannocchia, G., Rawlings, J.B., 2003\.  Disturbance models for offset-free model-predictive control.  AIChE journal 49, 426–437.
  * Peitz and Klus (2019) Peitz, S., Klus, S., 2019.  Koopman operator-based model reduction for switched-system control of PDEs.  Automatica 106, 184–191.
  * de la Peña and Christofides (2008) de la Peña, D.M., Christofides, P.D., 2008\.  Lyapunov-based model predictive control of nonlinear systems subject to data losses.  IEEE Transactions on Automatic Control 53, 2076–2089.
  * Preuss (1982) Preuss, H.P., 1982.  Perfect steady-state tracking and disturbance rejection by constant state feedback.  International Journal of Control 35, 75–94.
  * Proctor et al. (2018) Proctor, J.L., Brunton, S.L., Kutz, J.N., 2018.  Generalizing Koopman theory to allow for inputs and control.  SIAM Journal on Applied Dynamical Systems 17, 909–930.
  * Rawlings et al. (2017) Rawlings, J.B., Mayne, D.Q., Diehl, M., 2017.  Model predictive control: Theory, computation, and design. volume 2.  Nob Hill Publishing Madison, WI.
  * Rowley et al. (2009) Rowley, C.W., Mezić, I., Bagheri, S., Schlatter, P., Henningson, D.S., et al., 2009.  Spectral analysis of nonlinear flows.  Journal of fluid mechanics 641, 115–127.
  * Son et al. (2020a) Son, S.H., Oh, T.H., Kim, J.W., Lee, J.M., 2020a.  Move blocked model predictive control with improved optimality using semi-explicit approach for applying time-varying blocking structure.  Journal of Process Control 92, 50–61.
  * Son et al. (2020b) Son, S.H., Park, B.J., Oh, T.H., Kim, J.W., Lee, J.M., 2020b.  Move blocked model predictive control with guaranteed stability and improved optimality using linear interpolation of base sequences.  International Journal of Control , 1–13.
  * Sootla et al. (2018) Sootla, A., Mauroy, A., Ernst, D., 2018.  Optimal control formulation of pulse-based control using Koopman operator.  Automatica 91, 217–224.
  * Wang (2009) Wang, L., 2009.  Model predictive control system design and implementation using MATLAB®.  Springer Science & Business Media.
  * Williams et al. (2016) Williams, M.O., Hemati, M.S., Dawson, S.T., Kevrekidis, I.G., Rowley, C.W., 2016.  Extending data-driven Koopman analysis to actuated systems.  IFAC-PapersOnLine 49, 704–709.
  * Williams et al. (2015) Williams, M.O., Kevrekidis, I.G., Rowley, C.W., 2015.  A data–driven approximation of the Koopman operator: Extending dynamic mode decomposition.  Journal of Nonlinear Science 25, 1307–1346.


# Evaluating the Significance of Outdoor Advertising from Driver’s Perspective
Using Computer Vision

1st Zuzana Černeková  Faculty of Mathematics Physics and informatics
Comenius University
Bratislava, Slovakia
<EMAIL_ADDRESS>
0000-0002-7617-4192  2nd Zuzana Berger Haladová  Faculty of Mathematics
Physics and informatics
Comenius University
Bratislava, Slovakia
<EMAIL_ADDRESS>
0000-0002-5947-8063  3rd Ján Špirka  Faculty of Mathematics Physics and
informatics
Comenius University
Bratislava, Slovakia
<EMAIL_ADDRESS>4th Viktor Kocur  Faculty of Mathematics Physics and
informatics
Comenius University
Bratislava, Slovakia
<EMAIL_ADDRESS>
0000-0001-8752-2685

###### Abstract

Outdoor advertising, such as roadside billboards, plays a significant role in
marketing campaigns but can also be a distraction for drivers, potentially
leading to accidents. In this study, we propose a pipeline for evaluating the
significance of roadside billboards in videos captured from a driver’s
perspective. We have collected and annotated a new BillboardLamac dataset,
comprising eight videos captured by drivers driving through a predefined path
wearing eye-tracking devices. The dataset includes annotations of billboards,
including 154 unique IDs and 155 thousand bounding boxes, as well as eye
fixation data. We evaluate various object tracking methods in combination with
a YOLOv8 detector to identify billboard advertisements with the best approach
achieving 38.5 HOTA on BillboardLamac. Additionally, we train a random forest
classifier to classify billboards into three classes based on the length of
driver fixations achieving 75.8% test accuracy. An analysis of the trained
classifier reveals that the duration of billboard visibility, its saliency,
and size are the most influential features when assessing billboard
significance.

###### Index Terms:

neural networks, object detectors, saliency, eyetracking, billboards
detection, advertisement classification

††publicationid: pubid: 979-8-3503-4353-3/23/$31.00 ©2023 IEEE

##  I Introduction

Acquisition advertising is frequently employed as part of a marketing strategy
to increase brand awareness or support ongoing campaigns. Businesses utilize
various advertising platforms such as television, radio, and social networks.
One common outdoor advertising method is roadside advertising (billboards),
positioned along roads and highways to attract the attention of drivers and
other commuters. Additionally, advertisements can be found on buildings,
fences, bus stops, and even vehicles, creating a visual presence in road
traffic. Unfortunately, these ads are often considered visual pollution on the
roads.

The impact of advertising on individuals depends on various factors, including
the type of advertising, strategic placement, and the content itself.
Billboards undoubtedly capture the attention of drivers, benefiting marketing
efforts. However, it is important to acknowledge that even a momentary loss of
driver attention can result in traffic accidents with severe consequences.
Several studies, including our work, deal with the issue of advertisement
safety.

In this paper, we have created a pipeline for the estimation of billboard
advertisement’s significance. The pipeline consists of billboard detection and
tracking and classification of the billboard significance which is evaluated
by data collected in our user study utilizing eyetracking technology.

The paper is organized as follows. In the second section, we discuss the
previous works on billboard advertisements’ characteristics and significance.
In the third section, we describe our dataset BillboardLamac. The fourth
section is dedicated to the billboard tracking methods description and
evaluation. The fifth section discusses the classification of the billboard’s
significance.

##  II Billboards

The following section explores the previous works on the effectiveness of
electronic billboards and traditional billboards, considering factors such as
gaze duration, placement, and content impact on driver attention and road
safety.

Several studies indicate that electronic billboards achieve longer gaze
duration compared to traditional non-illuminated advertising [1, 2, 3]. This
type of advertisement is the most dangerous due to the average duration of
fixations [4, 5].

The placement of advertising significantly influences its reach.
Advertisements located on the driver’s side of the road recorded more
fixations than those on the opposite side [6]. The assumption that the height
of the advertisement placement prolongs fixation duration was negatively
confirmed in studies [6, 7]. The authors found that drivers looked more at
advertisements at road level than those elevated three or more meters above
the ground. An exception occurs when drivers deliberately seek out
advertisements, as they tend to look higher. Furthermore, some studies
considered the angle of the advertising surface [8], the complexity of the
traffic section, the type of environment, and its complexity [9, 6]. It was
found that all these elements to some extent influence the visibility of
advertisements.

The content of the advertisement is another important characteristic. Studies
[10, 11] confirm the assumption that advertisements with longer text achieve
longer fixations. This also applies to advertisements with a sexual undertone
[12] and advertisements featuring human beings [13]. Moreover, the more
complex the advertisement, the longer the driver’s gaze on it [14]. The
content of the advertisement that elicits negative emotions in drivers was
also examined [15]. It was found that such content also leads to longer
fixations and, at the same time, reduces speed. Additionally, questionnaires
filled out by drivers after their drive showed that negative advertisements
are the most memorable.

##  III Dataset

To test our method for billboard detection, tracking, and saliency evaluation
we have created a new BillboardLamac dataset which consists of 8 videos
captured by an eyetracker worn by 8 different drivers on the same round trip
in Lamač, Bratislava, Slovakia. A map of the round trip taken by the
participants of the user study is shown in Fig. 2. We have annotated the
billboards in the videos, for each billboard on the roundtrip we have a unique
ID (which is consistent through all the videos), and the bounding box for each
frame where it is visible. There are 145 unique billboards annotated with 155
912 bounding boxes in 8 videos of an average length of 24 minutes. An example
of an image from the dataset with displayed annotations can be seen in Fig. 1.
The dataset is available on demand and includes videos from the eyetracker, a
text file with bounding boxes and billboards’ IDs with the number of
corresponding frames of the video, and the text file with positions of the eye
fixations for each frame of the video.

Figure 1: A frame from the dataset BillboardLamac with annotated billboards
(ID 67 and 68) and their corresponding bounding boxes. Figure 2: Map of the
round trip taken by the participants of the user study displayed in Google
Maps.

###  III-A Eyetracking

To capture salient areas in the images and videos, we utilized the Tobii Pro
Glasses 3, a wearable eye-tracking device providing robust and precise gaze
data collection, all while allowing users unrestricted mobility. These glasses
feature 16 illuminators and 4 eye cameras integrated into scratch-resistant
lenses, ensuring optimal eye-tracking performance and an unobstructed view for
the wearer. Additionally, the scene camera has a wide field of view (106° H:
95°, V: 63°). An example of the fixations captured by Tobii Pro Glasses
displayed on the frame from a video captured in the user study can be seen in
Fig. 3

Figure 3: Fixations captured by Tobii Pro Glasses displayed on a frame from a
video captured in the user study.

##  IV Billboard tracking

To achieve information on the saliency of different billboards in the driving
scenario it is necessary to analyse them not only as individual frames but as
tracked instances. We decided to test YOLOv8large [16] in combination with
different trackers based on the SORT method [17].

###  IV-A Detector

The YOLOv8large [16] billboard detector was trained on Mappilary Vistas
Dataset [18]. Most of the images in Mapillary Vistas are captured from the
interior or rooftop of vehicles, showcasing traffic within urban areas with
numerous buildings, sidewalks, traffic signs, and various road objects,
including advertising surfaces. We excluded images without advertisements,
resulting in a reduced dataset of approximately 20,000 images. During the
initial training, we noticed higher error rates specifically for smaller
labeled advertisements. The smallest ads were labeled at such a distance that
it would be uncertain whether the driver was looking at the advertisement or
another object. Therefore, we created a second version of the dataset that
only included images where the advertising surface accounted for at least 0.1%
of the image, resulting in almost half the number of images compared to the
original dataset. To evaluate the detector we used the standard mAP50 metric.
We have evaluated the trained YOLOv8large on both Mapilary Vistas and our
dataset BillboardLamac and achieved 0.503 mAP50 and 0.656 mAP50 scores
respectively. The higher score on BillboardLamac was gained due to less
challenging billboards (in-occluded, bigger, etc.) in our dataset. Billboard’s
bounding boxes with corresponding confidences acquired by YOLOv8large on the
BillboardLamac dataset can be seen in Fig. 4.

Figure 4: Frame from our BillboardLamac dataset with 3 billboards detected by
the Yolov8large model [16] trained on a reduced version of the Mapillary
Vistas dataset [18].

.

###  IV-B Tracking

To track billboards, we have tested the methods OC-SORT [19], Deep OC-SORT
[20], StrongSORT [21], BoTSORT [22], and ByteTrack [23]. All methods are based
on the SORT (Simple Online and Realtime Tracking) method, which employs the
Kalman filter and Hungarian algorithm for object tracking [24]. SORT is a
tracking-by-detection framework for the problem of multiple object tracking
(MOT) where objects are detected in each frame and represented as bounding
boxes and the tracker associates them (see Fig. 5). In our approach, the input
video is processed frame by frame, and detections are obtained using the
trained YOLOv8 model and associated with references. Multiple objects and
classes can be tracked simultaneously. As shown in Fig. 5, the output of the
tracking method consists of detections with assigned references, describing
the trajectory of the tracked object.

Figure 5: Tracking by detection. Objects are first detected via bounding boxes
which are then associated into separate tracks.

To evaluate object tracking we decided to use the HOTA (Higher Order Tracking
Accuracy) metric [25], which consists of smaller components that assess three
aspects of tracking: localization, detection, and association. This
formulation ensures that both detection and association, unlike in many other
tracking metrics, are evenly balanced, and the final score lies somewhere
between them. The results for each method are recorded in Table I.

TABLE I: Comparison of the evaluated tracking methods. We report metrics proposed in [25]. Method | HOTA | DetPr | DetRe | DetA | AssPr | AssRe | AssA | Loc
---|---|---|---|---|---|---|---|---
OC-SORT [19] | 36.8 | 39.5 | 34.3 | 53.2 | 58.2 | 36.3 | 85.1 | 89.1
DeepOC-SORT [20] | 35.8 | 38.9 | 33.0 | 54.4 | 55.6 | 34.8 | 84.7 | 89.1
StrongSORT [21] | 38.5 | 41.4 | 35.8 | 54.6 | 60.6 | 42.7 | 72.9 | 88.7
BotSORT [22] | 37.0 | 39.5 | 34.7 | 53.0 | 58.4 | 37.5 | 82.6 | 89.2
ByteTrack [23] | 32.6 | 36.4 | 29.3 | 44.4 | 62.8 | 30.5 | 84.7 | 86.6

The StrongSORT method performed the best in the overall HOTA metric and in
most of the other metrics. However, in terms of association precision, the
ByteTrack method achieved the highest scores, although it obtained lower
values in almost all other metrics. For calculating the significance of
billboards, the StrongSORT method proved to be the most effective, achieving
the highest number of detections with relatively good reference assignments.

###  IV-C Billboard Saliency

If a billboard was detected in the frame, we analyzed the driver’s gaze. The
eyetracker captures up to two eye gaze fixation points per one video frame. In
cases where no fixation point was recorded for a frame, we estimated it by
extrapolating the positions of the two previous fixation points. If the
information from the previous frame was also not available we did not consider
the frame for further analysis. For frames where both a billboard was detected
and a fixation point was recorded, we measured their overlap. This provided
information on the number of frames where the driver looked at the billboard.
By multiplying the number of frames by the duration of one frame, we obtained
the gaze duration. Based on the viewing duration, we assessed the salience of
the billboard for each driver individually and then calculated the final
salience as median across the drivers.

A total of 145 billboards were identified along the chosen route. Table II
presents the time range (represented by the number of frames with the
billboard), the count of billboards in each category, and the average number
of drivers who viewed the billboard for each category. The categorization was
determined based on the outcomes of the StrongSORT tracking method.

On average, each billboard was viewed by three drivers, while 45
advertisements were not viewed by any driver. There were only two billboards
viewed by all drivers.

TABLE II: We split the detected billboards into 3 saliency categories based on the median gaze duration across the drivers. We also denote the average number of drivers who have looked at the billboards of the given cateogory. Category | Gaze duration | # of billboards | Avg. # of drivers
---|---|---|---
None | 0 ms | 45 | -
Short | 1-249 ms | 79 | 3
Long | 250+ ms | 21 | 3

##  V Classification of Billboard significance

We trained a Random Forest Classifier (RFC) to predict the significance of
billboards. Data were divided into a training and a testing set. The training
set contained 115 billboards and the test set had 30 billboards.

###  V-A Feature extraction

We extracted a total of 7 features for each billboard to be used as input into
the classifier. Each billboard was at most detected in 8 sequences. Not all of
the billboards were detected in all videos, since billboards may have been
hidden by cars or be present at such angles that they could not be seen when
the driver was in the most distant lane of traffic. We calculated several
features listed below for each billboard. To combine the data for each
billboard from different driving sessions we have calculated the average over
the sessions where the billboard was visible.

####  V-A1 Duration of visibility

We included a feature indicating how many frames the ad has been visible based
on the output of the tracker. We expect that the longer an ad has been
visible, the more likely it is to be noticed by the driver.

####  V-A2 Advertisement localization

Two features were computed from the localization of the advertisement. The
first one determined in which part of the image the ad was located the
longest. The image was divided into three sections: left, right, and center
with ratios of 40:20:40. As a second feature, we also calculated the Euclidean
distance of the center point of the ad from the center point of the image. We
expect the position of the billboard to have a significant effect on the
driver’s attention.

####  V-A3 Billboard size

Since the size of the advertisement increased as the vehicle approached the
advertisement, we calculated the size of the advertisement as the average size
over the entire sequence where the advertisement was visible. We also included
a separate feature with the size calculated as the average for the 10 frames
where the detected billboard bounding boxes were the largest.

####  V-A4 Saliency map

Figure 6: Example of overlaying the saliency map estimated by H2DS [26] with a
frame from the experiment.

An existing neural network HD2S [26] was used to generate saliency maps for
each frame in the videos. The network outputs a map indicating the saliency of
each frame. As depicted in Fig. 6, the saliency map effectively highlights the
region containing the billboard. We evaluated the generated saliency map
against the ground truth saliency obtained using the eyetracker with common
metrics [27] resulting in AUC-Judd accuracy of $0.8964$, NSS of $1.7823$, SIM
score of $0.3332$, and CC score of $0.4367$.

We included two features extracted from the saliency maps. One feature was
calculated as the average of the saliency values within the detected area of
the billboard and the second feature was calculated as the ratio of average
saliency values in the billboard and the whole frame. Because of the many low
or zero values in a frame that significantly reduce the average, all of the
mentioned averages were calculated only for points with values greater than 50
(from range 0-255). To obtain these features for a given billboard we averaged
them over the whole sequence of frames where the billboards were visible.

###  V-B Classification results

We used a random forest classifier for classification. To optimize the
training hyperparameters, we used 5-fold cross-validation. The trained
classifier used 100 trees with a maximum depth of 2.

|  | Predicted Class
---|---|---
|  | 0 ms | 0-250 ms | 250+ ms
True Class | 0 ms | 6 | 3 | 0
0 - 250 ms | 1 | 15 | 0
250+ ms | 0 | 3 | 1
TABLE III: The confusion matrix of the trained random forest classifier on the
testing data.

The resulting classifier accurately predicted the class for billboard
saliency, and if incorrect, it was typically misclassified as a neighboring
class. The overall classification performance achieved an accuracy of 75.8% on
the test set. The confusion matrix is shown in Table III.

Figure 7: Feature importance for the trained random forest classifier
calculated based on feature permutation (Permutation) and mean decrease in
impurity (MDI).

To analyse the impact of the features used for classification we performed
feature importance analysis based on feature permutation and mean decrease in
impurity. The resulting importance values can be seen in Fig. 7.

Among the various features considered, the duration of billboard visibility
emerged as the most influential factor. This finding indicates that the length
of time a billboard remains visible to drivers greatly impacts its perceived
significance. Additionally, features derived from the saliency map, obtained
using the HD2S method [26], demonstrated a notable impact on billboard
classification. This suggests that the visual prominence of billboards, as
captured by saliency analysis, plays a crucial role in determining their
effectiveness. Furthermore, we observed that the relative size of the
billboard in the ten frames where it occupied the largest area also
contributed significantly to the classification outcome.

By identifying these key features, our analysis provides valuable insights
into the factors that contribute to the perceived significance of billboards.
This information can guide future billboard design and placement strategies,
enabling advertisers to optimize their campaigns for maximum impact and driver
engagement.

##  VI Conclusion

We have presented a pipeline for evaluating the significance of roadside
billboards in videos captured from a driver’s perspective. Our approach does
not require the use of eye-tracking devices during the actual pipeline
execution.

Through the development and annotation of the BillboardLamac dataset,
consisting of eight videos captured by drivers along a predefined path, we
were able to train and evaluate our pipeline effectively. By evaluating
various object tracking methods with the YOLOv8 detector, we achieved good
performance in identifying billboard advertisements, with our best approach
achieving a 38.5 HOTA score on BillboardLamac.

Moreover, we employed a random forest classifier to classify billboards into
three distinct classes based on the length of driver fixations. This
classifier demonstrated strong performance, achieving a test accuracy of
75.8%. Notably, our analysis of the trained classifier revealed that the
duration of billboard visibility, saliency features derived from a neural
network, and the relative size of the billboard were the most influential
factors when assessing billboard significance.

The findings from our feature importance analysis provide valuable insights
for advertisers and marketers. Understanding the impact of factors such as
visibility duration, visual saliency, and size allows for informed decisions
regarding billboard design and placement strategies. By optimizing these
factors, advertisers can enhance the effectiveness of their campaigns in
capturing driver attention and engagement.

Overall, our pipeline offers a practical solution for evaluating the
significance of roadside billboards, contributing to safer driving experiences
and more effective marketing campaigns. Future work may involve expanding the
dataset, incorporating additional features, and exploring other machine-
learning techniques to further enhance the accuracy and applicability of the
pipeline.

## Acknowledgment

This publication is the result of support under the Operational Program
Integrated Infrastructure for the project: Advancing University Capacity and
Competence in Research, Development and Innovation (ACCORD,
ITMS2014+:313021X329), co-financed by the European Regional Development Fund.
The work presented in this paper was carried out in the framework of the
TERAIS project, a Horizon-Widera-2021 program of the European Union under the
Grant agreement number 101079338.

## References

  * [1] O. Oviedo-Trespalacios, V. Truelove, B. Watson, and J. A. Hinton, “The impact of road advertising signs on driver behaviour and implications for road safety: A critical systematic review,” _Transportation Research Part A: Policy and Practice_ , vol. 122, pp. 85–98, 2019. [Online]. Available: https://www.sciencedirect.com/science/article/pii/S0965856418310632
  * [2] D. Beijer, A. Smiley, and M. Eizenman, “Observed driver glance behavior at roadside advertising signs,” _Transportation Research Record_ , vol. 1899, no. 1, pp. 96–103, 2004.
  * [3] R. Brome, M. Awad, and N. M. Moacdieh, “Roadside digital billboard advertisements: Effects of static, transitioning, and animated designs on drivers’ performance and attention,” _Transportation research part F: traffic psychology and behaviour_ , vol. 83, pp. 226–237, 2021.
  * [4] K. Yellappan, Y. Ghani, M. Musa, M. Siam, and C. Tan, “Exposure and perception on distraction towards roadside digital advertisements,” 2016.
  * [5] A. Smiley, B. Persaud, G. Bahar, C. Mollett, C. Lyon, T. Smahel, and W. L. Kelman, “Traffic safety evaluation of video advertising signs,” _Transportation research record_ , vol. 1937, no. 1, pp. 105–112, 2005.
  * [6] M. Costa, L. Bonetti, V. Vignali, A. Bichicchi, C. Lantieri, and A. Simone, “Driver’s visual attention to different categories of roadside advertising signs,” _Applied ergonomics_ , vol. 78, pp. 127–136, 2019.
  * [7] D. Crundall, E. Van Loon, and G. Underwood, “Attraction and distraction of attention with roadside advertisements,” _Accident Analysis & Prevention_, vol. 38, no. 4, pp. 671–677, 2006.
  * [8] M. Zalesinska, “The impact of the luminance, size and location of led billboards on drivers’ visual performance—laboratory tests,” _Accident Analysis & Prevention_, vol. 117, pp. 439–448, 2018.
  * [9] K. Mollu, J. Cornu, K. Brijs, A. Pirdavani, and T. Brijs, “Driving simulator study on the influence of digital illuminated billboards near pedestrian crossings,” _Transportation research part F: traffic psychology and behaviour_ , vol. 59, pp. 45–56, 2018.
  * [10] J. Harasimczuk, N. E. Maliszewski, A. Olejniczak-Serowiec, and A. Tarnowski, “Are longer advertising slogans more dangerous? the influence of the length of ad slogans on drivers’ attention and motor behavior,” _Current Psychology_ , vol. 40, pp. 429–441, 2021.
  * [11] L. Meuleners, P. Roberts, and M. Fraser, “Identifying the distracting aspects of electronic advertising billboards: A driving simulation study,” _Accident Analysis & Prevention_, vol. 145, p. 105710, 2020.
  * [12] N. Maliszewski, A. Olejniczak-Serowiec, and J. Harasimczuk, “Influence of sexual appeal in roadside advertising on drivers’ attention and driving behavior,” _PloS one_ , vol. 14, no. 5, p. e0216919, 2019.
  * [13] A. Tarnowski, A. Olejniczak-Serowiec, and A. Marszalec, “Roadside advertising and the distraction of driver’s attention,” in _MATEC Web of Conferences_ , vol. 122. EDP Sciences, 2017, p. 03010.
  * [14] H. Marciano _et al._ , “The effect of billboard design specifications on driving: a pilot study,” _Accident Analysis & Prevention_, vol. 104, pp. 174–184, 2017.
  * [15] M. Chan and A. Singhal, “The emotional side of cognitive distraction: Implications for road safety,” _Accident Analysis & Prevention_, vol. 50, pp. 147–154, 2013.
  * [16] G. Jocher, A. Chaurasia, and J. Qiu, “YOLO by Ultralytics,” Jan. 2023. [Online]. Available: https://github.com/ultralytics/ultralytics
  * [17] M. Broström, “Real-time multi-object, segmentation and pose tracking using Yolov8 with DeepOCSORT and LightMBN.” [Online]. Available: https://github.com/mikel-brostrom/yolov8_tracking
  * [18] G. Neuhold, T. Ollmann, S. R. Bulò, and P. Kontschieder, “The mapillary vistas dataset for semantic understanding of street scenes,” _2017 IEEE International Conference on Computer Vision (ICCV)_ , pp. 5000–5009, 2017.
  * [19] J. Cao, J. Pang, X. Weng, R. Khirodkar, and K. Kitani, “Observation-centric sort: Rethinking sort for robust multi-object tracking,” in _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , 2023, pp. 9686–9696.
  * [20] G. Maggiolino, A. Ahmad, J. Cao, and K. Kitani, “Deep oc-sort: Multi-pedestrian tracking by adaptive re-identification,” _arXiv preprint arXiv:2302.11813_ , 2023.
  * [21] Y. Du, Z. Zhao, Y. Song, Y. Zhao, F. Su, T. Gong, and H. Meng, “Strongsort: Make deepsort great again,” _IEEE Transactions on Multimedia_ , pp. 1–14, 2023.
  * [22] N. Aharon, R. Orfaig, and B.-Z. Bobrovsky, “Bot-sort: Robust associations multi-pedestrian tracking,” _arXiv preprint arXiv:2206.14651_ , 2022.
  * [23] Y. Zhang, P. Sun, Y. Jiang, D. Yu, F. Weng, Z. Yuan, P. Luo, W. Liu, and X. Wang, “Bytetrack: Multi-object tracking by associating every detection box,” in _Computer Vision–ECCV 2022: 17th European Conference, Tel Aviv, Israel, October 23–27, 2022, Proceedings, Part XXII_. Springer, 2022, pp. 1–21.
  * [24] A. Bewley, Z. Ge, L. Ott, F. Ramos, and B. Upcroft, “Simple online and realtime tracking,” in _2016 IEEE International Conference on Image Processing (ICIP)_ , 2016, pp. 3464–3468.
  * [25] J. Luiten, A. Osep, P. Dendorfer, P. Torr, A. Geiger, L. Leal-Taixé, and B. Leibe, “Hota: A higher order metric for evaluating multi-object tracking,” _International journal of computer vision_ , vol. 129, pp. 548–578, 2021.
  * [26] G. Bellitto, F. Proietto Salanitri, S. Palazzo, F. Rundo, D. Giordano, and C. Spampinato, “Hierarchical domain-adapted feature learning for video saliency prediction,” _International Journal of Computer Vision_ , vol. 129, pp. 3216–3232, 2021.
  * [27] Z. Bylinskii, T. Judd, A. Oliva, A. Torralba, and F. Durand, “What do different evaluation metrics tell us about saliency models?” _IEEE transactions on pattern analysis and machine intelligence_ , vol. 41, no. 3, pp. 740–757, 2018.