# Low-rank approximated Kalman–Bucy filters using Oja’s principal component
flow for linear time-invariant systems

Daiki Tsuzuki1 and Kentaro Ohki1 *This work was supported by JSPS KAKENHI
Grant Number JP19K03619 and 23H01432.1Department of Applied Mathematics and
Physics, Graduate School of Informatics, Kyoto University, Yoshida-Honmachi,
Sakyo-ku, Kyoto 606-8501, Japan<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

The Kalman–Bucy filter is widely used in various applications. However, the
filter becomes computationally complex under large-scale systems. To address
this problem, a low-rank approximated Kalman–Bucy filter consisting of Oja’s
principal component flow and a low-dimensional Riccati differential equation
was proposed. However, the estimation error was established only for linear
time-invariant systems with a symmetric system matrix. This study removes
restrictions on the symmetricity of the system matrix and reveals the
equilibrium points of the Oja flow and their stability for general real square
matrices. In addition, the attraction domain for a set of stable equilibrium
points is estimated. Based on these results, we demonstrate that the low-rank
approximated Kalman–Bucy filter has a bounded estimation error covariance
matrix when the system is controllable and observable.

##  I Introduction

The Kalman/Kalman–Bucy filter [1, 2] is widely used in various fields such as
aerospace engineering [3, 4], weather forecasting [5, 6, 7], and quality
management [8, 9], for recursive state estimation from noisy and partial
observations. However, the computational burden of the Riccati equation
renders using the filter in high-dimensional systems limits the application of
the filter. Various methods, including model reduction [10, 11] and
approximations of Riccati differential equations [12, 13, 14, 15], remain
challenging.

Bonnabel et al. [14] proposed a low-rank approximated Kalman–Bucy (LRKB)
filter consisting of Oja’s principal component flow (Oja flow) [16] and low-
dimensional Riccati differential equations. In this filter, a low-rank
approximated covariance matrix is used to reduce the computational cost of the
Riccati equation. Although applicable to linear time-varying and nonlinear
systems as an extended Kalman–Bucy filter, no conditions for bounded
estimation error were provided. Besides, the estimation error can be quite
worse when the system satisfies a certain condition. Yamada et al. [15]
subsequently modified the low-rank filter to ensure a bounded estimation error
if a linear time-invariant (LTI) system is detectable and has a symmetric
system $A$ matrix. The restriction to symmetric matrices results from analysis
and use of the Oja flow to extract the principal or minor components that
focus on positive definite matrices [17, 18, 19, 20, 21]. Furthermore, the
equilibrium points and convergence properties have yet to be resolved for the
general $A$ matrix [22, 23]. Therefore, an analysis of the Oja flow is
critical for ensuring a bounded estimate error covariance using the modified
LRKB filter for a general $A$ matrix.

In this study, we removed the restriction on the symmetric system $A$ matrix
to ensure that the modified LRKB filter proposed in [15] exhibits satisfactory
performance. The contributions of this study are as follows:

  1. 1.

We show the explicit equilibrium points of Oja flow for general real square
matrices (Prop. 2) and their local stability (Theorem 1). Although disjoint
equilibrium sets exists, only one of them is stable. We also show a sufficient
condition for the domain of attraction of the stable equilibrium set (Theorem
2).

  2. 2.

The conditions for the bounded error covariance matrix of the LRKB filter are
derived (Theorem 3).

The remainder of this paper is organized as follows. We review the Kalman–Bucy
and LRKB filters in Section II. We show our contribution to the analysis of
the Oja flow in Section III. The equilibrium points of the Oja flow and their
stability are presented, and the domain of attraction for the stable
equilibrium points is derived. Based on the results of the Oja flow, Section
IV provides the conditions for the estimation error covariance matrix to be
remain finite. We present the conclusion of this paper in Section V.

#### Notation

The sets of real and complex numbers are denoted by $\mathbb{R}$ and
$\mathbb{C}$, respectively. The sets of $n\times m$ real and complex matrices
are described by $\mathbb{R}^{n\times m}$ and $\mathbb{C}^{n\times m}$,
respectively. The $n\times n$ identity matrix is denoted by $I_{n}$ and
$n\times m$ zero matrix is denoted by $O_{n,m}$. For $n\times m$ matrix $A$,
$A^{\top}$ and $A^{\dagger}$ are the transposed and Hermitian conjugates of
$A$, respectively. For a real symmetric matrix $A$, $A>0(A\geq 0)$ indicates
that $A$ is positive (semi-)definite. For a square matrix
$A\in\mathbb{C}^{n\times n}$, its eigenvalues are denoted by
$\lambda_{i}(A)\quad(\mathrm{Re}(\lambda_{1}(A))\geq\dots\geq\mathrm{Re}(\lambda_{n}(A)))$,
and the corresponding eigenvectors with norm 1 are denoted by
$\bm{\psi}_{i}(A)$, $i=1,\dots,n$. If the argument is clear, we can simply
write $\lambda_{i}$ and $\bm{\psi}_{i}$. For the degenerated eigenvalues, we
always use independent eigenvectors. The Stiefel manifold is denoted as
$\mathrm{St}(r,n):=\\{X\in{\mathbb{R}^{n\times r}}|{X}^{\top}X=I_{r}\\}$.

##  II Preliminary: Kalman–Bucy filter and its low-rank approximation

###  II-A Kalman–Bucy filter

As stated in Introduction, the introduced LRKB filter can be used in linear
time-varying or nonlinear systems. However, in this study we focus only on LTI
systems .

$\displaystyle\frac{dx(t)}{dt}=$ $\displaystyle Ax(t)+Gw(t),$ (1)
$\displaystyle y(t)=$ $\displaystyle Cx(t)+Hv(t),$ (2)

where $x(t)\in\mathbb{R}^{n}$ is the state variable, $w(t)\in\mathbb{R}^{n}$
is the process noise, $y(t)\in\mathbb{R}^{p}$ is the observation, and
$v(t)\in\mathbb{R}^{p}$ is the observation noise. We assume that $w(t)$ and
$v(t)$ are independent white Gaussian noises with a mean of $0$ and identity
covariance matrices, respectively. System matrices $A,G,C$, and $H$ are
appropriately defined real matrices, and $HH^{\top}$ is a positive-definite
matrix. Throughout this study, we assume that $r\in\\{1,\dots,n-1\\}$ exists
such that $\mathrm{Re}(\lambda_{r}(A))>\mathrm{Re}(\lambda_{r+1}(A))$. The
Kalman–Bucy filter for the system (1)-(2) is given as follows:

$\displaystyle\frac{d\hat{x}(t)}{dt}=$
$\displaystyle(A-\hat{P}(t){C}^{\top}{(HH^{\top})}^{-1}C)\hat{x}(t)$
$\displaystyle+\hat{P}(t){C}^{\top}{(HH^{\top})}^{-1}y(t),$ (3)
$\displaystyle\frac{d\hat{P}(t)}{dt}=$ $\displaystyle
A\hat{P}(t)+\hat{P}(t){A}^{\top}+GG^{\top}$
$\displaystyle-\hat{P}(t){C}^{\top}{(HH^{\top})}^{-1}C\hat{P}(t)$ (4)

with $\hat{x}(0)=\mathbb{E}[x(0)]$ and
$\hat{P}(0)=\mathbb{E}[(x(0)-\hat{x}(0))(x(0)-\hat{x}(0))^{\top}]$, where
$\hat{x}(t)\in\mathbb{R}^{n}$ is the estimate of $x(t)$ and its estimation
error covariance matrix is $\hat{P}(t)\in\mathbb{R}^{n\times n}$. The
following proposition provides the steady-state solution to Equation (4).

######  Proposition 1 ([24]).

If the system $(A,G,C)$ is controllable (stabilizable) and observable
(detectable), the following algebraic Riccati equation

$\displaystyle
A\hat{P}_{s}+\hat{P}_{s}A^{\top}+GG^{\top}-\hat{P}_{s}C^{\top}{(HH^{\top})}^{-1}C\hat{P}_{s}=O_{n,n}$

has a unique positive (semi-)definite solution $\hat{P}_{s}$ and
$A-\hat{P}_{s}C^{\top}{(HH^{\top})}^{-1}C$ is stable. Furthermore, the
solution $\hat{P}(t)$ to Equation (4) converges to $\hat{P}_{s}$.

###  II-B Low-Rank approximated Kalman–Bucy filter

Implementing the Kalman–Bucy filter requires real-time numerical computation
or calculation of the steady-state solution of the algebraic Riccati equation
offline; both are difficult if the system’s dimensions are high. To address
this problem, we considered a modified LRKB filter [15] hereafter referred to
as the LRKB filter. This concept, originating from [14], can substitute a
positive (semi-)definite matrix $\hat{P}(t)$ with
$\tilde{P}(t)=U(t)\tilde{R}(t)U(t)^{\top}$, where $U(t)\in\mathrm{St}(r,n)$
and $\tilde{R}(t)>0$. The time evolution equations for $U(t)$ and
$\tilde{R}(t)$ are as follows:

$\displaystyle\varepsilon\frac{dU(t)}{dt}=$
$\displaystyle(I_{n}-U(t){U}(t)^{\top})AU(t),\ U(0)^{\top}U(0)=I_{r},$ (5)
$\displaystyle\frac{d\tilde{R}(t)}{dt}=$ $\displaystyle
A_{U(t)}\tilde{R}(t)+\tilde{R}(t)A_{U(t)}^{\top}+G_{U(t)}G_{U(t)}^{\top}$
$\displaystyle-\tilde{R}(t)C_{U(t)}^{\top}{(HH^{\top})}^{-1}C_{U(t)}\tilde{R}(t),\
\tilde{R}(0)>0,$ (6)

where $A_{U}\coloneqq U^{\top}AU$, $C_{U}\coloneqq CU$ and $G_{U}\coloneqq
U^{\top}G$. $\varepsilon\in(0,1]$ is a small parameter that is used to adjust
the convergence speed of Equation (5). Equation (5) is known as the Oja flow
[17] and usually takes $\varepsilon=1$. Provided that it follows an Oja flow,
$U(t)\in\mathrm{St}(r,n)$ holds at any time [25]. Equation (6) is an
$r$-dimensional Riccati differential equation. The differential equation for
the estimated value $\tilde{x}(t)$ is as follows:

$\displaystyle\frac{d\tilde{x}(t)}{dt}=$
$\displaystyle(A-\tilde{P}(t){C}^{\top}{(HH^{\top})}^{-1}C)\tilde{x}(t)$
$\displaystyle+{\tilde{P}}(t){C}^{\top}{(HH^{\top})}^{-1}y(t).$ (7)

The differential equation for the error covariance matrix
$\tilde{V}(t)\coloneqq\mathbb{E}[(\tilde{x}(t)-x(t)){(\tilde{x}(t)-x(t))}^{\top}]$
is as follows:

$\displaystyle\frac{d\tilde{V}(t)}{dt}=$
$\displaystyle(A-\tilde{P}(t){C}^{\top}{(HH^{\top})}^{-1}C)\tilde{V}(t)$
$\displaystyle+\tilde{V}(t){(A-\tilde{P}(t){C}^{\top}{(HH^{\top})}^{-1}C)}^{\top}$
$\displaystyle+GG^{\top}+\tilde{P}(t){C}^{\top}{(HH^{\top})}^{-1}C\tilde{P}(t).$
(8)

Bonnabel et al. [14] used $\lambda_{n}(GG^{\top})I_{r}$ instead of
$G_{U(t)}G_{U(t)}^{\top}$ in Equation (6). However, when $GG^{\top}$ is
singular, implying $\lambda_{n}(GG^{\top})=0$, the estimation error of
Bonnabel’s LRKB filter can become large if $(A,G)$ is controllable. To
overcome this problem, Yamada et al. [15] modified Bonnabel’s filter,
providing conditions for bounded estimation error for LTI systems with
symmetric $A$. The analysis was restrictive; therefore, this study presents
the general results of [15].

##  III Equilibrium points and convergence of the Oja flow

First, we detail all the equilibrium points of Oja flow (5). In addition to
their local stability, this section presents the domains of attraction of the
stable equilibrium points.

###  III-A Equilibrium points of the Oja flow

We present several lemmas to examine the structure of the equilibrium points
of the Oja flow (5). At any time $t$, we set $U(t)=\Psi
K(t)\in\mathrm{St}(r,n)$, where
$\Psi:=[\bm{\psi}_{1}(A),\dots,\bm{\psi}_{n}(A)]$ and
$K(t)\in\mathbb{C}^{n\times r}$. The differential equation for $K(t)$ is then
given by the following expression:

$\displaystyle\varepsilon\frac{dK(t)}{dt}=\Lambda
K(t)-K(t){K(t)}^{\dagger}{\Psi}^{\dagger}\Psi\Lambda K(t),$ (9)

where $\Lambda={\Psi}^{-1}A\Psi$ is a diagonal matrix or Jordan normal form.
Next, we investigate the equilibrium points of Eq. (9). The following lemma
holds true:

######  Lemma 1.

Let $K\in\mathbb{C}^{n\times r}$ be $K=\Psi^{-1}U$, where
$U\in\mathrm{St}(r,n)$. Then, $K$ has the form

$\displaystyle K=P\begin{bmatrix}K_{r}\\\ O_{n-r,n}\end{bmatrix},$ (10)

where $P\in\mathbb{C}^{n\times n}$ is a unitary matrix, and
$K_{r}\in\mathbb{C}^{r\times r}$ is a regular matrix.

###### Proof.

Because the rank of $U\in\mathrm{St}(r,n)$ is $r$, the rank of $K$ is
$\mathrm{rank}(K)=\mathrm{rank}({\Psi}^{-1}U)=r$. Using singular value
decomposition, this statement holds. ∎

For a unitary matrix $P\in\mathbb{C}^{n\times n}$, let
$\Sigma^{P}:=\begin{bmatrix}\Sigma_{11}^{P}&\Sigma_{12}^{P}\\\
\Sigma_{21}^{P}&\Sigma_{22}^{P}\end{bmatrix}:={P}^{\dagger}{\Psi}^{\dagger}\Psi
P$, where $\Sigma_{11}^{P}\in\mathbb{C}^{r\times r}$,
$\Sigma_{12}^{P}=(\Sigma_{21}^{P})^{\dagger}\in\mathbb{C}^{(n-r)\times r}$,
and $\Sigma_{22}^{P}\in\mathbb{C}^{(n-r)\times(n-r)}$. If $P=I_{n}$, the
superscript is omitted. Then, the following lemma holds for $K_{r}$.

######  Lemma 2.

Let $K_{r}$ be of Lemma 1. Then, $K_{r}$ is represented by a unitary matrix
$W_{0}\in\mathbb{C}^{r\times r}$ such that

$\displaystyle K_{r}=(\Sigma_{11}^{P})^{-1/2}W_{0}.$

For any orthogonal matrix $W\in\mathbb{R}^{r\times r}$, the following equation
also satisfies (10):

$\displaystyle K_{r}=(\Sigma_{11}^{P})^{-1/2}W_{0}W.$

###### Proof.

From polar decomposition, a unitary matrix $W_{0}$ exists such that
$K_{r}=(K_{r}K_{r}^{\dagger})^{1/2}W_{0}$. Because
$U^{\top}U=U^{\dagger}U=I_{r}$, then
$\begin{bmatrix}K_{r}^{\dagger}&&O_{n,n-r}\end{bmatrix}{P}^{\dagger}{\Psi}^{\dagger}\Psi
P\begin{bmatrix}{K_{r}}\\\ O_{n,n-r}\end{bmatrix}=I_{r}$. Therefore,
$K_{r}^{\dagger}\Sigma_{11}^{P}K_{r}=I_{r}$ and
$K_{r}=(\Sigma_{11}^{P})^{-1/2}W_{0}$. Because
$U\in\mathrm{St}(r,n)\subset\mathbb{R}^{n\times r}$,
$K_{r}=(\Sigma_{11}^{P})^{-1/2}W_{0}W$ for any orthogonal matrix
$W\in\mathbb{R}^{r\times r}$ always satisfies $U^{\top}U=I_{r}$. Therefore,
the statement holds. ∎

Using Lemmas 1 and 2, the following proposition holds for the equilibrium
points of the Oja flow (5) with a real square matrix $A$.

######  Proposition 2.

The equilibrium point is given as $\bar{U}_{P}=\Psi P\begin{bmatrix}K_{r}\\\
O_{n-r,r}\end{bmatrix}\in\mathrm{St}(r,n)$, where $P\in\mathbb{R}^{n\times n}$
is a permutation matrix and $K_{r}\in\mathbb{C}^{r\times r}$ is of the form in
Lemma 2. Furthermore, any $\bar{U}_{P}^{\prime}$ of
$\mathcal{U}_{P}:=\\{\bar{U}_{P}W\ |\ W\in\mathbb{R}^{r\times r},\
W^{\top}W=I_{r}\\}$ is also an equilibrium point.

We write $\bar{U}$ and $\mathcal{U}$ if $P=I_{n}$.

###### Proof.

From Lemma 1, $\bar{U}_{P}=\Psi P\begin{bmatrix}K_{r}\\\
O_{n,n-r}\end{bmatrix}$, Let $\Lambda_{P}:=P^{\dagger}\Lambda P$. By using
Equation (5) and Lemma 2, we have the following equation:

$\displaystyle 0=$
$\displaystyle(I_{n}-\bar{U}_{P}{\bar{U}_{P}}^{\top})A\bar{U}_{P}$
$\displaystyle=$ $\displaystyle\Psi
P\left(\Lambda_{P}-\begin{bmatrix}K_{r}K_{r}^{\dagger}&O_{r,n-r}\\\
O_{n-r,r}&O_{n-r,n-r}\end{bmatrix}\Sigma^{P}\Lambda_{P}\right)\begin{bmatrix}K_{r}\\\
O_{n,n-r}\end{bmatrix}.$

Because $\Psi P$ is a nonsingular matrix, the following equation holds:

$\displaystyle 0$
$\displaystyle=\left(\Lambda_{P}-\begin{bmatrix}I_{r}&({\Sigma_{11}^{P}})^{-1}{\Sigma_{12}^{P}}\\\
O_{n-r,r}&O_{n-r,n-r}\end{bmatrix}\Lambda_{P}\right)\begin{bmatrix}K_{r}\\\
O_{n,n-r}\end{bmatrix}.$

The equality holds if $\Lambda_{P}$ is a block diagonal matrix. To satisfy the
block diagonality of $\Lambda_{P}$, a unitary matrix $P$ should be the product
of permutation matrix $P_{\rm per}$ and block diagonal unitary matrix $P_{\rm
BD}=\mathrm{block}\mbox{-}\mathrm{diag}[P_{r},P_{\perp}]$, where
$P_{r}\in\mathbb{C}^{r\times r}$ and
$P_{\perp}\in\mathbb{C}^{(n-r)\times(n-r)}$ are unitary matrices. We can
rewrite $P_{r}K_{r}$ with $K_{r}$; thus, matrix $P$ is a permutation matrix
without loss of generality. From Lemma 2, if $\bar{U}_{P}$ is also an
equilibrium point, then for any orthogonal matrix $W\in\mathbb{R}^{r\times
r}$, $\bar{U}_{P}W$ is also an equilibrium point. The proof is completed. ∎

######  Remark 1.

In [14], a stable equilibrium set for a general $A$ matrix was characterized
as the symmetric part of $A$.

$\displaystyle\\{[\bm{\psi}_{1}^{\rm s},\dots,\bm{\psi}_{r}^{\rm
s}]F_{r}\in\mathrm{St}(r,n)\ |\ F_{r}\in\mathbb{R}^{r\times r},\
F_{r}^{\top}F_{r}=I_{r}\\},$

where $\bm{\psi}^{\rm s}_{i}:=\bm{\psi}_{i}(A+A^{\top})\in\mathbb{R}^{n}$,
$i=1,\dots,r$, are the first $r$ unit eigenvectors of the eigenvalues
$\lambda_{i}(A+A^{\top})$, $i=1,\dots,n$. Proposition 2 provides an
alternative and direct characterization of the equilibrium sets of the Oja
flow. Proposition 2 does not deny the existence of other invariant sets such
as limit cycles.

Any $\bar{U}_{P}\in\mathcal{U}_{P}$ retains certain eigenvalues of $A$.

######  Proposition 3.

Let $P$ be any permutation matrix and $\bar{U}_{P}\in\mathcal{U}_{P}$. Then
$\lambda_{i}({\bar{U}_{P}}^{\top}A\bar{U}_{P})=\lambda_{\mathcal{I}_{P}(i)}(A)$
for $i\in\\{1,\dots,r\\}$, where $\mathcal{I}_{P}$ is a permutation related to
$P$ and preserves the order
$\mathrm{Re}(\lambda_{\mathcal{I}_{P}(1)}(A))\geq\dots\geq\mathrm{Re}(\lambda_{\mathcal{I}_{P}(r)}(A))$.

###### Proof.

Let ${\Psi}_{rp}\coloneqq\Psi P\begin{bmatrix}I_{r}\\\
O_{n-r,r}\end{bmatrix}\in\mathbb{C}^{n\times r}$,
$\Lambda_{rp}\coloneqq\begin{bmatrix}I_{r}&&O_{r,n-r}\end{bmatrix}P^{\top}\Lambda
P\begin{bmatrix}I_{r}\\\ O_{n-r,r}\end{bmatrix}\in\mathbb{C}^{r\times r}$.
Then,
$\bar{U}_{P}^{\top}\bar{U}_{P}=K_{r}^{\dagger}{\Psi}_{rp}^{\dagger}{\Psi}_{rp}K_{r}=I_{r}$.
Because $K_{r}$ is a nonsingular matrix, we have the following equation:

$\displaystyle\lambda_{i}(\bar{U}_{P}^{\top}A\bar{U}_{P})$
$\displaystyle=\lambda_{i}(K_{r}^{\dagger}{\Psi}_{rp}^{\dagger}{\Psi}_{rp}K_{r}K_{r}^{-1}\Lambda_{rp}K_{r})$
$\displaystyle=\lambda_{i}(K_{r}^{-1}\Lambda_{rp}K_{r}).$

Because the similarity transformation preserves the eigenvalues,
$\lambda_{i}({\bar{U}_{P}}^{\top}A\bar{U}_{P})=\lambda_{\mathcal{I}_{P}(i)}(A)$
holds. ∎

The preservation of the eigenvalues of $A$ is crucial for analyzing the
controllability and observability of the reduced system matrices in Props. 5
and 6.

###  III-B Convergence and domain of attraction of the Oja flow

In this section, we analyze the local stability of each $\mathcal{U}_{P}$
defined in Prop. 2 and the domain of attraction of the Oja flow. Additionally,
we set $K(t)=\begin{bmatrix}K_{r}(t)\\\
K_{\perp}(t)\end{bmatrix}=\Psi^{-1}U(t),$ where
$K_{r}(t)\in\mathbb{C}^{r\times r},K_{\perp}(t)\in\mathbb{C}^{(n-r)\times r}$.

######  Theorem 1.

The set $\mathcal{U}$ is locally asymptotically stable. Any equilibrium set
$\mathcal{U}_{P}\neq\mathcal{U}$ is unstable.

###### Proof.

We consider the perturbations from the equilibrium point $\Psi
P\begin{bmatrix}K_{r}\\\ O_{n-r,r}\end{bmatrix}$ given by $\Psi
P\begin{bmatrix}\delta K_{r}(t)\\\ \delta K_{\perp}(t)\end{bmatrix}$. We
define $\mathrm{block}\mbox{-}\mathrm{diag}[\Lambda_{rp},\Lambda_{\perp
p}]:=P^{\top}\Lambda P$. Because $\lim_{t\to\infty}\delta
K_{\perp}(t)=O_{n-r,r}$ is sufficient for this statement, we focus only on
$\frac{d}{dt}\delta K_{\perp}(t)$. Ignoring higher-order terms and using
$\varepsilon=1$ for simplicity, the differential equation for $\delta
K_{\perp}(t)$ is as follows:

$\displaystyle\frac{d\delta K_{\perp}(t)}{dt}=$ $\displaystyle\Lambda_{\perp
p}\delta K_{\perp}(t)-\delta K_{\perp}(t)K_{r}^{-1}\Lambda_{rp}K_{r}.$

The vectorization of the matrices yields

$\displaystyle\frac{d\mathrm{vec}(\delta K_{\perp}(t))}{dt}$ $\displaystyle=$
$\displaystyle(I_{r}\otimes\Lambda_{\perp
p}-(K_{r}^{\top}\Lambda_{rp}^{\top}K_{r}^{-\top})\otimes
I_{n-r})\mathrm{vec}(\delta K_{\perp}(t)),$

where $\otimes$ denotes the Kronecker product. Because
$\mathrm{Re}(\lambda_{r}(A))>\mathrm{Re}(\lambda_{r+1}(A))$ and $P=I_{n}$,
$\delta K_{\perp}(t)$ converges zero matrix and consequently, $K_{r}+\delta
K_{r}(t)$ converges to $K_{r}W$, where $W\in\mathbb{R}^{r\times r}$ is an
orthogonal matrix. Therefore, $\mathcal{U}$ is asymptotically stable. Note
that if the permutation $P$ is of the form
$P=\mathrm{block}\mbox{-}\mathrm{diag}[P_{r},\ P_{\perp}]$, where
$P_{r}\in\mathbb{R}^{r\times r}$ and
$P_{\perp}\in\mathbb{R}^{(n-r)\times(n-r)}$, then
$\mathcal{U}=\mathcal{U}_{P}$. However, any element in
$\mathcal{U}_{P}\neq\mathcal{U}$ leaves $\mathcal{U}_{P}$. Therefore, any set
$\mathcal{U}_{P}\neq\mathcal{U}$ is unstable. ∎

Theorem 1 does not show the explicit domain of attraction of $\mathcal{U}$.
However, if we consider $(n,r)=(2,1)$, the global behavior of the Oja flow is
given as follows:

######  Proposition 4.

Suppose $(n,r)=(2,1)$. If
$K(0)\neq\begin{bmatrix}0,&\pm\Sigma_{22}^{-1/2}\end{bmatrix}^{\top}$, then
the solution $\Psi K(t)$ of (9) converges to
$\mathcal{U}=\\{\pm\Sigma_{11}^{-1/2}\psi_{1}(A)\\}$.

###### Proof.

For simplicity, we put $\varepsilon=1$. For $A\in\mathbb{R}^{2\times 2}$,
$\mathrm{Re}(\lambda_{1}(A))>\mathrm{Re}(\lambda_{2}(A))$ implies
$\lambda_{i}$ are real. Therefore, $\Sigma=\Sigma^{\top}\in\mathbb{R}^{2\times
2}$. Here, $K(t)$ moves along an ellipse $\\{K\in\mathbb{R}^{2}\ |\
K^{\top}\Sigma K=1\\}$. From (9) and $K(t)^{\top}\Sigma
K(t)=\Sigma_{11}K_{r}(t)^{2}+\Sigma_{22}K_{\perp}(t)^{2}+2\Sigma_{12}K_{r}(t)K_{\perp}(t)=1$,

$\displaystyle\frac{d}{dt}K(t)=$
$\displaystyle\frac{\delta\lambda}{2}\left\\{\begin{bmatrix}0\\\
-2K_{\perp}(t)\end{bmatrix}+\gamma(K_{r}(t),K_{\perp}(t))\begin{bmatrix}K_{r}(t)\\\
K_{\perp}(t)\end{bmatrix}\right\\},$ (11)

where $\delta\lambda:=\lambda_{1}(A)-\lambda_{2}(A)$ and
$\gamma(K_{r},K_{\perp}):=(1+\Sigma_{2}K_{\perp}^{2}-\Sigma_{1}K_{r}^{2})$.
Equation (11) contains four equilibrium points:
$(\bar{K}_{r},\bar{K}_{\perp})=(\pm\Sigma_{1}^{-1/2},0)$ and
$(0,\pm\Sigma_{2}^{-1/2})$. Now we consider the following four cases: (i) If
$\gamma(K_{r}(t),K_{\perp}(t))>2$, then $|K_{r}(t)|$ and $|K_{\perp}(t)|$
increase. (ii) If $\gamma(K_{r}(t),K_{\perp}(t))\in(0,2]$, then $|K_{r}(t)|$
increases and $|K_{\perp}(t)|$ decreases. (iii) If
$\gamma(K_{r}(t),K_{\perp}(t))=0$, then $K_{r}(t)$ does not change and
$|K_{\perp}(t)|$ decreases or does not change if $K_{\perp}(t)=0$. (iv) If
$\gamma(K_{r}(t),K_{\perp}(t))<0$, then $|K_{r}(t)|$ and $|K_{\perp}(t)|$
decrease. In (i), the rate of increase of $|K_{r}(t)|$ is larger than that of
$K_{\perp}(t)$. Then, Equation (11) enters case (ii) and becomes case (iii).
After case (iii), Equation (11) does not return to case (ii); thus,
$|K_{\perp}(t)|$ converges to $0$. From the constraint $K(t)^{\top}\Sigma
K(t)=1$, this phenomenon implies that $K_{r}(t)$ converges to
$\pm\Sigma_{11}^{-1/2}$.

∎

Unlike Proposition 4, estimating the exact domain of attraction for a general
$(n,r)$ is difficult. We provide sufficient conditions for the domain of
attraction. Before proceeding, we introduce the following lemmas:

######  Lemma 3.

For any matrix $X,Y\in\mathbb{C}^{m\times n}$ and any real number $\alpha>0$,
the following inequalities hold:

$\displaystyle\pm(X^{\dagger}Y+Y^{\dagger}X)\leq\alpha
X^{\dagger}X+\frac{1}{\alpha}Y^{\dagger}Y.$

###### Proof.

For any $X,Y\in\mathbb{C}^{m\times n}$, we obtain
$\begin{bmatrix}X&Y\end{bmatrix}^{\dagger}\begin{bmatrix}X&Y\end{bmatrix}\geq
0$. By multiplying both sides of the inequality by $[\sqrt{\alpha}I_{n},\
\pm\frac{1}{\sqrt{\alpha}}I_{n}]$ and its transposition from the left and
right sides, we obtain Lemma 3. ∎

######  Lemma 4.

For any $A\in\mathbb{R}^{n\times n}$, a unitary matrix
$S\in\mathbb{C}^{n\times n}$ and upper triangle matrix
$Q\in\mathbb{C}^{n\times n}$ exist such that $S^{\dagger}AS$ becomes an upper
triangle and satisfies $S=\Psi Q$.

###### Proof.

From the QR decomposition of $\Psi$, $\Psi=ST$, where $T$ is the upper
triangle and regular matrix. Because $Q=T^{-1}$ is also an upper triangle,
this statement holds. ∎

If the matrix $A$ is decomposed using the unitary matrix $S$ satisfying Lemma
4, we can put $S^{\dagger}AS=\begin{bmatrix}L_{11}&&L_{12}\\\
O_{n-r,r}&&L_{22}\end{bmatrix}$, where $L_{11}$ is an upper triangular matrix
with eigenvalues $\lambda_{1}(A),\dots,\lambda_{r}(A)$, $L_{22}$ is an upper
triangular matrix with eigenvalues $\lambda_{r+1}(A),\dots,\lambda_{n}(A)$,
and $L_{12}$ is in $\mathbb{C}^{r\times(n-r)}$. Notice that
$\mathrm{St}(r,n)=\\{S\begin{bmatrix}F_{1}\\\
F_{2}\end{bmatrix}\in\mathbb{R}^{n\times r}\mid F_{1}\in\mathbb{C}^{r\times
r},F_{2}\in\mathbb{C}^{(n-r)\times r},\
F_{1}^{\dagger}F_{1}+F_{2}^{\dagger}F_{2}=I_{r}\\}$.

######  Theorem 2.

If
${\lambda}^{l_{1}}_{r}\coloneqq\lambda_{r}(L_{11}+L_{11}^{\dagger})>{\lambda}^{l_{2}}_{1}\coloneqq\lambda_{1}(L_{22}+L_{22}^{\dagger})$,
and $U(0)\in\mathcal{V}_{\beta}:=\\{S\begin{bmatrix}F_{1}\\\
F_{2}\end{bmatrix}\in\mathrm{St}(r,n)\mid F_{1}\in\mathbb{C}^{r\times
r},F_{2}\in\mathbb{C}^{(n-r)\times r},\ 0\leq F_{2}^{\dagger}F_{2}<\beta
I_{r}\\}$, then the solution to Equation (5) converges to an element of
$\mathcal{U}$. Here,
$\beta=(1+4\ell_{\max}/({\lambda}^{l_{1}}_{r}-{\lambda}^{l_{2}}_{1})^{2})^{-1}$
and $\ell_{\max}$ is the maximum eigenvalue of $L_{12}^{\dagger}L_{12}$.

###### Proof.

By using Lemma 4, $U(t)$ can be expressed as $U(t)=\Psi
Q\begin{bmatrix}F_{1}(t)^{\top}&F_{2}(t)^{\top}\end{bmatrix}^{\top}$. We show
the initial domain is such that $\lim_{t\to\infty}F_{2}(t)=O_{(n-r),r}$,
implying $U(t)$ converges to $\mathcal{U}$. The differential equation of
$\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t)]$ is

$\displaystyle\frac{d}{dt}\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t)]$
$\displaystyle=$
$\displaystyle\mathrm{Tr}[{F_{2}(t)}^{\dagger}(L_{22}+L_{22}^{\dagger})F_{2}(t)\\{I_{r}-{F_{2}(t)}^{\dagger}F_{2}(t)\\}]$
$\displaystyle-\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t)F_{1}(t)^{\dagger}(L_{11}+L_{11}^{\dagger})F_{1}(t)]$
$\displaystyle-\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t)\\{{F_{1}(t)}^{\dagger}L_{12}{F_{2}(t)}+{F_{2}(t)}^{\dagger}L_{12}^{\dagger}{F_{1}(t)}\\}].$

By applying Lemma 3 with $\alpha\in(0,\delta\lambda)$, where
$\delta\lambda:={\lambda}^{l_{1}}_{r}-{\lambda}^{l_{2}}_{1}$, and
${F_{1}(t)}^{\dagger}{F_{1}(t)}=I_{r}-{F_{2}(t)}^{\dagger}F_{2}(t)$ to the
above equation, the following inequality holds:

$\displaystyle\frac{d}{dt}\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t)]$
$\displaystyle\leq$
$\displaystyle-\delta\lambda\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t)\\{I_{r}-{F_{2}(t)}^{\dagger}F_{2}(t)\\}]$
$\displaystyle+\frac{1}{\alpha}\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t){F_{2}(t)}^{\dagger}L_{12}^{\dagger}L_{12}{F_{2}(t)}]$
$\displaystyle+\alpha\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t)\\{I_{r}-{F_{2}(t)}^{\dagger}F_{2}(t)\\}]$
$\displaystyle\leq$
$\displaystyle-(\delta\lambda-\alpha)\mathrm{Tr}[{F_{2}(t)}^{\dagger}F_{2}(t)]$
$\displaystyle+\left(\delta\lambda-\alpha+\frac{\ell_{\max}}{\alpha}\right)\mathrm{Tr}[({F_{2}(t)}^{\dagger}F_{2}(t))^{2}].$
(12)

We consider the following Riccati differential equation:

$\displaystyle\frac{d}{dt}Z(t)=$
$\displaystyle-(\delta\lambda-\alpha)Z(t)+\left(\delta\lambda-\alpha+\frac{\ell_{\max}}{\alpha}\right)Z(t)^{2},$
(13)

where $Z(t)=Z(t)^{\dagger}\in\mathbb{C}^{r\times r}$ and
$Z(0)=F_{2}(0)^{\dagger}F_{2}(0)$. The solution of Equation (13) converges to
$O_{r,r}$ if
$Z(0)<\frac{\alpha(\delta\lambda-\alpha)}{\alpha(\delta\lambda-\alpha)+\ell_{\max}}I_{r}$.
Since $\mathrm{Tr}[F_{2}(t)^{\dagger}F_{2}(t)]\leq\mathrm{Tr}[Z(t)]$ for
$t\geq 0$, if
$F_{2}(0)^{\dagger}F_{2}(0)<\frac{\alpha(\delta\lambda-\alpha)}{\alpha(\delta\lambda-\alpha)+\ell_{\max}}I_{n-r}$,
then $F_{2}(t)$ converges to $O_{n-r,r}$. Thus, $U(t)$ converges to
$\mathcal{U}$.

Next, we estimate the maximum domain of attraction. Let

$\displaystyle
f(\alpha):=\frac{\alpha(\delta\lambda-\alpha)}{\alpha(\delta\lambda-\alpha)+\ell_{\max}},\
\alpha\in(0,\delta\lambda).$

Then, $f$ is a concave, positive function over $(0,\delta\lambda)$ and is
maximum when $\alpha_{\max}=\delta\lambda/2$. Thus,
$f(\alpha_{\max})=\frac{\delta\lambda^{2}}{\delta\lambda^{2}+4\ell_{\max}}=\beta$,
thereby completes the proof.

∎

By Bendixson’s inequality, ${\lambda}^{l_{1}}_{r}>{\lambda}^{l_{2}}_{1}$
implies $\mathrm{Re}(\lambda_{r}(A))>\mathrm{Re}(\lambda_{r+1}(A))$. From
numerical study as in [26], if
$\mathrm{Re}(\lambda_{r}(A))>\mathrm{Re}(\lambda_{r+1}(A))$ is satisfied, then
solution $U(t)$ to Eq. (5) can converge to $\mathcal{U}$ for almost initial
matrices. In the future, we will study an appropriate Lyapunov function to
prove this type of general convergence results. Theorem 2 is similar to
Proposition 4 if $\ell_{\max}=0$. This condition is satisfied when $A$ is a
normal matrix:

######  Corollary 1.

Let $A$ be the normal matrix. If $U(0)\in\mathcal{V}_{1}$, then the solution
to Eq. (5) converges to an element of $\mathcal{U}$.

###### Proof.

For a normal matrix $A\in\mathbb{R}^{n\times n}$, $S$ diagonalizes $A$.
Therefore, $\ell_{max}=0$. This statement follows from Theorem 2.

∎

The closure of $\mathcal{V}_{1}$ is equal to $\mathrm{St}(r,n)$; therefore,
for any normal matrix $A$, almost all the elements of $\mathrm{St}(r,n)$ are
in the domain of attraction. For example, almost all samples from the uniform
distribution of $\mathrm{St}(r,n)$ satisfy Corollary 1. Theorem 2.2.1 of [27]
can be used to generate such samples.

##  IV Boundedness of estimation error

To ensure a bounded estimation error for the LRKB filter, we establish the
conditions for the steady-state solution of Eq. (8) as in [15] that dealt with
symmetric $A$. We use $\varepsilon\to 0$ to simplify the discussion, and
consider only the steady-state element $\bar{U}\in\mathcal{U}$ discussed.
First, we begin the analysis by studying the controllability and observability
using the following equality for any permutation matrix
$P\in\mathbb{R}^{n\times n}$:

$\displaystyle\bar{U}_{P}\bar{U}_{P}^{\top}A\bar{U}_{P}=A\bar{U}_{P}\quad\forall\bar{U}_{P}\in\mathcal{U}_{P}.$
(14)

This equality is proven in the proof of Prop. 2.

######  Proposition 5.

We assume that $\bar{U}\in\mathcal{U}$. If $(C,A)$ is observable, then
$(C_{\bar{U}},A_{\bar{U}})$ is also observable.

###### Proof.

If $(C_{\bar{U}},A_{\bar{U}})$ is unobservable, then a nonzero vector
$\bm{v}\in{\mathbb{C}^{r}}$ and $\lambda_{i}$ exist such that

$\displaystyle\bar{U}^{\top}A\bar{U}\bm{v}=\lambda_{i}\bm{v},\quad
C\bar{U}\bm{v}=0.$

From Prop. 3, $\lambda_{k}=\lambda_{k}(A)$, $k=1,\dots,r$. From Eq. (14), we
have

$\displaystyle\bar{U}\bar{U}^{\top}A\bar{U}\bm{v}=A\bar{U}\bm{v}=\lambda_{i}\bar{U}\bm{v},\quad
C\bar{U}\bm{v}=0.$

Because the column vectors of $\bar{U}$ are linearly independent,
$\bar{U}\bm{v}\neq 0$. This result contradicts the observability of $(C,A)$.
Therefore, if $(C,A)$ is observable, then $(C_{\bar{U}},A_{\bar{U}})$ is also
observable. ∎

######  Proposition 6.

We assume that $\bar{U}\in\mathcal{U}$. If $(A,G)$ is controllable, then
$(A_{\bar{U}},G_{\bar{U}})$ is also controllable.

###### Proof.

If $(A,G)$ is controllable, then we have

$\displaystyle\mathrm{rank}\begin{bmatrix}G&A-\lambda_{i}I_{n}\end{bmatrix}=n$

for any $\lambda_{i}=\lambda_{i}(A)$, $i=1,\dots,n$. Next, we consider the
rank of the following matrix;

$\displaystyle\bar{U}^{\top}\begin{bmatrix}G&A-\lambda_{i}I_{n}\end{bmatrix}\begin{bmatrix}I_{n}&O_{n,r}\\\
O_{n,n}&\bar{U}\end{bmatrix}.$

Since $\bar{U}$ is full rank, the following equation holds:

$\displaystyle\mathrm{rank}\begin{bmatrix}\bar{U}^{\top}G&\bar{U}^{\top}A\bar{U}-\lambda_{i}I_{r}\end{bmatrix}=r\quad\forall\lambda_{i},\
i=1,\dots,n.$

This result implies that $(A_{\bar{U}},G_{\bar{U}})$ is controllable.
Therefore, if $(A,G)$ is controllable, then $(A_{\bar{U}},G_{\bar{U}})$ is
also controllable. ∎

######  Remark 2.

If $(A,G,C)$ is stabilizable and detectable, then
$(A_{\bar{U}},G_{\bar{U}},C_{\bar{U}})$ are also stabilizable and detectable.
The proof follows the procedure in Propositions 5 and 6.

If Propositions 5 and 6 hold, we can apply Prop. 1 to Eq. (6). Consequently,
$\tilde{R}(t)$ converges to $\tilde{R}_{\bar{U}}>0$ and
$A_{\bar{U}}-\tilde{R}_{\bar{U}}C_{\bar{U}}^{\top}{(HH^{\top})}^{-1}C_{\bar{U}}$
is stable. However, $\tilde{R}(t)$ can converge to a different solution
depending on the choice of $\bar{U}$. Even in such cases, the following
proposition ensures the uniqueness of
$\tilde{P}_{r}:=\lim_{t\to\infty}\tilde{P}(t)$.

######  Proposition 7.

Suppose $(A,G,C)$ is controllable and observable. Let
$U_{1},U_{2}\in\mathcal{U}$. Subsequently,
$U_{1}\tilde{R}_{U_{1}}U_{1}^{\top}=U_{2}\tilde{R}_{U_{2}}U_{2}^{\top}=\tilde{P}_{r}$
holds.

The proof follows the proof of Prop. 4 in [15].

Next, we present a proposition regarding the stability of the estimation error
system. We then provide the conditions for the boundedness of the estimation
error.

######  Proposition 8.

We assume that $\bar{U}\in\mathcal{U}$. If the eigenvalues of
$A_{\bar{U}}-\tilde{R}_{\bar{U}}C_{\bar{U}}^{\top}{(HH^{\top})}^{-1}C_{\bar{U}}$
are denoted by $\sigma_{i}\ (i=1,\dots,r)$, then the eigenvalues of
$A-\tilde{P}_{r}{C}^{\top}{(HH^{\top})}^{-1}C$ are
$\sigma_{1},\dots,\sigma_{r},\lambda_{r+1}(A),\dots,\lambda_{n}(A)$.

###### Proof.

Let $\bm{u}_{i}$ $(i=1,\dots,r)$ be the eigenvectors of
$A_{\bar{U}}-\tilde{R}_{\bar{U}}C_{\bar{U}}^{\top}{(HH^{\top})}^{-1}C_{\bar{U}}$
and $\bm{\psi}_{i}=\bm{\psi}_{i}(A)$. Note that from the proof of Prop. 5,
$\bm{\psi}_{i}$ $(i=1,\dots,r)$ is described by $\bar{U}\bm{v}_{i}$, where
$\bm{v}_{i}$ is the eigenvector of $\bar{U}^{\top}A\bar{U}$. From Eq. (14),

$\displaystyle\bar{U}(\bar{U}^{\top}A\bar{U}-\tilde{R}_{\bar{U}}C_{\bar{U}}^{\top}{(HH^{\top})}^{-1}C_{\bar{U}})\bm{u}_{k}$
$\displaystyle=$
$\displaystyle(A-\tilde{P}_{r}C^{\top}{(HH^{\top})}^{-1}C)\bar{U}\bm{u}_{k}$
$\displaystyle=$ $\displaystyle\sigma_{k}\bar{U}\bm{u}_{k}.$

By $\bar{U}\bm{u}_{k}\neq 0$, $\sigma_{1},\dots,\sigma_{r}$ are eigenvalues of
$A-\tilde{P}_{r}C^{\top}{(HH^{\top})}^{-1}C$. Moreover, by using
$\bar{U}\bar{U}^{\top}\bm{\psi}_{j}=\bm{\psi}_{j}\ (j=1,\dots,r)$ and $k>r$,
$(I_{n}-\bar{U}\bar{U}^{\top})\bm{\psi}_{k}\neq 0$,
$(I_{n}-\bar{U}\bar{U}^{\top})\tilde{P}_{r}=O_{n,n}$,
$(I_{n}-\bar{U}\bar{U}^{\top})A(I_{n}-\bar{U}\bar{U}^{\top})=A-\bar{U}\bar{U}^{\top}A$,
we can conclude that the following:

$\displaystyle\bm{\psi}_{k}^{\dagger}(I_{n}-\bar{U}\bar{U}^{\top})(A-\tilde{P}_{r}C^{\top}{(HH^{\top})}^{-1}C-\lambda_{k}I_{n})$
$\displaystyle\quad\quad\times(I_{n}-\bar{U}\bar{U}^{\top})\bm{\psi}_{k}$
$\displaystyle=$
$\displaystyle\bm{\psi}_{k}^{\dagger}(I_{n}-\bar{U}\bar{U}^{\top})(A{\bm{\psi}_{k}}-\lambda_{k}{\bm{\psi}_{k}})=0.$

This result implies that $\lambda_{r+1}(A),\dots,\lambda_{n}(A)$ are
eigenvalues of $A-\tilde{P}_{r}C^{\top}{(HH^{\top})}^{-1}C$. Thus, Proposition
8 has been proved. ∎

######  Remark 3.

From Propositions 1, 5, and 6, $\sigma_{1},\dots,\sigma_{r}$ are stable
eigenvalues. Thus, if we choose $r$ such that
$\mathrm{Re}(\lambda_{r+1}(A))<0$, then the LRKB filter is stable.

From Proposition 8, we derive the condition for the necessary rank $r$ for the
bounded estimation error using the LRKB filter.

######  Theorem 3.

Suppose system $(A,G,C)$ is controllable and observable. If $A$ has
$r^{\prime}$ unstable eigenvalues and
$A-\tilde{P}_{r}C^{\top}{(HH^{\top})}^{-1}C$ is stable, then $r\geq
r^{\prime}$.

###### Proof.

From Prop. 8, if $r<r^{\prime}$, then
$\lambda_{r^{\prime}}(A-\tilde{P}_{r}C^{\top}{(HH^{\top})}^{-1}C)=\lambda_{r^{\prime}}(A)$,
which means $A-\tilde{P}_{r}C^{\top}{(HH^{\top})}^{-1}C$ has an unstable
eigenvalue; this contradicts to the assumption. Therefore, $r\geq r^{\prime}$.
∎

From Theorems 1 and 3, we conclude that selecting the rank $r$ is sufficient
when the real parts of the eigenvalues of $A$ differ. This difference must be
greater than the number of unstable eigenvalues of $A$. The number of unstable
eigenvalues can be found, e.g., using the Routh table. However, finding a
suitable rank for linear time-varying or nonlinear systems remains
challenging.

##  V Conclusion

In this study, we analyzed the properties of the Oja flow and a modified LRKB
filter with a general matrix. We removed the restriction on symmetric matrices
and investigated the Oja flow for general matrices. Based on the results, we
revealed the rank $r$ condition for stable estimation.

This study investigated continuous-time LTI systems. In the future, the
applicability of the modified LRKB filter should be analyzed for discrete-time
LTI, linear time-varying, and nonlinear systems.

## References

  * [1] R. E. Kalman, “A New Approach to Linear Filtering and Prediction Problems,” Journal of Basic Engineering, vol. 82, pp. 35–45, 03 1960.
  * [2] R. E. Kalman and R. S. Bucy, “New Results in Linear Filtering and Prediction Theory,” Journal of Basic Engineering, vol. 83, pp. 95–108, 03 1961.
  * [3] B. O. Teixeira, M. A. Santillo, R. S. Erwin, and D. S. Bernstein, “Spacecraft tracking using sampled-data Kalman filters,” IEEE Control Systems Magazine, vol. 28, no. 4, pp. 78–94, 2008.
  * [4] M. E. Pittelkau, “Kalman Filtering for Spacecraft System Alignment Calibration,” Journal of Guidance, Control, and Dynamics, vol. 24, no. 6, pp. 1187–1195, 2001.
  * [5] C. Federico and B. Massimiliano, “Wind speed and wind energy forecast through Kalman filtering of Numerical Weather Prediction model output,” Applied Energy, vol. 99, pp. 154–166, 2012.
  * [6] M. Kunii, T. Miyoshi, and E. Kalnay, “Estimating the Impact of Real Observations in Regional Numerical Weather Prediction Using an Ensemble Kalman Filter,” Monthly Weather Review, vol. 140, no. 6, pp. 1975 – 1987, 2012.
  * [7] Y. Che, X. Peng, L. Delle Monache, T. Kawaguchi, and F. Xiao, “A wind power forecasting system based on the weather research and forecasting model and Kalman filtering over a wind-farm in Japan,” Journal of Renewable and Sustainable Energy, vol. 8, p. 013302, 01 2016.
  * [8] M. Senel, K. Chintalapudi, D. Lal, A. Keshavarzian, and E. J. Coyle, “A Kalman filter based link quality estimation scheme for wireless sensor networks,” in IEEE GLOBECOM 2007-IEEE Global Telecommunications Conference, pp. 875–880, IEEE, 2007.
  * [9] X. Lai, T. Yang, Z. Wang, and P. Chen, “IoT Implementation of Kalman Filter to Improve Accuracy of Air Quality Monitoring and Prediction,” Applied Sciences, vol. 9, no. 9, 2019.
  * [10] S. Gugercin and A. C. Antoulas, “A Survey of Model Reduction by Balanced Truncation and Some New Results,” International Journal of Control, vol. 77, no. 8, pp. 748–766, 2004.
  * [11] H. Shingin and Y. Ohta, “Reduction of state variables based on regulation and filtering performances,” Automatica, vol. 92, pp. 143–154, 2018.
  * [12] G. Kirsten and V. Simoncini, “Order Reduction Methods for Solving Large-Scale Differential Matrix Riccati Equations,” SIAM Journal on Scientific Computing, vol. 42, no. 4, pp. A2182–A2205, 2020.
  * [13] T. Stillfjord, “Adaptive high-order splitting schemes for large-scale differential Riccati equations,” Numerical Algorithms, vol. 78, no. 4, pp. 1129–1151, 2018.
  * [14] S. Bonnabel and R. Sepulchre, The Geometry of Low-Rank Kalman Filters, pp. 53–68.  Springer Berlin Heidelberg, 2013.
  * [15] S. Yamada and K. Ohki, “On a New Low-Rank Kalman-Bucy Filter and its Convergence Property,” Proceedings of the ISCIE International Symposium on Stochastic Systems Theory and its Applications, vol. 2021, pp. 16–20, 2021.
  * [16] E. Oja, “Simplified neuron model as a principal component analyzer,” Journal of Mathematical Biology, vol. 15, pp. 267–273, 1982.
  * [17] E. Oja and J. Karhunen, “On stochastic approximation of the eigenvectors and eigenvalues of the expectation of a random matrix,” Journal of Mathematical Analysis and Applications, vol. 106, no. 1, pp. 69–84, 1985.
  * [18] L. Xu, “Least mean square error reconstruction principle for self-organizing neural-nets,” Neural Networks, vol. 6, no. 5, pp. 627–648, 1993.
  * [19] T. Chen and S.-i. Amari, “Unified stabilization approach to principal and minor components extraction algorithms,” Neural Networks, vol. 14, no. 10, pp. 1377–1387, 2001.
  * [20] S. Yoshizawa, U. Helmke, and K. Starkov, “Convergence analysis for principal component flows,” International Journal of Applied Mathematics and Computer Science, vol. 11, pp. 223–236, 2001.
  * [21] P.-A. Absil, R. Mahony, and R. Sepulchre, Optimization Algorithms on Matrix Manifolds.  Princeton: Princeton University Press, 2008.
  * [22] J.-G. Liu and Z. Liu, “Convergence of Oja’s online principal component flow,” arXiv preprint arXiv:2202.11308v1, 2022.
  * [23] S. Yoshizawa, “Dynamical systems for eigenvalue problems of axisymmetric matrices with positive eigenvalues,” arXiv preprint arXiv:2307.09635v1, 2023.
  * [24] R. S. Bucy and P. D. Joseph, Filtering for Stochastic Processes with Applications to Guidance.  Chelsea Pub Co, 2 ed., 1987.
  * [25] W.-Y. Yan, U. Helmke, and J. B. Moore, “Global analysis of Oja’s flow for neural networks,” IEEE Transactions on Neural Networks, vol. 5, no. 5, pp. 674–683, 1994.
  * [26] S. Yamada and K. Ohki, “Comparison of Estimation Error between Two Different Low-Rank Kalman-Bucy Filters,” in 2021 60th Annual Conference of the Society of Instrument and Control Engineers of Japan (SICE), Sept. 2021.
  * [27] Y. Chikuse, Statistics on Special Manifolds, vol. 174 of Lecture Notes in Statistics.  Springer New York, 2003.


#  A framework for power line inspection tasks with multi-robot systems from
signal temporal logic specifications

Giuseppe Silano1,3, Davide Liuzza2, Luigi Iannelli3, and Martin Saska1
1Giuseppe Silano and Martin Saska are with the Faculty of Electrical
Engineering, Czech Technical University in Prague, Czech Republic, email:
{giuseppe.silano<EMAIL_ADDRESS>Liuzza is with the ENEA
Fusion and Nuclear Safety Department, Italy, email:
<EMAIL_ADDRESS>Silano and Luigi Iannelli are with the
Department of Engineering, University of Sannio, Italy, email:
{giuseppe.silano<EMAIL_ADDRESS>work was partially funded
by the European Union’s Horizon 2020 research and innovation programme AERIAL-
CORE under grant agreement no. 871479, and by ECSEL Joint Undertaking research
and innovation programme COMP4DRONES under grant agreement no. 826610.

###### Abstract

Inspection of power line infrastructures must be periodically conducted by
electric companies in order to ensure reliable electric power distribution.
Research efforts are focused on automating the power line inspection process
by looking for strategies that satisfy different requirements expressed in
terms of potential damage and faults detection. This problem comes up with the
need of safe planning and control techniques for autonomous robots to perform
visual inspection tasks. Such an application becomes even more interesting and
of critical importance when considering a multi-robot extension.

In this paper, we propose to compute feasible and constrained trajectories for
a fleet of quad-rotors leveraging on Signal Temporal Logic (STL)
specifications. The planner allows to formulate rather complex missions
avoiding obstacles and forbidden areas along the path. Simulations results
achieved in MATLAB show the effectiveness of the proposed approach leading the
way to experimental tests on the hardware.

Electric companies invest significantly on the inspection and preemptive
maintenance of the power line infrastructure. The most common strategy is to
perform aerial inspection of the power line corridor, at regular intervals.
The traditional (and the most common) approach to inspection uses a manned
helicopter, equipped with multiple sensors, e.g., light and radar sensors,
visual, infrared and ultra-violet cameras, mounted on gyroscope stabilized
gimbals, and an expert crew, for recording and documenting the relevant data
captured from these sensors. This data is later manually examined to detect
potential faults and damage on different power line components (e.g., cables,
towers, insulators). This process is not only extremely time consuming and
dangerous for human operators, but also very expensive ($1,500 for one hour
flight) and prone to human error [1].

In the last two decades, multiple complementary research directions have been
explored for automating the task of visual inspection. One key direction has
been on developing Unmanned Aerial Vehicles (UAVs) capable of inspecting the
power line corridor and the towers [1]. UAVs can execute two different levels
of inspections, i.e., macro and micro, depending on the wing type. For
example, a fixed-wing UAV cannot hover itself, thus it is used to provide a
macro-level inspection which can be surveying power lines in a large area. On
the other hand, a rotary-wing UAV can perform a micro-level inspection such as
checking for mechanical failures on transmission components. However, using
UAVs to achieve these tasks is particularly challenging due to their limited
battery capacity and to the fact that operating close to lines may be
impractical for the strong electromagnetic field and the several obstacles
(e.g., branches, markers ball) [1].

For all such reasons, having UAVs capable of interpreting high-level, possibly
vague, tasks specifications, and nevertheless plan and execute appropriate
actions for a particular context in which the system is operating, is more and
more important. Symbolic control proposes to fulfill this need by
automatically designing feedback controllers that lead to the satisfaction of
formal logic specifications. Temporal-logic (TL) based motion planning
techniques can be used for this purpose. In particular, Signal Temporal Logic
(STL) provides formal high-level languages that can describe planning
objectives more complex than well-suited point-to-point navigation algorithms
[2]. The task specification is given as a temporal logic formula w.r.t. the
discretized abstraction of the robot motion modelled as a finite transition
system [3]. Then, a high-level discrete plan is found by off-the-shelf model-
checking algorithms [4] given the finite transition system and the task
specification. This discrete plan is then implemented through a corresponding
low-level hybrid controller.

In this paper, we propose a framework to encode missions for a fleet of quad-
rotors as STL specifications. Then, using motion primitives for quad-rotor
trajectory generation [5], we construct an optimization problem to generate
optimal strategies that satisfy such specifications. The proposed approach
generates feasible dynamic trajectories accounting for velocity and
accelerations constraints of the vehicles that satisfy missions limitations,
too. Simulations results in MATLAB show the effectiveness of the proposed
approach.

Let us consider a continuous-time dynamical system $\pazocal{H}$ and its
discretized time version $x^{+}=f(x,u)$, where $x,x^{+}\in
X\subset\mathbb{R}^{n}$ are the current and next state of the system,
respectively, $u\in U\subset\mathbb{R}^{m}$ is the control input and $f\colon
X\times U\rightarrow X$ is differentiable in both arguments. The system’s
initial state is denoted by $x_{0}$ and takes values from some initial set
$X_{0}\subset\mathbb{R}^{n}$. Let $T_{s}=t_{k+1}-t_{k}$ for all
$k\in\mathbb{N}_{\geq 0}$ be the sampling period and $T\in\mathbb{R}_{\geq 0}$
be a trajectory duration, we can write the time interval as
$t=(0,T_{s},2T_{s},\dots,NT_{s})$ with $NT_{s}=T$. Therefore, given an initial
state $x_{0}$ and a finite control input sequence
$\mathbf{u}=(u_{0},u_{1},\dots,u_{N-1})$, a trajectory of the system is the
unique sequence of states $\mathbf{x}=(x_{0},x_{1},\dots,x_{N})$ with
$x_{k+1}=f(x_{k},u_{k})$.

The trajectory generator is designed to make the closed loop system satisfying
a specification expressed in STL. STL is a logic that allows the succinct and
unambiguous specification of a wide variety of desired system behaviors over
time, such as “The quad-rotor reaches the goal within 10 time units while
always avoiding obstacles”. Formally, let
$M=\\{\mu_{1},\mu_{2},\dots,\mu_{L}\\}$ be a set of real-valued functions of
the state $\mu_{k}\colon X\rightarrow\mathbb{R}$. For each $\mu_{k}$ define
the predicate $p_{k}\coloneqq\mu_{k}(x)\geq 0$. Set
$AP\coloneqq\\{p_{1},p_{2},\dots,p_{L}\\}$. Thus each predicate defines a set
over the system state space, namely $p_{k}$ defines $\\{x\in X|\mu_{k}(x)\geq
0\\}$. Let $I\subset\mathbb{R}$ denote a nonsingleton interval, $\top$ the
Boolean True, $p$ a predicate, $\neg$ and $\wedge$ the Boolean negation and
AND operators, respectively, and $\pazocal{U}$ the Until temporal operator. An
STL formula $\varphi$ is built recursively from the predicates using the
grammar
$\varphi\coloneqq\top|p|\neg\varphi|\varphi_{1}\wedge\varphi_{2}|\varphi_{1}\pazocal{U}_{I}\varphi_{2}$,
where $p$ is a predicate and $\varphi_{1}$ and $\varphi_{2}$ are STL formulas.
Informally, $\varphi_{1}\pazocal{U}_{I}\varphi_{2}$ means that $\varphi_{2}$
must hold at some point in $I$, and until then, $\varphi_{1}$ must hold
without interruption. The disjunction ($\vee$), implication
($\Longrightarrow$), Always ($\square$) and Eventually ($\lozenge$) operators
can be also defined. Formally, the pointwise semantics of an STL formula
$\varphi$ defines what it means for a system trajectory $\mathbf{x}$ to
satisfy $\varphi$ [6].

Designing a controller that satisfies a STL formula $\varphi$ is not always
enough. In a dynamic environment, where the system must react to new
unforeseen events, it is useful to have a margin of maneuverability. That is,
it is useful to control the system such that to maximize a degree of
satisfaction of the formula. The latter can be defined and computed using the
robust semantics of TL [7, 8].

Starting from the definition of a robust STL formula, indicated as
$\rho_{\varphi}(\mathbf{x},t)$, we can compute control inputs $\mathbf{u}$ by
maximizing the robustness over the set of finite state sequences $\mathbf{x}$.
The obtained sequence $\mathbf{u}^{\star}$ is valid if
$\rho_{\varphi}(\mathbf{x}^{\star},t)$ is positive, where $\mathbf{x}^{\star}$
and $\mathbf{u}^{\star}$ obey to the dynamical system $\pazocal{H}$. Thus, the
goal is to find a provably correct control scheme for a robot system which
makes it meeting a control objective $\varphi$ expressed in temporal logic.
So, let $\epsilon>0$ be a desired minimum robustness, we solve the following
problem

$\begin{split}&\operatorname*{maximize}_{\mathbf{u}\in
U^{N-1}}\;\;{\rho_{\varphi}(\mathbf{x})}\\\
&\quad\text{s.t.}\leavevmode\nobreak\ \qquad x_{k+1}=f(x_{k},u_{k}),\forall
k=\\{0,1,\dots,N-1\\},\\\ &\qquad\;\;\qquad x_{k}\in X,u_{k}\in U,\forall
k=\\{0,1,\dots,N\\},\\\
&\qquad\;\;\qquad\rho_{\varphi}(\mathbf{x})\geq\epsilon\end{split}.$ (1)

The same formalization applies if $f$, $x$ and $u$ represent the overall
dynamics, state and input of a multi-robot system. Because $\rho_{\varphi}$
uses the non-differentiable functions $\max$ and $\min$ [7, 8], solving (1)
can be achieved by using Mixed-Integer Programming solvers, non-smooth
optimizers or stochastic heuristics [6]. However, as recently shown in [9], it
is more efficient and reliable to approximate the non-differentiable objective
$\rho_{\varphi}$ by a smooth (infinitely differentiable) function
$\tilde{\rho}_{\varphi}$ and to solve the resulting optimization problem using
Sequence Quadratic Programming.

To come up with a trajectory that satisfies the vehicle constraints, the
motion primitives defined in [5] were considered. Such a method allows to
obtain a rapid generation and feasibility verification of motion primitives
for quad-rotors. The results of the modified optimization problem is a set of
splines [5, eq. (22)] whose parameters can be tuned to achieve a desired
motion fixing a combination of position, velocity, and acceleration at the
start and end points. Therefore, the objective function can be reformulated
replacing $\rho_{\varphi}(\mathbf{x})$ with its smooth version accounting for
the mathematical formulation of the trajectory generator $L$, i.e.,
$\tilde{\rho}_{\varphi}\left(L(\mathbf{x})\right)$. To track the obtained
references, any of the quad-rotor control algorithms within the literature
[10] can be used.

We show the effectiveness of the proposed approach through two numerical
examples: (i) a multi-drone reach and avoid problem in a constrained
environment (eq. (2a), https://youtu.be/SHThWKYw1v8.), and (ii) a multi-
mission example where several drones have to fly one of two missions in the
same environment (see eq. (2b), https://youtu.be/uqYq5xmCckY.).

$\begin{split}&\varphi_{\mathrm{safe}}^{i,j}=\square_{[0,T]}\varphi_{\mathrm{dist}}^{i,j}=\square_{[0,T]}\left(\lVert{{}^{W}\mathbf{r}}^{i}-{{}^{W}\mathbf{r}}^{j}\rVert\geq\delta_{\mathrm{min}}\right),\\\
&\varphi_{\mathrm{mra}}=\wedge_{k=1}^{N}\varphi_{\mathrm{ra}}^{k}\bigwedge\wedge_{k=1}^{N}\left(\wedge_{i\neq
j}\varphi_{\mathrm{safe}}^{i,j}\right),\end{split}$

(2a)

$\begin{split}\varphi_{\mathrm{pli}}^{N}&=\wedge_{k=1}^{N}\Bigl{(}\square_{[0,T]}\left(\varphi_{\mathrm{dist}}^{k}\wedge\varphi_{\mathrm{ws}}^{k}\right)\Bigr{)}\bigwedge\lozenge_{[0,T]}\biggl{(}\left(\wedge_{k=1}^{\nicefrac{{N}}{{2}}}\varphi_{\mathrm{pole1}}^{k}\wedge\varphi_{\mathrm{pole4}}^{k}\right)\\\
&\bigwedge\left(\wedge_{k=\nicefrac{{N}}{{2}}+1}^{N}\varphi_{\mathrm{pole2}}^{k}\wedge\varphi_{\mathrm{pole3}}^{k}\right)\pazocal{U}_{[0,\nicefrac{{T}}{{2}}]}\left(\wedge_{k=\nicefrac{{N}}{{2}}+1}^{N}\varphi_{\mathrm{pole3}}^{k}\wedge\varphi_{\mathrm{pole2}}^{k}\right)\biggr{)}\\\
\end{split}.$

(2b)

## References

  * [1] H. Baik and J. Valenzuela, “Unmanned Aircraft System Path Planning for Visually Inspecting Electric Transmission Towers,” _Journal of Intelligent & Robotic Systems_, vol. 95, pp. 1097–1111, 2018.
  * [2] S. Karaman and E. Frazzoli, “Smapling-based algorithms for optimal motion planning,” _The International Journal of Robotics Research_ , vol. 7, no. 30, pp. 846–894, 2011.
  * [3] M. Kloetzer and C. Belta, “Automatic Deployment of Distributed Teams of Robots From Temporal Logic Motion Specifications,” _IEEE Transactions on Robotics_ , vol. 26, no. 1, pp. 48–61, 2010.
  * [4] A. Bathia, M. R. Maly, L. E. Kabraki, and M. Y. Vardi, “Motion planning with complex goals,” _IEEE Robotics and Automation Magazine_ , vol. 3, no. 18, pp. 55–64, 2011.
  * [5] M. W. Mueller, M. Hehn, and R. D’Andrea, “A Computationally Efficient Motion Primitive for Quadrocopter Trajectory Generation,” _IEEE Transactions on Robotics_ , vol. 31, no. 6, pp. 1294–1310, 2015.
  * [6] V. Raman, A. Donzé, M. Maasoumy, R. M. Murray, A. Sangiovanni-Vincentelli, and S. A. Seshia, “Model predictive control with signal temporal logic specifications,” in _53rd IEEE Conference on Decision and Control_ , 2014, pp. 81–87.
  * [7] G. Fainekos and G. Pappas, “Robustness of temporal logic specifications for continuous-time signals,” _Theoretical Computer Science_ , vol. 410, pp. 4262–4291, 2009.
  * [8] A. Donzé and O. Maler, “Robust satisfaction of temporal logic over real-valued signals,” in _Formal Modeling and Analysis of Timed Systems_ , K. Chatterjee and T. A. Henzinger, Eds., 2010, pp. 92–106.
  * [9] Y. V. Pant, H. Abbas, and R. Mangharam, “Smooth operator: Control using the smooth robustness of temporal logic,” in _2017 IEEE Conference on Control Technology and Applications_ , 2017, pp. 1235–1240.
  * [10] T. P. Nascimento and M. Saska, “Position and attitude control of multi-rotor aerial vehicles: A survey,” _Annual Reviews in Control_ , vol. 48, pp. 129–146, 2019.


11institutetext: College of Computer Science and Technology, Zhejiang
University, Hangzhou City, Zhejiang Province, China 22institutetext: ZJU-
Hangzhou Global Scientific and Technological Innovation Center, Hangzhou City,
Zhejiang Province, China 22email:
gaosde@zju.edu.cn,fuyichao@zju.edu.cn,lk2017<EMAIL_ADDRESS>
# Contrastive Knowledge Amalgamation for Unsupervised Image Classification

Shangde Gao Contributed equally.1122 Yichao Fu ⋆ 11 Ke Liu⋆ 11 Yuqiang Han
1122

###### Abstract

Knowledge amalgamation (KA) aims to learn a compact student model to handle
the joint objective from multiple teacher models that are are specialized for
their own tasks respectively. Current methods focus on coarsely aligning
teachers and students in the common representation space, making it difficult
for the student to learn the proper decision boundaries from a set of
heterogeneous teachers. Besides, the KL divergence in previous works only
minimizes the probability distribution difference between teachers and the
student, ignoring the intrinsic characteristics of teachers. Therefore, we
propose a novel Contrastive Knowledge Amalgamation (CKA) framework, which
introduces contrastive losses and an alignment loss to achieve intra-class
cohesion and inter-class separation. Contrastive losses intra- and inter-
models are designed to widen the distance between representations of different
classes. The alignment loss is introduced to minimize the sample-level
distribution differences of teacher-student models in the common
representation space. Furthermore, the student learns heterogeneous
unsupervised classification tasks through soft targets efficiently and
flexibly in the task-level amalgamation. Extensive experiments on benchmarks
demonstrate the generalization capability of CKA in the amalgamation of
specific task as well as multiple tasks. Comprehensive ablation studies
provide a further insight into our CKA.

###### Keywords:

Knowledge amalgamation Contrastive learning.

##  1 Introduction

Reusing pre-trained models to get lite ones for reducing computation costs of
training a new one from scratch has been a trending research topic in recent
years [8, 7]. Knowledge Distillation (KD) methods [6] train a light-weight
target model (the _“student"_ model) by learning from a well-trained
cumbersome model (the _“teacher"_ model), which improves the performance of
students with any architectures compared to the models trained from scratch.
Knowledge Amalgamation (KA) [13, 11] aims to train a versatile student model
by transferring knowledge from multiple pre-trained teachers. The above method
requires mapping teachers and student to a common representation space. The
student learn similar intermediate features through the aggregated cues from
the pre-trained teachers. Further, by integrating probability knowledge from
pre-trained teachers using KL divergence, student can predict the joint of
teachers’ label sets.

However, complex optimization designs are required for heterogeneous teachers
in previous works. Besides, direct application of previous KA methods to
downstream tasks causes severe performance degradation because of domain
shifts, additional noise, as well as information loss in feature projections.
Moreover, due to the imperfection of pre-trained teachers and absence of human
annotation, the supervision signals for students are confused.

In this work, we endeavor to explore an efficient and effective KA scheme for
unsupervised image classification. We aim to transfer knowledge as much as
possible from pre-trained teachers who specialize in heterogeneous
unsupervised image classification tasks to a compact and versatile student.
For example, if one teacher classifies cars and the other classifies
airplanes, the student should be able to classify both cars and airplanes. To
achieve this, we first extend the contrastive learning paradigm to the
knowledge fusion environment, for two reasons. Firstly, CL can effectively
push positive sample pairs together and pull negative sample pairs apart
without the need for manual annotations. Additionally, different teacher and
student models are natural augmentation schemes, and their combination
significantly increases the number of positive and negative samples for
training the student. Secondly, supervised contrastive loss models have been
shown to outperform traditional cross-entropy losses[9]. Thus, they can be
effectively used in teacher pre-training to alleviate the incompleteness and
unreliability of supervising teacher models.

We propose a novel Contrastive Knowledge Amalgamation, refered to as CKA, by
implementing the CKA framework via DNNs for unsupervised classification.
Concretely, we first construct a common representation space based on the
shared multilayer perceptron (MLP), and design contrastive and alignment
losses to achieve intra-class cohesion and inter-class separation of samples.
As a way of unsupervised learning, the contrastive loss intra- and inter-
models aims to enlarge the distance between feature representations of
different sample categories and reduce the distance between feature
representations of the same sample category. Besides, alignment losses are
proposed to minimize the sample-level distribution difference between
different models. Apart from learning the teachers’ features, a soft target
distillation loss finally is designed to effectively and flexibly transfer
probability knowledge from pre-trained teachers to a student, enabling the
student to make inferences similar to or the same as the teachers’ during
task-level amalgamation.

The contributions of this work are summarized as follows:

  * •

We propose a novel model reuse paradigm to supervise the student model without
annotations, named CKA, which introduces contrastive losses and an alignment
loss to achieve intra-class cohesion and inter-class separation.

  * •

We design a soft target distillation loss to effectively transfer knowledge
from the pre-trained teachers to a student in the output probability space.

  * •

Extensive experiments on standard benchmarks demonstrate that CKA provides
more accurate supervision, and is generalizable for amalgamating heterogeneous
teachers.

##  2 Related Works

###  2.1 Knowledge Distillation & Knowledge Amalgamation

Knowledge distillation (KD) [6, 16] is a method of transferring knowledge from
one model to another. However, existing approaches are still performed under a
single teacher-student relationship with a sharing task, are not applicable to
multiple and heterogeneous teachers. Knowledge amalgamation (KA) aims to
acquire a compact student model capable of handling the comprehensive joint
objective of multiple teacher models, each specialized in its own task. There
are two kinds of approaches: (1) _Homogeneous_ KA, where all teachers and
students have identical network architectures [13]. (2) _Heterogeneous_ KA,
where each teacher has different architecture and specializes in its own class
set [14, 11]. Among these, [14] matches the outputs of students to the
corresponding teachers, while [11] aligns the features of students and
teachers in a shared latent space by minimizing the maximum mean discrepancy.
However, when facing with the imperfect teachers with unreliable supervisions,
previous studies suffer from conflicting supervisions in the student training
process, which significantly harms the performance of the student model. To
the best of our knowledge, it is the first time to explore the CKA paradigm
for unsupervised classification tasks.

###  2.2 Contrastive Learning

Contrastive Learning is an unsupervised learning method where supervision is
automatically generated from the data. Currently, contrastive learning (CL)
has achieved state-of-the-art performance in representation learning [4, 1,
2]. SimCLR [1] proposes the proposal by performing data augmentation on the
raw input data and mapping it to a feature space, constructing a contrastive
loss (i.e., InfoNCE loss) to maximize the similarity between positive pairs
and minimize the similarity between negative pairs. BYOL [4] and SimSiam [2]
extend the work by designing their losses to measure the similarity between
positive samples, effectively eliminating the need for negative samples.
However, all of these approaches are tailored for single-model and single-
task. Our method extends the concept of CL to a knowledge amalgamation
environment, designing intra-and inter- model contrastive losses to explore
the model-agnostic semantic similarity and further apply them to downstream
unsupervised multi-classification tasks.

##  3 Problem Formulation

We define the problem of knowledge amalgamation as follows. Assume that we are
given $\mathcal{T}=\left\\{\mathcal{T}_{t}\right\\}_{t=1}^{N}$ well pre-
trained teachers, where each teacher $\mathcal{T}_{t}$ specializes a distinct
classification task, i.e., a set of full labeled classes
$\mathcal{T}_{t}=(\mathcal{D}_{t};\mathcal{Y}^{t})$. Our proposal is to learn
a versatile student with an unlabeled dataset
$\mathcal{D}=\bigcup_{t=1}^{N}{\mathcal{D}_{t}}$, which is able to perform
predictions over the comprehensive class set of distinct-task teachers,
$\mathcal{Y}=\bigcup_{t=1}^{N}{\mathcal{Y}^{t}}$. In our KA setting, $N$ tasks
$\mathcal{T}=\left\\{\mathcal{T}_{t}\right\\}_{t=1}^{N}$ can be built for
either the same or cross dataset. Without loss of generality, we assume that
for any two tasks $\mathcal{T}_{i},\mathcal{T}_{j}\in\mathcal{T}$, their
specialties are totally disjoint, i.e.,
$\mathcal{Y}^{i}\cap\mathcal{Y}^{j}=\oslash$.

##  4 Approach

This work is aimed to build a contrastive knowledge amalgamation framework,
and implement it by DNNs for unsupervised image classification. Knowledge
amalgamation is particularly challenging when teacher-student structures are
heterogeneous and data annotation is not available.

Figure 1: The overflow of contrastive knowledge amalgamation.

To tackle the difficulty, we first leverage the distance between feature
representations of the samples, and introduce contrastive and alignment losses
to achieve intra-class coherence and inter-class separation of the feature
representations. Additionally, we design a soft-target distillation loss to
effectively transfer the soft-target probability knowledge from pre-trained
teachers to the student. The overview of the proposed CKA is shown in Figure.
1, in which the knowledge of pre-trained teachers is fixed. By training the
student model in downstream tasks, the student is capable of making inferences
that are similar or identical to those of their teachers.

###  4.1 Margin-based Intra- and Inter-model Contrast

As there are no annotated data available, we novely use contrastive learning
(CL) to construct supervision for guiding the student. CL aims to maximize the
similarities of positive pairs while minimizing those of negative ones [1].
The characteristics of pairs can be defined by different criteria. Motivated
by this, we develop two types of contrastive losses, including edge-based
student-internal contrast (intra-model contrast) and distance-based teacher-
student model contrast (inter-models contrast), to increase the distance
between different sample class feature representations and decrease the
distance between the same sample class feature representations. The overall
schematic is shown in Figure. 2.

Figure 2: Illustrations of intra- and inter-model contrast loss via the
teacher-student pair.

####  4.1.1 Margin-based Intra-model Contrast

To begin with, we describe the standard contrastive loss term, following the
most popular setups of SimCLR [1], which is defined as:

$\mathcal{L}(\tilde{x},\hat{x})=-\log\frac{e^{{s}(\tilde{z},\hat{z})/\tau}}{e^{{s}(\tilde{z},\hat{z})/\tau+\sum_{\bar{z}\in\Lambda^{-}}e^{{s}\left(\tilde{z},\bar{z})\right)/\tau}}}$
(1)

Here, by way of randomized data augmentation Aug($\cdot$), two different views
$\tilde{x}$ and $\hat{x}$ for the input sample $x$ are generated. The two
images are then fed into an encoder network $\mathcal{E}(x)$, followed by a
two-layer nonlinear projection head MLP $h(\cdot)$, yielding a pair of
$L_{2}$-normalized positive embeddings $\hat{z}=h(\mathcal{E}(\hat{x}))$ and
$\tilde{z}=h(\mathcal{E}(\tilde{x}))$. $\bar{z}\in\Lambda^{-}$ represents the
negative sample in a mini-batch. $s(\cdot,\cdot)$ declares the _cosine
similarity_ for measuring the relationship between embedding pair $\tilde{z}$
and $\hat{z}$ (resp. $\tilde{x}$ and $\hat{x}$), formulated as:

$s\left(\tilde{z},\hat{z}\right)=\frac{\left(\tilde{z}\right)\left(\hat{z}\right)^{\top}}{\left\|\tilde{z}\right\|\cdot\left\|\hat{z}\right\|}$
(2)

To prevent the loss from being dominated by easy negatives (different class
samples with little similarity), a constant margin $\alpha$ is introduced that
only negative pairs with similarity larger than $\alpha$ contribute to the
contrastive loss in Eqn. 1. Formally, the margin-based intra-model contrastive
loss for training the student model is denoted as:

$\mathcal{L}_{intra}=\left(1-s\left(\tilde{z},\hat{z}\right)\right)+\sum_{\bar{z}\in\Lambda^{-}}\left(s\left(\tilde{z},\bar{z}\right)-\alpha\right)$
(3)

####  4.1.2 Distance-based inter-model contrast

For inter-model contrast, data across models are embedded as point
distributions in high-dimensional vector spaces. To measure the inter-model
distance between those two point distributions, we model two _metric measure
spaces_ (mm-spaces) $\mathcal{X}=\left(X,d_{X},\mu\right)\text{ and
}\mathcal{Y}=\left(Y,d_{Y},\nu\right)$, where data $X$ (resp. $Y$) is a
complete separable set endowed with a distance $d_{X}$ and a positive Borel
measure $\mu\in\mathcal{M}+(X)$. Those two mm-spaces are considered up to
isometry (denoted $\mathcal{X}\sim\mathcal{Y}$), meaning that there is a
bijection $\psi:\operatorname{spt}(\mu)\rightarrow\operatorname{spt}(\nu)$
(where $\operatorname{spt}(\mu)$ is the support of $\mu$) such that
$d_{X}(x,y)=d_{Y}(\psi(x),\psi(y))$ and $\psi_{\sharp}\mu=\nu$. Here
$\psi_{\sharp}$ is the push-forward operator.

Specifically, let $\mu\in\mathcal{P}\left(\mathbb{R}^{p}\right)$ and
$\nu\in\mathcal{P}\left(\mathbb{R}^{q}\right)$ with $p\neq q$ to be discrete
measures on mm-spaces with $\mu=\sum_{i=1}^{n}a_{i}\delta_{x_{i}}$ (here
$\delta_{x_{i}}$ is the mass at $x_{i}$) and
$\nu=\sum_{i=1}^{n}b_{j}\delta_{y_{j}}$ of supports $X$ and $Y$, where
$a\in\Sigma_{n}$ and $b\in\Sigma_{m}$ are simplex histograms. The distance
$\mathcal{D}$ between those points is defined as:

$\footnotesize\small\mathcal{D}\left(\mathcal{X},\mathcal{Y}\right)^{q}=\sum_{i,j,k,l}\left|d_{X}\left(x_{i},x_{k}\right)-d_{Y}\left(y_{j},y_{l}\right)\right|^{q}\pi_{i,j}\pi_{k,l}$
(4)

Here
$d_{X}(x_{i},x_{k}):\mathbb{R}^{p}\times\mathbb{R}^{p}\rightarrow\mathbb{R}_{+}$,
measures the euclidean distance between sample points $x_{i}$ and $x_{k}$ in
$\mu$. The intuition underpinning the definition of this distance is that
there exists a fuzzy correspondence map $\pi\in\mathcal{P}(X\times Y)$ between
the points of the distributions, which tends to associate pairs of points with
similar distances within each pair: the more similar $d_{X}(x_{i},x_{k})$ is
to $d_{Y}(y_{j},y_{l})$, the stronger the transport coefficients $\pi_{i,j}$
and $\pi_{k,l}$ are. From a semantic perspective, by simultaneously learning
the model structures of both the teacher and student, this distance can
measure the similarity between samples, reducing the distance between the
feature representations of similar sample classes and increasing the distance
between feature representations of dissimilar sample classes.

Given a mini-batch size $B$ of feature maps $P\in\mathbb{R}^{B\times c\times
h\times w}$ and $Q_{t}\in\mathbb{R}^{B\times c\times h\times w}$ extracted
from the student encoder and $t$-th teacher encoder, where $c$, $h$, and $w$
denote the number of channel, height and width of the feature maps
respectively. For simplicity, we omit the superscripts and subscripts and
denote the feature maps of two different models as $P$ and $Q$. The distance
metric on $P$ and $Q$ is designed firstly to guide the contrast across
different models, i.e., inter-model contrast. To this end, we first reshape
$P$ and $Q$ to $\mathbb{R}^{B\times m}$, i.e.,
$P=\left[p_{1},p_{2},\ldots,p_{B}\right]$ and
$Q=\left[q_{1},q_{2},\ldots,q_{B}\right]$, where $m=c\times h\times w$ is the
feature vectors. The transport map $\pi^{p},\pi^{q}\in\mathbb{R}^{B\times B}$
for $P$ and $Q$ can be derived by:

$\pi_{i,j}=\frac{e^{-d(p_{i},p_{j})}}{\sum_{j=1}^{N}e^{-d(p_{i},p_{j})}}$ (5)

where $d(\cdot,\cdot)$ is the mm-space distance between two instances $p_{i}$
and $p_{j}$. Unless stated otherwise, euclidean distance is used in our
experiments.

As for any $k$-th row vector in $\pi^{p}$ and $\pi^{q}$, $\pi^{p}_{k}$ and
$\pi^{q}_{k}$ can be termed as positive pairs because they both semantically
illustrate the distance of $k$-th sample and others in the mini-batch $N$,
regardless of the model representation. Our distance-based inter-model
contrastive loss, discovering fine-gained sample similarity matching between
the student and each teacher, can be defined as:

$\mathcal{L}_{inter}=\sum_{t=1}^{N}\left(\left(1-s\left(\pi^{p},\pi^{q}_{t\text{+}}\right)\right)+\sum_{\pi^{q}_{t\text{-}}\in\Lambda^{-}}s\left(\pi^{p},\pi^{q}_{t\text{-}}\right)\right)$
(6)

where $\pi^{q}_{t\text{-}}$ and $\pi^{q}_{t\text{+}}$ denote the distance-
based negative and positive pairs.

###  4.2 Common Feature Alignment

To enable a student to mimic the aggregated hints from heterogeneous teachers,
a shared multilayer perceptron (MLP) is designed for mapping all features to a
common latent space. Specifically, a $1\times 1$ kernel convolution is added
after the backbone network of each model separately, thereby unifying the
outputs of different models into the same channel, which is taken to be the
input of MLP and set to 256 in our implementation.

As represented in CFL [11], we adopt the Maximum Mean Discrepancy (MMD) to
measure the discrepancy between the output features of the student and that of
teachers in the unit ball of a reproducing kernel Hilbert space [3]. Take a
teacher-student pair as an example, we extract the mini-batch common space
features with the designed shared MLP and represent them as
$f_{\mathcal{S}},f_{\mathcal{T}}\in\mathbb{R}^{B\times d}$, of which $d$
denotes the output dimension of the MLP and is set to 128 in our
implementation. An empirical $l_{2}$ _norm_ approximation to the MMD distance
of $f_{\mathcal{S}}$ and $f_{\mathcal{T}}$ is computed as follow:

$\mathrm{MMD}=\frac{1}{B}\left\|\sum_{i=1}^{B}\phi\left(f_{\mathcal{T}}^{i}\right)-\sum_{j=1}^{B}\phi\left(f_{\mathcal{S}}^{j}\right)\right\|_{2}^{2}$
(7)

where $\phi$ is an explicit mapping function. The extension of multi-kernel
formulation of MMD can then be defined as:

$\displaystyle\mathrm{MMD}^{2}[K,f_{\mathcal{S}},f_{\mathcal{T}}]=$
$\displaystyle
K\left(f_{\mathcal{S}},f_{\mathcal{S}}\right)-2K\left(f_{\mathcal{T}}^{i},f_{\mathcal{S}}^{j}\right)+$
(8) $\displaystyle K\left(f_{\mathcal{T}},f_{\mathcal{T}}\right)$

$K$ is defined as the convex combination of $m$ PSD kernel:

$\mathcal{K}=\left\\{K=\sum_{u=1}^{m}\sigma_{u}K_{u}:\sum_{u=1}^{m}\sigma_{u}=1,\sigma_{u}\geq
0,\forall u\right\\}$ (9)

here $\mathcal{K}$ denotes the multi-prototypical kernel set. The constraints
on coefficients $\left\\{\sigma_{u}\right\\}$ are imposed to guarantee that
the derived multi-kernel $K$ is characteristic.

The process of aligning each teacher and student is equivalent to minimizing
the MMD distance between them. This can achieve intra-class cohesion of
similar samples. We aggregate all such MMDs between $N$ pairs of teachers and
students, and the overall alignment loss $\mathcal{L}_{align}$ in the shared
MLP can be written as:

$\mathcal{L}_{align}=\sum_{i=1}^{N}\mathrm{MMD}\left(f_{\mathcal{S}},f_{\mathcal{T}_{t}}\right)$
(10)

###  4.3 Soft-target Distillation

Apart from learning the teacher’s features, the student is also expected to
produce identical or similar inferences as the teachers do. We thus also take
the teachers’ predictions by feeding unlabelled input samples to them and then
supervise the student’s training. As there is no annotation available for each
instance $x$ in the target dataset $\mathcal{D}_{s}$, the predictions of pre-
trained teachers can be constructed as supervision for guiding the student,
named as _soft-target distillation_.

Table 1: Statistics of datasets used in this paper. Dataset | Images | Categories | Train/Test
---|---|---|---
CUB-200-2011 | 11,788 | 200 | 5,994/5,794
Stanford Dogs | 20,580 | 120 | 12,000/8,580
Stanford Cars | 16,185 | 196 | 8,144/8,041
FGVC-Aircraft | 102,000 | 102 | 6,667/3,333

Specifically, we first feed $x$ into each $T_{i}$ to obtain the golden label
probability distribution $\Phi(x;T_{i})$ in the _softmax layer_ , and then
concatenate them together for training the student by minimizing the KL-
divergence between their probability distribution:

$\mathcal{L}_{std}=\sum_{x\in\mathcal{D}_{s}}\operatorname{KL}(\Phi(x,S)\|\Phi(x,T))$
(11)

where $\Phi(x,S)$ and $\Phi(x,T)$ denote the _softmax_ probability
distribution of the student and that of the concatenated teachers for input
$x$, respectively.

Considering the weighted sum of contrastive losses (including
$\mathcal{L}_{intra}$, and $\mathcal{L}_{inter}$), alignment loss
$\mathcal{L}_{align}$ and soft-target distillation loss $\mathcal{L}_{std}$
together, the total training objective of our CKA can be described as:

$\small\mathcal{L}=\lambda_{intra}\mathcal{L}_{intra}+\lambda_{inter}\mathcal{L}_{inter}+\lambda_{a}\mathcal{L}_{align}+\lambda_{d}\mathcal{L}_{std}$
(12)

##  5 Experiments

In this section, we evaluate the proposed method on standard benchmarks and
compare the results with the recent state of the arts. We also conduct
ablation studies to validate the effect of the major components.

###  5.1 Experiments Setup

##### Datasets

We evaluate our proposed CKA on four widely used benchmarks, i.e.,
CUB-200-2011 [15], Stanford Cars [10], Stanford Dogs [9], and FGVC-Aircraft
[12]. The detailed statistics are summarized in Table 1.

##### Implementation Details

We adopt the _resnet_ family [5] including _resnet_ -18, _resnet_ -34, and
_resnet_ -50, as our model samples. Besides, all the teachers are first pre-
trained as [10] and fine-tuned to heterogeneous tasks. To construct
heterogeneous tasks on the given datasets, we split all the categories into
non-overlapping parts of equal size to train the teachers. The trained teacher
model weights are frozen during the student training process. In student
training phrase, data augmentation is performed via Random ResizedCrop, Random
ColorJitter, Random HorizontalFlip, and Random GaussianBlur while in testing,
Center Crop is used. During training, the learning rate is set to 0.0005, and
the cosine decay is used; the weight decay is set to 0.0005; Adam is used as
the optimizer, and the batch size is set to 64; a total of 100 epochs are
trained. All experiments are completed with GPUs of RTX 2080 Ti 11GB and CPUs
of Intel. There are several hyper-parameters involved in our method, including
$\alpha$ in Eqn. 3, set to 0.4, for alleviating the dominance of negative
sample pairs; $\lambda_{intra}$, $\lambda_{inter}$, $\lambda_{a}$ and
$\lambda_{d}$ for the final CKA loss in Eqn. 12, are set to $\lambda_{intra}$
= $\lambda_{inter}$ = $\lambda_{d}$ = 1 and $\lambda_{a}$ = 10.

Table 2: Comparison of different methods on comprehensive classification tasks. Best results are shown in bold. Method | Size | Dogs | Cars | CUB | Aircraft | Average
---|---|---|---|---|---|---
Supervised | 163M |  83.62 $\pm$ 0.00 |  89.64 $\pm$ 0.00 |  72.68 $\pm$ 0.00 |  82.78 $\pm$ 0.00 | 82.14
Teacher1 | 130M |  66.64 $\pm$ 0.00 |  70.33 $\pm$ 0.00 |  65.37 $\pm$ 0.00 |  63.01 $\pm$ 0.00 | 66.80
Teacher2 | 240M |  72.03 $\pm$ 0.00 |  87.85 $\pm$ 0.00 |  66.12 $\pm$ 0.00 |  81.12 $\pm$ 0.00 | 76.60
Ensemble | 370M |  73.90 $\pm$ 0.22 |  77.08 $\pm$ 0.64 |  68.25 $\pm$ 0.00 |  75.76 $\pm$ 0.00 | 73.38
Vanilla KD | 240M |  76.16 $\pm$ 0.60 |  80.39 $\pm$ 0.31 |  69.94 $\pm$ 0.79 |  78.00 $\pm$ 0.01 | 76.06
CFL | 240M |  76.23 $\pm$ 0.26 |  81.12 $\pm$ 0.21 |  70.67 $\pm$ 0.97 |  79.98 $\pm$ 0.22 | 76.86
CKA-Intra | 240M |  78.89 $\pm$ 0.59 |  82.33 $\pm$ 0.31 |  71.07 $\pm$ 0.04 |  79.02 $\pm$ 0.21 | 77.71
CKA-Inter | 240M |  79.72 $\pm$ 0.60 | 82.95 $\pm$ 1.20 | 71.49 $\pm$ 0.25 |  80.45 $\pm$ 0.51 | 78.46
CKA | 240M | 79.76 $\pm$ 0.09 |  82.88 $\pm$ 0.21 |  71.32 $\pm$ 0.55 | 80.78 $\pm$ 0.08 | 78.45

##### Compared Methods

We implement various baselines to evaluate the effectiveness of our proposal,
which are categorized as: (1) _Original Teacher_ : The teacher models are used
independently for prediction. We set the probabilities of classes out of the
teacher specialty to zeros. (2) _Ensemble_ : The output logits of teachers are
directly concatenated for predictions over the union label set. (3) Vanilla KD
[6]: The student is trained to mimic the soft targets produced by logits
combination of all teacher models, via minimizing the vanilla KL-divergence
objective. (4) CFL [11]: CFL first maps the hidden representations of the
student and the teachers into a common feature space. The student is trained
by aligning the mapped features to that of the teachers, with supplemental
supervision from the logits combination. We also include a supervised learning
method, which trains the student with labeled data for a better understanding
of the performance. We compare the average accuracy of each method in three
random experiments.

###  5.2 Quantitative Analysis

We compare our proposed method CKA with SOTA on above-mentioned classification
datasets. The experiment results and corresponding model size are listed in
Table 2. Our findings are: (1) Simple baselines can be seriously affected by
incomplete datasets and annotations, showing that it is necessary to conduct
amalgamation. (2) CFL cannot achieves consistent improvements on comprehensive
tasks, demonstrating the instability of supervision based on simple feature
alignments. (3) Our proposed CKA and its variants outperform the previous
baseline models on all the datasets, and the average accuracy of CKA-Inter is
achieves a 1.60 points gain over the best performing baseline model. On the
FGVC-Aircraft dataset, the knowledge consolidation accuracy of CKA reached
80.78% without label information, approaching that of supervised learning
methods. We attribute this success to the fact that CKA provides the student
with natural semantic relevance estimated on the sample set based on
contrastive losses, and the intra-class cohesion and inter-class separation
methods effectively transfer feature-level knowledge. Furthermore, supervisory
contradictions from incomplete teachers are avoided by soft labels at task-
level amalgamation. These promising results indicate that our CKA framework
produces better supervisions for training the student model, yields great
potentials for model reusing.

###  5.3 Ablation Study

We conduct ablation studies to investigate the contribution of the contrastive
losses and soft-target distillation loss described in our proposed approach.

For margin-based intra-model contrastive loss, we compare the performances by
turning them on and off. For inter-model loss between teacher-student pairs,
on the other hand, we define three different distances in Eqn. 4, including
euclidean distance, cosine and MMD distance. We summarize the comparative
results in Table 3, where we observe that the CKA-Inter with MMD distance
yields better performance than others. Moreover, CKA and its variants also
improve with a large room over KD and CFL, validating the complement of
contrastive losses and flexibility of soft-target loss.

Table 3: Ablation analysis of CKA. Remove modules lead to deteriorated performance Method |  Vanilla KDkd | CFL | CKA |  W/O Inter-model loss | W Inter-model loss |  W/O Intra-model loss
---|---|---|---|---|---|---
Euclidean | Cosine | MMD
Cars | 80.22 | 81.12 | 82.88 | 80.04 | 82.33 | 82.95 | 83.21 | 82.33
Aircraft | 78.00 | 79.98 | 80.78 | 77.97 | 80.21 | 80.45 | 81.42 | 79.02
Table 4: Result of merging heterogeneous teachers with different architectures is demonstrated by Stanford Dogs dataset.  Teachers |  $\mathcal{T}_{1}$: _restnet-18_ |  $\mathcal{T}_{2}$: _restnet-34_ |  $\mathcal{T}_{1}$: _restnet-50_ |  $\mathcal{T}_{2}$: _restnet-34_
---|---|---|---|---
Method | Vanilla KD | CFL | CKA | Vanilla KD | CFL | CKA
Student Net | _resnet-34_ | 80.67 | 81.09 | 82.54 | 80.54 | 81.23 | 82.08
_resnet-50_ | 82.04 | 82.25 | 83.18 | 82.62 | 84.55 | 85.21

###  5.4 Results in Challenging Settings

####  5.4.1 CKA with Heterogeneous Teachers

We further consider merging knowledge from heterogeneous teachers with
different structures. Specifically, we random select two different _resnet_
architectures as the teachers, respectively. The results are listed in Table
4. We find that while a larger student tends to perform better, indicating
that the wider and larger the model, the more complete the knowledge can be
learned. Our CKA achieves the best results on the Stanford Dogs, showing its
effectiveness for heterogeneous teachers.

Table 5:  Results of merging from teacher models with different knowledge domains and in a cross-dataset scenario of Stanford Cars and FGVC-Aircraft. Method |  $\mathcal{T}_{1}$: Stanford Cars |  $\mathcal{T}_{2}$: FGVC-Aircraft | Merge
---|---|---|---
Supervised |  89.64 $\pm$ 0.00 |  82.78 $\pm$ 0.00 |  86.90 $\pm$ 0.00
Teacher1 |  89.64 $\pm$ 0.00 | — | —
Teacher2 | — |  78.00 $\pm$ 0.00 | —
Ensemble | — | — |  82.08 $\pm$ 0.54
Vanilla KD |  85.26 $\pm$ 0.25 |  80.31 $\pm$ 0.85 |  83.76 $\pm$ 0.59
CFL |  87.99 $\pm$ 0.73 |  84.22 $\pm$ 0.48 |  86.76 $\pm$ 0.46
CKA-Intra |  88.95 $\pm$ 0.00 |  84.93 $\pm$ 0.58 |  87.75 $\pm$ 0.60
CKA-Inter | 89.48 $\pm$ 0.31 |  84.91 $\pm$ 0.59 |  88.11 $\pm$ 0.79
CKA |  89.28 $\pm$ 0.75 | 85.78 $\pm$ 0.07 | 88.21 $\pm$ 0.50

####  5.4.2 CKA with Heterogeneous Teachers for Cross-Dataset

Specifically, we pretrain distinct-task teacher models on different datasets
separately and then train a student to perform classification over the union
label set of both datasets. The results of merging knowledge from two combined
datasets, Stanford Cars and FGVC-Aircraft are listed in Table 5. _resnet-34_
is adopted for training student in the cross-dataset setting. Our CKA still
outperforms previous baseline models in this settings. Interestingly, we find
that the performance of CKA is superior to all baselines and even to the
results of supervision. We speculate that the reason is that the correlation
between classes in different datasets is weak and the data classification
categories are complex, which is prone to confusion by label supervision
alone. In contrast, our CKA uses contrast loss to compute the distance between
samples, which is more robust and discriminative.

##  6 Conclusion

In this paper, we explore knowledge amalgamation for unsupervised
classification tasks for promoting better model reuse. We present a principled
framework CKA, in which contrastive losses and alignment loss are designed to
enlarge the distance between feature representations of samples from different
categories and decrease that of samples from the same categories, as a self-
supervised way to guide the student to learn discriminative features. Besides,
we present a soft-target distillation loss to efficiently and flexibly
transfer the dark knowledge in the task-level amalgamation. Experiments on
several benchmarks demonstrate our CKA can substantially outperform strong
baselines. More extensive investigations show that CKA is generalizable for
challenging settings, including merging knowledge from heterogeneous teachers,
or even cross-dataset teachers.

## References

  * [1] Chen, T., Kornblith, S., Norouzi, M., Hinton, G.: A simple framework for contrastive learning of visual representations. In: International conference on machine learning. pp. 1597–1607. PMLR, Virtual (2020)
  * [2] Chen, X., He, K.: Exploring simple siamese representation learning. In: Proceedings of the IEEE/CVF conference on computer vision and pattern recognition. pp. 15750–15758. Computer Vision Foundation / IEEE, Virtual (2021)
  * [3] Gretton, A., Borgwardt, K.M., Rasch, M.J., Schölkopf, B., Smola, A.: A kernel two-sample test. The Journal of Machine Learning Research 13(1), 723–773 (2012)
  * [4] Grill, J.B., Strub, F., Altché, F., Tallec, C., Richemond, P., Buchatskaya, E., Doersch, C., Avila Pires, B., Guo, Z., Gheshlaghi Azar, M., et al.: Bootstrap your own latent-a new approach to self-supervised learning. Advances in neural information processing systems 33, 21271–21284 (2020)
  * [5] He, K., Zhang, X., Ren, S., Sun, J.: Deep residual learning for image recognition. In: Proceedings of the IEEE conference on computer vision and pattern recognition. pp. 770–778. IEEE Computer Society, Las Vegas, NV, USA (2016)
  * [6] Hinton, G., Vinyals, O., Dean, J.: Distilling the knowledge in a neural network. arXiv preprint arXiv:1503.02531 (2015)
  * [7] Jiang, W., Synovic, N., Jajal, P., Schorlemmer, T.R., Tewari, A., Pareek, B., Thiruvathukal, G.K., Davis, J.C.: Ptmtorrent: A dataset for mining open-source pre-trained model packages. arXiv preprint arXiv:2303.08934 (2023)
  * [8] Jiang, W., Synovic, N., Sethi, R., Indarapu, A., Hyatt, M., Schorlemmer, T.R., Thiruvathukal, G.K., Davis, J.C.: An empirical study of artifacts and security risks in the pre-trained model supply chain. In: Proceedings of the 2022 ACM Workshop on Software Supply Chain Offensive Research and Ecosystem Defenses. pp. 105–114. ACM, Los Angeles, U.S.A (2022)
  * [9] Khosla, A., Jayadevaprakash, N., Yao, B., Li, F.F.: Novel dataset for fine-grained image categorization: Stanford dogs. In: Proc. CVPR workshop on fine-grained visual categorization (FGVC). vol. 2. Citeseer, Colorado Springs, CO (2011)
  * [10] Krause, J., Stark, M., Deng, J., Fei-Fei, L.: 3d object representations for fine-grained categorization. In: 2013 IEEE International Conference on Computer Vision Workshops. pp. 554–561 (2013). https://doi.org/10.1109/ICCVW.2013.77
  * [11] Luo, S., Wang, X., Fang, G., Hu, Y., Tao, D., Song, M.: Knowledge amalgamation from heterogeneous networks by common feature learning. In: 28th Proceedings of the International Joint Conference on Artificial Intelligence. ijcai.org, Macao, China (2019)
  * [12] Maji, S., Rahtu, E., Kannala, J., Blaschko, M., Vedaldi, A.: Fine-grained visual classification of aircraft. arXiv preprint arXiv:1306.5151 (2013)
  * [13] Shen, C., Wang, X., Song, J., Sun, L., Song, M.: Amalgamating knowledge towards comprehensive classification. In: Proceedings of the AAAI Conference on Artificial Intelligence. pp. 3068–3075. AAAI Press, Honolulu, Hawaii, USA (2019)
  * [14] Vongkulbhisal, J., Vinayavekhin, P., Visentini-Scarzanella, M.: Unifying heterogeneous classifiers with distillation. In: Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition. pp. 3175–3184. Computer Vision Foundation / IEEE, Long Beach, CA, USA (2019)
  * [15] Wah, C., Branson, S., Welinder, P., Perona, P., Belongie, S.: The caltech-ucsd birds-200-2011 dataset (2011)
  * [16] Zhao, B., Cui, Q., Song, R., Qiu, Y., Liang, J.: Decoupled knowledge distillation. In: Proceedings of the IEEE/CVF Conference on computer vision and pattern recognition. pp. 11953–11962. Computer Vision Foundation / IEEE, Virtual (2022)