# Discrete Linear Canonical Transform on Graphs: Uncertainty Principle and
Sampling

Yu Zhanga,b,c Bing-Zhao Lia,b<EMAIL_ADDRESS>aSchool of Mathematics
and Statistics, Beijing Institute of Technology, Beijing 100081, China
bBeijing Key Laboratory on MCAACI, Beijing Institute of Technology, Beijing
100081, China  cDepartment of Mechanical Engineering, Keio University,
Yokohama 223-8522, Japan

###### Abstract

With an increasing influx of classical signal processing methodologies into
the field of graph signal processing, approaches grounded in discrete linear
canonical transform have found application in graph signals. In this paper, we
initially propose the uncertainty principle of the graph linear canonical
transform (GLCT), which is based on a class of graph signals maximally
concentrated in both vertex and graph spectral domains. Subsequently,
leveraging the uncertainty principle, we establish conditions for recovering
bandlimited signals of the GLCT from a subset of samples, thereby formulating
the sampling theory for the GLCT. We elucidate interesting connections between
the uncertainty principle and sampling. Further, by employing sampling set
selection and experimental design sampling strategies, we introduce optimal
sampling operators in the GLCT domain. Finally, we evaluate the performance of
our methods through simulations and numerical experiments across applications.

###### keywords:

Graph signal processing, graph Fourier transform, graph linear canonical
transform, uncertainty principle, sampling theory.

##  1 Introduction

Advancements in information and communication technologies have led to the
integration of various types of data in fields such as social networks,
recommendation systems, medical image processing, and bioinformatics. Unlike
traditional time series or digital images, these datasets often exhibit
complex and irregular structures, posing challenges for conventional signal
processing tools [1, 2]. Over the past decade, significant progress has been
made in the development of tools for analyzing signals defined on graphs,
giving rise to the field of graph signal processing (GSP) [1, 2, 3, 4, 5, 6].

GSP relies on graph topology, extending traditional discrete signal processing
to signals characterized by the complex and irregular underlying structures
represented by corresponding graphs. Emphasizing the relationship,
interaction, and effects between signals and their graph structures in various
real-world scenarios, GSP aims to build a theoretical framework for high-
dimensional and large-scale data analysis tasks from a signal processing
perspective. This includes graph transforms [1, 2, 3, 4], frequency analysis
[7, 8], filtering [9, 10], sampling and interpolation [11, 12, 13, 14, 15, 16,
17], reconstruction and recovery [18, 19, 20], uncertainty principles on
graphs [12, 21, 22], fast computation algorithms [23, 14], and more. In GSP,
spectral analysis plays a crucial role, inspired by harmonic analysis, using
the graph Laplacian operator as its core theoretical foundation, naturally
extending concepts such as frequency and filter banks to the graph domain [1,
3]. On the other hand, inspired by algebraic methods, the multiplication of
the graph signal with the adjacency matrix of the underlying graph yields the
fundamental shift operation for graph signals [2, 4]. These two approaches
bring complementary perspectives and tools, collectively fostering attention
and development in the GSP field. Although the adjacency matrix does not
directly represent signal variations as the Laplacian operator does, its
eigenvectors still reveal the smoothness and variation characteristics of
signals on the graph [7]. In this paper, we adopt an algebraic method
framework to broaden its scope.

A crucial tool in this framework is the graph Fourier transform (GFT), which
expands signals into eigenvectors of the adjacency matrix, defining the
spectrum through corresponding eigenvalues. It is a direct generalization of
the discrete Fourier transform to graphs [4]. Although the GFT can perform
various operations on real-world graph signals, it cannot capture the
transformation process from the vertex domain to the spectral domain and
cannot handle graph signals with chirp-like characteristics [24, 25]. In
recent years, approaches like windowed graph Fourier transform [8], graph
fractional Fourier transform (GFRFT) [24], windowed graph fractional Fourier
transform [26], and the directed graph fractional Fourier transform [27] have
been proposed to address issues with these methods. However, these approaches
still suffer from insufficient degrees of freedom, lack of flexibility, and
underutilization of parameters. In our previous study [28], we extended the
linear canonical transform (LCT) [29] to graphs, introducing the graph linear
canonical transform (GLCT). As the generalization of the Fourier transform and
fractional Fourier transform, the GLCT, with its three free parameters,
exhibits strong adaptability and flexibility in signal processing. It
demonstrates clear advantages in handling non-stationary signals and non-
integer-order cases. We prove that it satisfies all expected properties,
unifying GFT, GFRFT, and graph scale transform, and provide two examples of
GLCT.

In this paper, we further explore the uncertainty principle in the GLCT
domain. The Heisenberg uncertainty principle in continuous-time signals [30],
highlights the fundamental trade-off between the signal’s spread in time and
its spread in frequency. Specifically, continuous-time signals cannot be
perfectly localized in both time and frequency domains. In [21], the
uncertainty principle for graph signals was first introduced using geodesic
distance, aiming to establish a connection between the signal spread at the
vertices of the graph and its spectrum spread defined by GFT in the graph
spectral domain. However, a key difference between graph signals and time
signals is that, while time or frequency has a well-defined distance concept,
graphs are not metric spaces. Additionally, the vertices of a graph can
represent more complex signals, such as those in Hilbert spaces [31]. To
overcome the complexity of defining distances, M. Tsitsvero, and others [12]
creatively introduced prolate spheroidal wave functions [32, 33] into graph
signals, defining a new graph uncertainty principle using the energy
percentage. We leverage this approach to derive the uncertainty principle in
the GLCT domain. Leveraging this approach, we delineate the uncertainty
principle in the GLCT domain, which proves to be broader in scope compared to
those in the GFT domain. As the parameters of the GLCT undergo modifications,
the scope of the uncertainty principle correspondingly shifts. This dynamic
interplay furnishes pivotal insights for the formulation of the sampling
theorem in the GLCT domain.

Furthermore, based on the proposed uncertainty principles of GLCT and its
localization properties, we consider sampling and reconstruction in the GLCT
domain. In [11], it was demonstrated that bandlimited graph signals can be
perfectly recovered with a sufficient number of samples. There is an
intriguing connection between sampling theory and uncertainty principles in
the GLCT domain, leading to conditions for recovering bandlimited signals from
subsets of their values. These conditions are found to be equivalent to those
derived in [11, 12, 13, 34]. We present several signal recovery algorithms and
sampling strategies aimed at finding the optimal sampling set. Building on the
relationship between the GFT, the GFRFT and, the GLCT, this paper establishes
a framework for graph signal sampling based on prior information in the GLCT
domain. The sampling process can be viewed as a filtering operation in the
GLCT domain, and the choice of filters allows the reconstruction of graph
signals. In practice, our introduced sampling approaches in the GLCT domain
not only address the limitations of the GFT in handling signals with chirp-
like features but also demonstrates robustness in processing a variety of
real-world signals. By judicious selection of optimal GLCT parameters, we can
enhance performance in signal reconstruction, filtering, classification, and
clustering tasks.

Our contributions are summarized as follows:

  * 1.

Firstly, based on conditions for perfect localization in the vertex and graph
spectral domains, we propose an uncertainty principle in the GLCT domain and
discuss its intriguing connection with sampling.

  * 2.

Secondly, by defining bandlimited signals in the GLCT domain, we propose
conditions for sampling and perfect recovery in the GLCT domain. We present a
method for designing and selecting optimal sampling operators from qualified
sampling operators, comparing various sampling strategies.

  * 3.

Finally, through simulation experiments, we compare our proposed sampling
framework with GFT and GFRFT approaches. The results demonstrate that our
sampling framework produces minimal errors and showcases its competitive
performance in the application of semi-supervised classification for online
blogs and clustering of IEEE 118 bus test cases.

This paper is structured as follows. Section 2 provides a brief overview of
our prior work in [28], the graph uncertainty principle, and the graph
sampling framework. In Section 3, based on conditions for perfect
localization, we propose uncertainty principles in the GLCT domain and discuss
their relationship with sampling. In Section 4, we begin by defining
$\mathbf{M}$-bandlimited signals, introduce conditions for sampling and
perfect recovery in the GLCT domain, and select the optimal sampling operator
from qualified sampling operators. We present various sampling strategies and
compare them. Section 5 evaluates the classification performance of GLCT
sampling on online blogs and clustering of IEEE 118 bus test cases, comparing
it with GFT and GFRFT sampling. Finally, Section 6 concludes this paper.

##  2 Preliminaries

In this section, we briefly review concepts in discrete signal processing on
GLCT [28] and introduce the uncertainty principle of graph signals [21] and
the basic theory of sampling [11, 12, 13].

###  2.1 Graph Linear Canonical Transform

In discrete signal processing on graphs, the high-dimensional structure of a
signal is represented by a graph $\mathcal{G}=(\mathcal{V},\mathbf{A}),$ where
$\mathcal{V}=\\{v_{0},\dots,v_{N-1}\\}$ denotes the set of vertices and
$\mathbf{A}\in\mathbb{C}^{N\times N}$ denotes the weighted adjacency matrix.
The graph signal is defined as a mapping that maps the vertex $v_{n}$ to a
signal coefficient $x_{n}\in\mathbb{C}$, which can also be written as a
complex vector

$\mathbf{x}=\begin{bmatrix}x_{0}&x_{1}&\ldots&x_{N-1}\end{bmatrix}^{\top}\in\mathbb{C}^{N}.$

The GFT can be defined through the eigendecomposition of either the adjacency
matrix or the Laplacian matrix [1, 2, 3, 4]. Given that the discrete Fourier
transform is equivalent to the GFT in the context of cyclic graphs (directed
graphs) [2, 4], we opt for the adjacency matrix to broaden the scope of
application and derive the equivalent transformation of the discrete linear
canonical transform (DLCT) [35] on graphs. Hence, the adjacency matrix

$\mathbf{A=V\Lambda_{A}V}^{-1}=\sum^{N-1}_{i=0}\lambda_{\mathbf{A}_{i}}u_{i}u^{\ast}_{i},$
(1)

where the columns $u_{0},u_{1},\dots,u_{N-1}$ of $\mathbf{V}$ are the
eigenvectors of $\mathbf{A}$, and the eigenvalues in the corresponding
diagonal matrix $\mathbf{\Lambda_{A}}\in\mathbb{C}^{N\times N}$ are
$\lambda_{\mathbf{A}_{0}},\dots,\lambda_{\mathbf{A}_{N-1}}$. Consequently, the
GFT is formulated as $\mathbf{\hat{x}}=\mathbf{V}^{-1}\mathbf{x}$.

Analogous to how the GFT extends the DFT, we propose the extension of the DLCT
[35] to GSP, introducing the GLCT. The definition of GLCT can be divided into
three stages: the GFRFT, graph scaling, and graph chirp multiplication (CM)
[28].

  1. 1.

GFRFT Stage: By performing eigendecomposition on the orthogonal and
diagonalized GFT matrix $\mathbf{V}^{-1}$, the GFRFT operator can be obtained
as [24]

$\mathbf{F^{\alpha}}=\mathbf{Q}\mathbf{\Lambda}^{\alpha}\mathbf{Q}^{\top},$
(2)

where the columns $q_{0},\dots,q_{N-1}$ of $\mathbf{Q}$ are the eigenvectors
of $\mathbf{V}^{-1}$, and the eigenvalues in the corresponding diagonal matrix
$\mathbf{\Lambda}$ are $\lambda_{0},\dots,\lambda_{N-1}$. And
$\mathbf{\Lambda}^{\alpha}=\mathrm{diag}(\lambda^{\alpha}_{i}),i=0,\dots,N-1$.

  2. 2.

Graph Scaling Stage: We consider the scaled adjacency matrix
$\mathbf{S}=\frac{1}{\beta}\mathbf{A}$ as the graph shift operator [9].
Similar to the GFT, by performing eigendecomposition on $\mathbf{S}$, we
obtain the scaled GFT matrix $\mathbf{V}^{-1}_{\beta}$, from which further
eigendecomposition yields

$\mathbf{V}^{-1}_{\beta}=\mathbf{Q_{\beta}\Lambda_{S}Q^{\top}_{\beta}},$ (3)

where $\beta$ is the scaling parameter, and $\mathbf{Q}_{\beta}$ and
$\mathbf{\Lambda_{S}}$ denote the eigenvectors of $\mathbf{V}^{-1}_{\beta}$
and its corresponding eigenvalues, respectively.

  3. 3.

Graph CM Stage: Similar to the definition of the discrete chirp-Fourier
transform [36], the CM is performed on the GFT matrix $\mathbf{V}^{-1}$, as

$\mathbf{F^{\xi}}=\mathbf{Q}\mathbf{\Lambda}^{\xi}\mathbf{Q}^{\top},$ (4)

where $\xi=lk+f$ is the CM parameter which is a linear function of $k$, with
$l$ and $f$ being constants, and the matrix $\mathbf{\Lambda}^{\xi}$ is
$\mathrm{diag}(\lambda^{\xi}_{k})=\mathrm{diag}(\exp\left(j\frac{\pi}{2}k\xi\right))$.

Based on the three stages described above, according to Eqs (2)-(4), we can
derive the eigendecomposition of the GLCT as depicted in Figure 1. Its
definition is as follows.

Figure 1: Eigendecomposition of the GLCT [28].

######  Definition 1.

Parallel to the GFT, the GLCT of $\mathbf{x}\in\mathbb{C}^{N}$ can be defined
as [28]

$\mathbf{\hat{x}}=\mathbf{O}^{\mathbf{M}}\mathbf{x}=\mathbf{\Lambda}^{\xi}\mathbf{Q}_{\beta}\mathbf{\Lambda}^{\alpha}\mathbf{Q}^{\top}\mathbf{x},$
(5)

where the operator $\mathbf{O}^{\mathbf{M}}$ is a operator of the GLCT, and
the matrix $\mathbf{M}$ entries are $[a,b;c,d]$, which is called an GLCT
parameter matrix, and $ad-bc=1,\ a,b,c,d\in\mathbb{R}$. The parameter
relations between $\mathbf{M}$ and $(\xi,\beta,\alpha)$ are

$\xi=\frac{ac+bd}{a^{2}+b^{2}},\ \ \beta=\sqrt{a^{2}+b^{2}},\ \
\alpha=\cos^{-1}\left(\frac{a}{\beta}\right)=\sin^{-1}\left(\frac{b}{\beta}\right).$
(6)

Furthermore, its inverse GLCT (IGLCT) is

$\mathbf{x}=\mathbf{O}^{-\mathbf{M}}\mathbf{\hat{x}}=\mathbf{Q}\mathbf{\Lambda}^{-\alpha}\mathbf{Q}^{\top}_{\beta}\mathbf{\Lambda}^{-\xi}\mathbf{\hat{x}}.$
(7)

######  Remark 1.

In particular, when $a=\cos\alpha,b=\sin\alpha,c=-\sin\alpha,$ and
$d=\cos\alpha$, the GLCT operator will be the GFRFT operator [24], when
$a=0,b=1,c=-1,$ and $d=0$, the GLCT operator will be the GFT operator.

######  Remark 2.

Traditional GFT-based spectral analysis, limited by the Laplacian matrix,
overlooks signals’ vertex-frequency geometry. The GLCT offers flexibility in
shaping spectral characteristics, enabling signal compression, expansion, and
rotation in the vertex-frequency plane to form different smoothness patterns.

###  2.2 Graph Uncertainty Principle

A fundamental property of signals is the Heisenberg uncertainty principle
[30], which states that there is a fundamental trade-off between the spread of
a signal in time and the spread of its spectrum in frequency. The traditional
uncertainty principle of time signals is

$\Delta t^{2}\Delta w^{2}\geq\frac{1}{4},$

where $\Delta t^{2}$ and $\Delta w^{2}$ denote the time and frequency spreads,
respectively, of the classical continuous-time signal $x(t)$ and its Fourier
domain counterpart $\hat{x}(w)$.

After the GFT was proposed, the signal uncertainty principle defined on
undirected connected graphs was first derived in [21]. The spreads of vector
$\mathbf{x}$ in the vertex domain and the GFT domain are respectively defined
as

$\begin{cases}\Delta
g^{2}:=\min_{u_{0}\in\mathcal{V}}\frac{1}{\left|\left|\mathbf{x}\right|\right|^{2}}\mathbf{x}^{\top}\mathbf{P}^{2}_{u_{0}}\mathbf{x},&\\\
\Delta
s^{2}:=\frac{1}{\left|\left|\mathbf{x}\right|\right|^{2}}\sum_{i}\lambda_{\mathbf{A}_{i}}\left|\hat{\mathbf{x}}_{i}\right|^{2},&\end{cases}$
(8)

where
$\mathbf{P}_{u_{0}}:=\mathrm{diag}(d(u_{0},v_{1}),d(u_{0},v_{2}),...,d(u_{0},v_{N}))$,
and $d(u,v)$ denotes the geodesic distance between nodes $u$ and $v$. Thus,
the graph uncertainty principle can be represented by the admissible region of
$\Delta g^{2}$ and $\Delta s^{2}$, i.e.,

$\Gamma_{u_{0}}=\left\\{(s,g):\Delta s^{2}(\mathbf{x})=s,\Delta
g^{2}_{u_{0}}(\mathbf{x})=g\right\\},$

where $\Delta
g^{2}_{u_{0}}\left(\mathbf{x}\right)=\frac{1}{\left|\left|\mathbf{x}\right|\right|^{2}}\mathbf{x}^{\top}\mathbf{P}^{2}_{u_{0}}\mathbf{x}$,
and the curve can be expressed as

$\gamma_{u_{0}}=\min_{\mathbf{x}}\Delta
g^{2}_{u_{0}}(\mathbf{x})~{}\text{subject to}~{}\Delta s^{2}(\mathbf{x})=s.$

The above uncertainty principle studies the trade-off between the graph signal
distribution and its spectral domain based on the specific definition of graph
distances. However, a key distinction between the graph signal and the time
signal is that time exists in a metric space with a well-defined distance
concept, whereas the graph is not a metric space, and, its vertices can be
assigned mathematical objects endowed with a richer structure, rather than
assigning real numbers to each vertex. An example of such mathematical objects
originates from Hilbert spaces [31]. In this case, defining distances between
vertices is not straightforward.

To overcome the limitations associated with graph distance definitions, this
paper adopts an alternative definition based on spread [12]. We derive an
uncertainty principle for the GLCT domain that does not require any additional
distance definitions.

This method characterizes temporal spread based on the percentage of energy
falling within the interval $[-T/2,T/2]$, defined as

$\frac{\int^{T/2}_{-T/2}\left|x\left(t\right)\right|^{2}\mathrm{d}t}{\int^{+\infty}_{-\infty}\left|x\left(t\right)\right|^{2}\mathrm{d}t}=\zeta^{2}.$
(9)

Similarly, the definition of frequency domain spread is given by

$\frac{\int^{W/2}_{-W/2}\left|\hat{x}\left(w\right)\right|^{2}\mathrm{d}w}{\int^{+\infty}_{-\infty}\left|\hat{x}\left(w\right)\right|^{2}\mathrm{d}w}=\eta^{2}.$
(10)

Using these two defined spreads, we can further study the uncertainty
principle in the GLCT domain in Section 3.

###  2.3 Sampling on Graph Signals

Sampling is one of the fundamental challenges in GSP. Its objective is to
identify conditions for recovering bandlimited graph signals from a subset of
values and to devise suitable sampling and recovery strategies [12].

Consider the sampling set
$\mathcal{S}=\left(\mathcal{S}_{0},\dots,\mathcal{S}_{|\mathcal{S}|-1}\right),\mathcal{S}_{i}\in\\{0,1,\dots,N-1\\}$,
its dimension coefficient is defined as $|\mathcal{S}|$. This coefficient is
the size of the sampling indexes for obtaining a sampled signal
$\mathbf{x}_{\mathcal{S}}\in\mathbb{C}^{|\mathcal{S}|}(|\mathcal{S}|<N)$ from
a bandlimited graph signal $\mathbf{x}\in\mathbb{C}^{N}$. The sampling
operator $\mathbf{D}$ is defined as a linear mapping from $\mathbb{C}^{N}$ to
$\mathbb{C}^{|\mathcal{S}|}$ [11] , which is expressed as

$\mathbf{D}_{i,j}=\left\\{\begin{array}[]{rl}1,&j=\mathcal{S}_{i},\\\
0,&\mathrm{otherwise}.\end{array}\right.$ (11)

We then recover $\mathbf{x}$ from $\mathbf{x}_{\mathcal{S}}$ with
interpolation operator $\mathbf{R}$, which is a linear mapping from
$\mathbb{C}^{|\mathcal{S}|}$ to $\mathbb{C}^{N}$.

$\displaystyle\begin{array}[]{rl}\mathrm{sampling:}&\mathbf{x}_{\mathcal{S}}=\mathbf{D}\mathbf{x}\in\mathbb{C}^{|\mathcal{S}|},\\\
\mathrm{interpolation:}&\mathbf{x}_{\mathcal{R}}=\mathbf{Rx}_{\mathcal{S}}=\mathbf{RDx}\in\mathbb{C}^{N},\end{array}$

where $\mathbf{x}_{\mathcal{R}}\in\mathbb{C}^{N}$ recovers $\mathbf{x}$ either
approximately or exactly.

##  3 Uncertainty Principle of the GLCT

In the realm of continuous-time signals, Landau, Pollak, and Slepian
extensively explored space-frequency analysis linked to projection operators.
Their groundbreaking work, spanning the 1960s, encompassed the uncertainty
principle, eigenvalue distributions, and prolate ellipsoid wave functions [32,
33]. By extending this uncertainty principle to graph signals [12], we
circumvent the need for any extra distance definition, mitigating the
drawbacks associated with graph distance definitions. Furthermore, we apply
this methodology to establish the uncertainty principle in the GLCT domain.

###  3.1 Localization in Vertex and Spectral Domains

Before defining the uncertainty principle, we first introduce two localization
(bandlimited) operators [37]. For a subset of vertices
$\mathcal{S}\subseteq\mathcal{V}$, the vertex-limiting operator is defined as
follows

$\mathbf{D}_{\mathcal{S}}=\mathrm{diag}\left(d_{1},\ldots,d_{i},\ldots,d_{N}\right),$
(12)

where, $d_{i}=1$, if $i\in\mathcal{S}$, and $d_{i}=0$ if
$i\notin\mathcal{S},i=1,\ldots,N$. For this operator, it is the corresponding
sampling operator $\mathbf{D}_{\mathcal{S}}$ in (11).

Similarly, given the unitary matrix $\mathbf{O}^{\mathbf{M}}$ used in (5) and
the index subset $\mathcal{F}\subseteq\mathcal{\hat{G}}$, where
$\mathcal{\hat{G}}=\left\\{{1,...,N}\right\\}$ represents the set of all graph
spectral indices, we introduce the spectral-limiting operator as

$\mathbf{B}^{\mathbf{M}}_{\mathcal{F}}=\mathbf{O}^{-\mathbf{M}}\mathbf{\Sigma}_{\mathcal{F}}\mathbf{O}^{\mathbf{M}},$
(13)

where, $\mathbf{\Sigma}_{\mathcal{F}}$ is a diagonal matrix whose definition
is consistent with $\mathbf{D}_{\mathcal{S}}$, i.e. $\mathbf{\Sigma}_{ii}=1$,
if $i\in\mathcal{F}$, and $\mathbf{\Sigma}_{ii}=0$ if $i\notin\mathcal{F}$.

According to the properties of diagonal matrices, it is evident that the
operators $\mathbf{D}_{\mathcal{S}}$ and
$\mathbf{B}^{\mathbf{M}}_{\mathcal{F}}$ are symmetric and positive semi-
definite, with a spectral norm of exactly 1 for both operators. To simplify
notation without loss of generality, we use $\mathbf{D}$ and
$\mathbf{B}^{\mathbf{M}}$ to represent these two operators, respectively.

Therefore, the vector $\mathbf{x}$ is perfectly localized (bandlimited) over
$\mathcal{S}$ if

$\mathbf{Dx=x},$ (14)

and perfectly localized (bandlimited) over $\mathcal{F}$ if

$\mathbf{B}^{\mathbf{M}}\mathbf{x=x}.$ (15)

Here we have a lemma [12] about perfect localization.

######  Lemma 1.

A vector $\mathbf{x}$ is perfectly localized over both the vertex set
$\mathcal{S}$ and the spectral set $\mathcal{F}$ if and only if

$\lambda_{\max}\left(\mathbf{B}^{\mathbf{M}}\mathbf{DB}^{\mathbf{M}}\right)=1.$
(16)

In such a case, $\mathbf{x}$ is the eigenvector associated with the unit
eigenvalue.

###### Proof.

The proof is reported in A. ∎

Equivalently, the perfect localization over sets $\mathcal{S}$ and
$\mathcal{F}$ can also be achieved if and only if

$\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{D}\right|\right|_{2}=\left|\left|\mathbf{DB}^{\mathbf{M}}\right|\right|_{2}=1.$
(17)

###  3.2 Uncertainty Principle of the GLCT

Using the energy expectations of equations (9) and (10), we can extend them to
the GLCT domain similar to [12, 22].

Suppose a vector $\mathbf{x}$, vertex subset $\mathcal{S}$ and graph spectral
subset $\mathcal{F}$, two subsets are used (12) and (13). The vectors
$\mathbf{Dx}$ and $\mathbf{B}^{\mathbf{M}}\mathbf{x}$ represent the projection
of $\mathbf{x}$ on the vertex set $\mathcal{S}$ and the spectral set
$\mathcal{F}$ respectively. Then, similar to equations (9) and (10),
$\zeta^{2}$ and $(\eta^{\mathbf{M}})^{2}$ are used to denote the energy
expectations falling within the sets $\mathcal{S}$ and $\mathcal{F}$
respectively, as follows

$\frac{\left|\left|\mathbf{Dx}\right|\right|^{2}_{2}}{\left|\left|\mathbf{x}\right|\right|^{2}_{2}}=\zeta^{2},\frac{\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|^{2}_{2}}{\left|\left|\mathbf{x}\right|\right|^{2}_{2}}=(\eta^{\mathbf{M}})^{2}.$
(18)

In this paper, we find closed form region boundaries for all admissible pairs
$\phi(\zeta,\eta^{\mathbf{M}})$ by generalizing the uncertainty principle of
the GFT to the case of the GLCT. In (18), the graph topology is captured by
the matrix $\mathbf{O}^{\mathbf{M}}$, which appears in the definition of the
GLCT in (5), inside the operator $\mathbf{B}^{\mathbf{M}}$. Our goal is to
establish a trade-off between $\zeta$ and $\eta^{\mathbf{M}}$ by localizing
operators $\mathbf{D}$ and $\mathbf{B}^{\mathbf{M}}$, and find a signal that
obtains all admissible pairs to prove the uncertainty principle.

We begin with the first uncertainty relation for the localization operators
$\mathbf{D}$ and $\mathbf{B^{M}}$, where localization depends on the maximal
eigenvalue $\lambda_{\max}$ of the operator $\mathbf{B^{M}DB^{M}}$ [22]. That
is,
$\left\|\mathbf{B^{M}DB^{M}}\right\|_{2}=\lambda_{\max}\left(\mathbf{B}^{\mathbf{M}}\mathbf{DB}^{\mathbf{M}}\right)<1$.
It can be thought that the admissible region $\Gamma$ provides clear bounds
for specifying this uncertainty relation on $\mathcal{G}$.

######  Lemma 2.

Suppose for a vector $\mathbf{x}$, there are
$\left|\left|\mathbf{x}\right|\right|_{2}=1,\left|\left|\mathbf{Dx}\right|\right|_{2}=\zeta,\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}=\eta^{\mathbf{M}}$,
and $\lambda_{\max}<1$ and if $\lambda_{\max}\leq\zeta^{2}\eta^{2}$, there is
an inequality

$\arccos\zeta+\arccos\eta^{\mathbf{M}}\geq\arccos\sqrt{\lambda_{\max}\left(\mathbf{B^{M}DB^{M}}\right)},$
(19)

and finally reached the upper bound

$\eta^{\mathbf{M}}\leq\zeta\sqrt{\lambda_{\max}}+\sqrt{\left(1-\zeta^{2}\right)\left(1-\lambda_{\max}\right)}.$
(20)

###### Proof.

The proof is reported in B. ∎

Lemma 2 provides a general restriction on the set
$\phi(\zeta,\eta^{\mathbf{M}})$ in the upper right corner of the unit square
$[\lambda_{\max},1]^{2}$. By simple generalization extension, we obtain
similar results for the other three corners. For a subset $\mathcal{S}$, we
denote its complement as $\bar{\mathcal{S}}$ such that
$\mathcal{V}=\mathcal{S}\cup\bar{\mathcal{S}}$ and
$\mathcal{S}\cap\bar{\mathcal{S}}=\emptyset$. Accordingly, we define the
vertex projection on $\bar{\mathcal{S}}$ as $\bar{\mathbf{D}}$. Similarly, the
projection on the complement $\bar{\mathcal{F}}$ is denoted by
$\bar{\mathbf{B}}^{\mathbf{M}}$. Considering all subdomains of the square
$[0,1]^{2}$, we obtain the uncertainty principle for the GLCT.

######  Theorem 1.

There exists a vector $\mathbf{x}$ such that
$\left|\left|\mathbf{x}\right|\right|_{2}=1,\left|\left|\mathbf{Dx}\right|\right|_{2}=\zeta,\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}=\eta^{\mathbf{M}}$,
thus the admissible pair $\phi(\zeta,\eta^{\mathbf{M}})\in\Gamma$ is obtained,
where

$\displaystyle\Gamma$
$\displaystyle=\biggl{\\{}\phi\left(\zeta,\eta^{\mathbf{M}}\right):$ (21)
$\displaystyle\arccos\zeta+\arccos\eta^{\mathbf{M}}\geq\arccos\sqrt{\lambda_{\max}\left(\mathbf{B^{M}DB^{M}}\right)},$
$\displaystyle\arccos\sqrt{1-\zeta^{2}}+\arccos\eta^{\mathbf{M}}\geq\arccos\sqrt{\lambda_{\max}\left(\mathbf{B^{M}}\bar{\mathbf{D}}\mathbf{B^{M}}\right)},$
$\displaystyle\arccos\zeta+\arccos\sqrt{1-(\eta^{\mathbf{M}})^{2}}\geq\arccos\sqrt{\lambda_{\max}\left(\bar{\mathbf{B}}^{\mathbf{M}}\mathbf{D}\bar{\mathbf{B}}^{\mathbf{M}}\right)},$
$\displaystyle\arccos\sqrt{1-\zeta^{2}}+\arccos\sqrt{1-(\eta^{\mathbf{M}})^{2}}\geq\arccos\sqrt{\lambda_{\max}\left(\bar{\mathbf{B}}^{\mathbf{M}}\bar{\mathbf{D}}\bar{\mathbf{B}}^{\mathbf{M}}\right)}.\biggr{\\}}$

###### Proof.

The proof follows a similar structure to that of Lemma 2, demonstrating the
inequalities for the remaining three corners. ∎

Fig. 2 illustrates the uncertainty principle of Theorem 1, highlighting some
noteworthy characteristics. The blue portion in the figure corresponds to the
range of $\Gamma$ in the GFT domain, while in the GLCT domain, the ranges in
the lower left and upper right respectively tend to $(0,0)$ and $(1,1)$ as the
parameter continuously changes, as indicated by the additional red portion. If
we set the inequalities given by Theorem 1 to equality, we can determine the
equations of the curve at the four corners of Fig. 2, namely the upper right,
upper left, lower right, and lower left corners. Notably, the upper right
corner of the admissible region $\Gamma$ specifies the domain where the
maximum energy concentration occurs in the vertex and graph spectral for
$\phi(\zeta,\eta^{\mathbf{M}})$. The equation for this curve is obtained

$\eta^{\mathbf{M}}=\sqrt{\lambda_{\max}\left(\mathbf{B^{M}DB^{M}}\right)}\zeta+\sqrt{1-\lambda_{\max}\left(\mathbf{B^{M}DB^{M}}\right)}\sqrt{1-\zeta^{2}}.$
(22)

When the sets $\mathcal{S}$ and $\mathcal{F}$ yield matrices $\mathbf{D}$ and
$\mathbf{B}^{\mathbf{M}}$ satisfying the perfect localization condition in
Lemma 1, the curve collapses to a single point, precisely the upper right
corner. Consequently, similar to [12, 22], any curve originating from a corner
of the admissible region $\Gamma$ in Fig. 2 can be mapped onto the
corresponding corner, provided that the associated operator meets the
criterion for perfect localization.

Figure 2: The admissible region $\Gamma$ in Theorem 1.

###  3.3 The Connection between the Uncertainty Principle and Sampling

For the proposed uncertainty principle and as described in [22], there are
intriguing connections with the sampling conditions. Examining the upper left
corner of the admissible region in Fig. 2, it is evident that if the signal is
perfectly localized in the graph spectral subset $\mathcal{F}$, then
$\left(\eta^{\mathbf{M}}\right)^{2}=1$. However, for the vertex domain subset
$\mathcal{S}$, we require $\zeta^{2}\neq 0$ to avoid perfect localization of
the signal on the complement of $\mathcal{S}$, making it unrecoverable. This
is analogous to the perfect recovery conditions given in Section 4.

If we allow for some energy dissipation in the graph spectral domain, i.e.,
$\left(\eta^{\mathbf{M}}\right)^{2}<1$, to ensure $\zeta^{2}\neq 0$, we need
to check

$1-\lambda_{\max}\left(\mathbf{B^{M}\bar{D}B^{M}}\right)>\left(\eta^{\mathbf{M}}\right)^{2}.$

This form of condition is highly useful for designing potential sampling
strategies. Specifically, the largest eigenvalue,
$\lambda_{\max}\left(\mathbf{B^{M}\bar{D}B^{M}}\right)$, of the operator
$\mathbf{B^{M}\bar{D}B^{M}}$ is equivalent to the square of the maximum
singular value, $\sigma^{2}_{\max}\left(\mathbf{\bar{D}B^{M}}\right)$, of the
operator $\mathbf{\bar{D}B^{M}}$. This equivalence suggests the use of
$\sigma^{2}_{\max}\left(\mathbf{\bar{D}B^{M}}\right)$ or
$\lambda_{\max}\left(\mathbf{B^{M}\bar{D}B^{M}}\right)$ as the objective
function for minimization, as exemplified in Table 1. Specifically,
$\sigma^{2}_{\max}\left(\mathbf{\bar{D}B^{M}}\right)$ also serves as a
criterion for perfect signal recovery, requiring that
$\sigma^{2}_{\max}\left(\mathbf{\bar{D}B^{M}}\right)<1$ always holds. These
aspects will be further explored in Section 4.

##  4 Sampling Theory of the GLCT

In this section, we present a novel sampling theory for the GLCT, drawing
inspiration from the established sampling theories of the GFT and the GFRFT
[11, 12, 13, 14, 15, 34]. Our approach to graph sampling, rooted in the graph
linear canonical basis and graph frequency-based sampling methods, tackles the
challenge of reconstructing bandlimited graph signals from their subsampled
counterparts.

An $\mathbf{M}$-bandlimited graph signal is characterized as a signal with
zero graph linear canonical coefficients aligned with the eigenvalues of
$\lambda_{|\mathcal{F}|-1}$. This definition implies that $f(i)=0$ for
$\lambda,i\geq|\mathcal{F}|$, where $\mathcal{F}$ denotes the set of indices
associated with non-zero graph linear canonical coefficients.

Specifically, the indices in $\mathcal{F}$ identify those components of the
graph linear canonical coefficients that contribute significantly to the
signal. We further elaborate on the concept of $\mathbf{M}$-bandlimited graph
signals in the subsequent discussion.

###  4.1 $\mathbf{M}$-Bandlimited Graph Signals

In the context of GSP, the notion of $\mathbf{M}$-bandlimited signals becomes
imperative, necessitating adherence to the previously introduced conditions of
perfect localization. Assuming an arbitrary signal $\mathbf{y}$ lacks
bandlimited characteristics, preprocessing is required in the form of
bandlimiting.

This involves the utilization of a graph spectral domain filter
$\mathbf{\Sigma}_{\mathcal{F}}$, defined analogously to Eq. (13). First, the
GLCT is applied to the signal $\mathbf{y}$, yielding
$\mathbf{\hat{y}}=\mathbf{O}^{\mathbf{M}}\mathbf{y}$. Subsequently,
$\mathbf{\hat{y}}$ undergoes filtering

$\mathbf{\hat{x}}=\mathbf{\Sigma}_{\mathcal{F}}\mathbf{\hat{y}}=\mathrm{diag}\left(1,...,1\right)\mathbf{\hat{y}}.$

Following this, an IGLCT is applied to $\mathbf{\hat{x}}$, resulting in
$\mathbf{x}=\mathbf{O}^{-\mathbf{M}}\mathbf{\hat{x}}$. Consolidating these
steps, we arrive at

$\mathbf{x}=\mathbf{O}^{-\mathbf{M}}\mathbf{\Sigma}_{\mathcal{F}}\mathbf{O}^{\mathbf{M}}\mathbf{y},$
(23)

where, it is discernible that
$\mathbf{O}^{-\mathbf{M}}\mathbf{\Sigma}_{\mathcal{F}}\mathbf{O}^{\mathbf{M}}=\mathbf{B^{M}}$.
Through this methodology, a $\mathbf{M}$-bandlimited signal $\mathbf{x}$ is
attained, satisfying the conditions of perfect localization as delineated in
Eqs. (14) and (15). Therefore, we introduce the following definition for the
$\mathbf{M}$-bandlimited set.

######  Definition 2.

The closed subspace of graph signals in $\mathbb{C}^{N}$ with a bandwidth of
at most $|\mathcal{F}|$ is denoted as
$\mathrm{BL}_{|\mathcal{F}|}(\mathbf{O}^{\mathbf{M}})$, where
$\mathbf{O}^{\mathbf{M}}$ is defined as in Eq. (5).

###  4.2 Sampling and Perfect Recovery

We initiate our exploration by addressing the fundamental challenge of
determining conditions and methodologies for the perfect recovery of
$\mathbf{x}$ from the sampled graph signal $\mathbf{x}_{\mathcal{S}}$. Drawing
upon certain definitions of graph signal sampling outlined in the
preliminaries, we present the following theorem, mirroring the principles
introduced in [11].

######  Theorem 2.

Let $\mathbf{D}\in\mathbb{C}^{|\mathcal{S}|\times N}$ denote the sampling
operator. The interpolation (recovery) operator
$\mathbf{R}\in\mathbb{C}^{N\times|\mathcal{S}|}$ is characterized by the
following conditions: (1) $\mathbf{R}$ spans the space
$\mathrm{BL}_{|\mathcal{F}|}(\mathbf{O}^{\mathbf{M}})$; (2) $\mathbf{RD}$
functions as a projection operator, achieving perfect recovery

$\mathbf{x}=\mathbf{RDx}=\mathbf{Rx}_{\mathcal{S}},\ \ \text{for all}\
\mathbf{x}\in\mathrm{BL}_{|\mathcal{F}|}\left(\mathbf{O}^{\mathbf{M}}\right),$
(24)

where $\mathbf{x}_{\mathcal{S}}$ denotes the sampled graph signal.

If the original signal $\mathbf{x}$ can be perfectly recovered, then

$\mathbf{x-RDx=x-R\left(I-\bar{D}\right)x=x-R}\left(\mathbf{I-\bar{D}}\mathbf{B}^{\mathbf{M}}\right)\mathbf{x}.$

Hence, the matrix $\mathbf{R=\left(I-\bar{D}B^{M}\right)}^{-1}$. The perfect
recoverability of $\mathbf{x}$ implies the existence of $\mathbf{R}$, which is
tantamount to the invertibility of $\mathbf{\left(I-\bar{D}B^{M}\right)}$.
Consequently, we infer that when
$\left|\left|\mathbf{\bar{D}B^{M}}\right|\right|<1$, the original signal can
be perfectly recovered, as reported in [12]. Herein, the maximum norm of
$\mathbf{\bar{D}B^{M}}$ being $1$ restricts our consideration to this
scenario, where perfect localization on the complement of the sampling set
prevails, rendering signal recovery unattainable. This is closely related to
the uncertainty principles discussed in Section 3. Furthermore, since
$\mathbf{x}\in\mathrm{BL}_{|\mathcal{F}|}\left(\mathbf{O}^{\mathbf{M}}\right)$,

$\mathbf{\left(I-\bar{D}B^{M}\right)x=\left(I-\bar{D}\right)x=Dx=DB^{M}x}.$

By analyzing the preceding equation, it is evident that
$\mathbf{\left(I-\bar{D}B^{M}\right)}$ is equivalent to $\mathbf{DB^{M}}$,
leading to the conclusion that $\mathbf{R=\left(DB^{M}\right)}^{-1}$. However,
in the case of a degenerate matrix, the rank of $\mathbf{DB^{M}}$ may be less
than full. In such instances, we employ the pseudo-inverse to define
$\mathbf{R=\left(DB^{M}\right)}^{{\dagger}}$.

We have the following theorem regarding sampling and perfect recovery.

######  Theorem 3.

For $\mathbf{D}$ such that

$\mathrm{rank}\left(\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|}\right)=|\mathcal{F}|,\quad|\mathcal{F}|<N,$
(25)

where $\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|}$ denotes the first
$|\mathcal{F}|$ columns of $\mathbf{DO}^{-\mathbf{M}}$. For all
$\mathbf{x}\in\mathrm{BL}_{|\mathcal{F}|}(\mathbf{O^{M}})$, perfect recovery,
$\mathbf{x=RDx}$, is achieved by choosing

$\mathbf{R=\left(DB^{M}\right)^{{\dagger}}}=\left(\mathbf{DO^{-M}\Sigma}_{\mathcal{F}}\mathbf{O^{M}}\right)^{{\dagger}}=\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\left(\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|}\right)^{{\dagger}},$

where let
$\mathbf{P}=(\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|})^{{\dagger}}$, thus
$\mathbf{P}\mathbf{D}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}=\mathbf{I}_{|\mathcal{F}|\times|\mathcal{F}|}$
is a $|\mathcal{F}|\times|\mathcal{F}|$ identity matrix.

###### Proof.

Since $\mathrm{rank}(\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|})=|\mathcal{F}|$
and
$\mathrm{rank}((\mathbf{PD}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|})=|\mathcal{F}|$,
$\mathrm{rank}(\mathbf{P})=|\mathcal{F}|$, i.e. $\mathbf{P}$ spans
$\mathbb{C}^{|\mathcal{F}|}$. Hence,
$\mathbf{R}=\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\mathbf{P}$ spans
$\mathrm{BL}_{|\mathcal{F}|}(\mathbf{O}^{\mathbf{M}})$.

Let $\mathbf{J=RD}$, we have

$\displaystyle\mathbf{J}^{2}$
$\displaystyle=\mathbf{RDRD}=\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\left(\mathbf{PD}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\right)\mathbf{PD}$
$\displaystyle=\left(\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\mathbf{I}_{|\mathcal{F}|\times|\mathcal{F}|}\mathbf{P}\right)\mathbf{D=RD}=\mathbf{J},$

thus $\mathbf{J=PD}$ is a projection operator.

Since $\mathbf{RD}$ is a projection operator, $\mathbf{RDx}$ is an
approximation of $\mathbf{x}$ in the space of
$\mathrm{BL}_{|\mathcal{F}|}\left(\mathbf{O}^{\mathbf{M}}\right)$. We achieve
perfect recovery when $\mathbf{x}$ is in the space of
$\mathrm{BL}_{|\mathcal{F}|}\left(\mathbf{O}^{\mathbf{M}}\right)$. ∎

######  Remark 3.

When $|\mathcal{S}|<|\mathcal{F}|$, it holds that
$\mathrm{rank}(\mathbf{PD}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|})\leq\mathrm{rank}(\mathbf{P})\leq|\mathcal{S}|<|\mathcal{F}|$,
thus rendering $\mathbf{PDO}^{-\mathbf{M}}_{|\mathcal{F}|}$ incapable of being
an identity matrix, and perfect recovery of the original signal is
unattainable. In the case of $|\mathcal{S}|=|\mathcal{F}|$, for
$\mathbf{PDO}^{-\mathbf{M}}_{|\mathcal{F}|}$ to be the identity matrix,
$\mathbf{P}$ must be the inverse of
$\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|}$. For
$|\mathcal{S}|>|\mathcal{F}|$, $\mathbf{P}$ is the pseudo-inverse of
$\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|}$, and in cases where
$|\mathcal{S}|\geq|\mathcal{F}|$, perfect recovery of the original signal is
possible. For simplicity, we consider the scenario
$|\mathcal{S}|=|\mathcal{F}|$.

###  4.3 Qualified Sampling Operator

The implications of Theorem 3 underscore that achieving perfect recovery may
be unattainable, even for $\mathbf{M}$-bandlimited graph signals, with
arbitrary sampling operators. A qualified sampling operator must, at a
minimum, select $|\mathcal{F}|$ linearly independent rows from within
$\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$. We will introduce the definition
of a qualified sampling operator, ensuring its adherence to the conditions
delineated in Theorem 3. To achieve this, we formally define the adjacency
matrix in the graph linear canonical domain, thus

$\mathbf{A^{M}=O^{-M}\Lambda_{A}O^{M}},$ (26)

where $\mathbf{A^{M}}$ lacks physical significance as it represents a complex
matrix, specifically the adjacency matrix. However, when
$\mathbf{M}=(0,1;-1,0)$, $\mathbf{O}^{\mathbf{M}}$ degenerates into
$\mathbf{V}^{-1}$, and in this case, $\mathbf{A^{M}}$ serves as the shift
operator.

For any sampling and recovery operators satisfying Theorem 3, when
$\mathbf{x}$ is $\mathbf{M}$-bandlimited, we have

$\mathbf{x=RDx=O}^{-\mathbf{M}}_{|\mathcal{F}|}\mathbf{Px}_{\mathcal{S}}=\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\mathbf{\hat{x}}_{|\mathcal{F}|},$

where $\mathbf{\hat{x}}_{|\mathcal{F}|}$ represents the first $|\mathcal{F}|$
values of $\mathbf{\hat{x}}$. Therefore, we derive

$\mathbf{x}_{\mathcal{S}}=\mathbf{P}^{-1}\mathbf{\hat{x}}_{|\mathcal{F}|}=\mathbf{P}^{-1}\mathbf{\hat{x}}_{\mathcal{S}}.$

As observed, the sampled signal $\mathbf{x}_{\mathcal{S}}$ and the frequency
signal $\mathbf{\hat{x}}_{|\mathcal{F}|}$ in the dual domain form a graph
linear canonical operator. Therefore, we obtain a new adjacency matrix,
leading to the following theorem similar to [11].

######  Theorem 4.

For $\mathbf{x}\in\mathrm{BL}_{|\mathcal{F}|}(\mathbf{O^{M}})$, where
$|\mathcal{S}|=|\mathcal{F}|$, and if $\mathbf{D}$ is a qualified sampling
operator, then

$\mathbf{A}^{\mathbf{M}}_{\mathcal{S}}=\mathbf{P^{-1}\Lambda_{A_{|\mathcal{F}|}}P}\in\mathbb{C}^{|\mathcal{F}|\times|\mathcal{F}|},$
(27)

where $\mathbf{P}=(\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|})^{-1}$, and
$\mathbf{\Lambda_{A_{|\mathcal{F}|}}}$ is the diagonal matrix of the first
$|\mathcal{F}|$ eigenvalues of $\mathbf{\Lambda_{A}}$.

In the graph linear canonical domain, $\mathbf{A}^{\mathbf{M}}_{\mathcal{S}}$
is obtained through $\mathbf{A}^{\mathbf{M}}$ sampling. The original graph
signal $\mathbf{x}$ and the sampled graph signal $\mathbf{x}_{\mathcal{S}}$
have equivalent frequency content after undergoing the corresponding GLCT. The
information preserved by $\mathbf{A}^{\mathbf{M}}_{\mathcal{S}}$ can be
expressed as follows

$\displaystyle\mathbf{x}_{\mathcal{S}}-\mathbf{A}^{\mathbf{M}}_{\mathcal{S}}\mathbf{x}_{\mathcal{S}}=$
$\displaystyle\mathbf{P^{-1}\hat{x}}_{\mathcal{S}}-\mathbf{P^{-1}\Lambda_{A_{|\mathcal{F}|}}PP^{-1}\hat{x}}_{\mathcal{S}}$
$\displaystyle=$
$\displaystyle\mathbf{P^{-1}\left(\hat{x}_{\mathcal{S}}-\Lambda_{A_{|\mathcal{F}|}}\hat{x}_{\mathcal{S}}\right)}$
$\displaystyle=$
$\displaystyle\mathbf{DO^{-M}_{|\mathcal{F}|}\left(I-\Lambda_{A_{|\mathcal{F}|}}\right)O^{M}_{|\mathcal{F}|}x}$
$\displaystyle=$ $\displaystyle\mathbf{D\left(x-A^{M}x\right)},$

where,
$\mathbf{x}\in\mathrm{BL}_{|\mathcal{F}|}\left(\mathbf{O}^{\mathbf{M}}\right)$,
and $\mathbf{\hat{x}}_{|\mathcal{F}|}=\mathbf{\hat{x}}_{\mathcal{S}}$.

By reordering and rearranging the eigenvalues in the matrix of the GLCT,
Theorem 4 applies to all graph signals that are $\mathbf{M}$-bandlimited in
the GLCT domain.

###  4.4 Optimal Sampling Operator

When sampling graph signals, aside from determining the appropriate sample
size, devising a strategy for selecting sampling node locations is crucial as
the positions of the sampled nodes play a pivotal role in the performance of
the reconstruction algorithm. Given the multiple choices of $|\mathcal{F}|$
linearly independent rows in $\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$, our
objective is to opt for an optimal set that minimizes the impact of noise. One
strategy involves selecting sampling positions to minimize the normalized MSE
(NMSE) [12, 34]. Consider the noise $e$ introduced during the sampling process

$\mathbf{x}_{\mathcal{S}}=\mathbf{Dx}+e,$

where $\mathbf{D}$ is qualified. The recovered signal
$\mathbf{x}_{\mathcal{R}}$ is then given by

$\mathbf{x}_{\mathcal{R}}=\mathbf{Rx}_{\mathcal{S}}=\mathbf{RDx}+\mathbf{R}e=\mathbf{x}+\mathbf{R}e.$

Consequently, the bound of the recovery error, obtained through the norm and
Cauchy-Schwarz inequality, is

$||\mathbf{x}_{\mathcal{R}}-\mathbf{x}||_{2}=||\mathbf{R}e||_{2}=||\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\mathbf{P}e||_{2}\leq||\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}||_{2}||\mathbf{P}||_{2}||e||_{2}.$

Given that $||\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}||_{2}$ and $||e||_{2}$
are fixed, the objective is to minimize the spectral norm of $\mathbf{P}$,
which is the pseudo-inverse of
$\mathbf{D}\mathbf{O}^{-\mathrm{M}}_{|\mathcal{F}|}$. The norm considered here
can be either the 2-norm or the Frobenius norm. The optimal norm is selected
based on optimization algorithms and numerical methods.

In the following, we provide several alternative strategies for selecting
sampling sets in the GLCT domain based on the optimal design approaches
proposed by [11, 12, 13, 14, 38].

  1. 1.

Minimization of the Frobenius Norm of $\mathbf{P}$ (MinFro)
In this strategy, the objective is to directly minimize the Frobenius norm of
the matrix $\mathbf{P}$, aiming to minimize the NMSE between the original and
recovery signals. The method involves selecting columns from the matrix
$\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$ to minimize the Frobenius norm of
the pseudo-inverse matrix $\mathbf{P}$,

$\mathbf{D}^{\mathrm{opt}}=\arg\min_{\mathbf{D}}~{}||\mathbf{P}||_{F}=\arg\min_{\mathbf{D}}~{}||\mathbf{\mathbf{\Sigma_{\mathcal{F}}O^{M}D}}||_{F}.$
(28)

  2. 2.

Maximization of the Volume of the Parallelepiped Formed With the Columns of
$\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$ (MaxVol)
This strategy, inspired by [12], aims to select a set $\mathcal{S}$ of columns
from the matrix $\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$ to maximize the
volume of the parallelepiped constructed with the chosen columns. Essentially,
this approach involves completing the matrix
$\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$. The volume can be computed as the
determinant of the matrix
$\mathbf{O}^{\mathbf{M}}_{|\mathcal{F}|}\mathbf{D}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$.

$\mathbf{D}^{\mathrm{opt}}=\arg\max_{\mathbf{D}}~{}\prod^{|\mathcal{F}|}_{i=1}\lambda_{i}\left(\mathbf{O}^{\mathbf{M}}_{|\mathcal{F}|}\mathbf{D}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\right).$
(29)

  3. 3.

Minimization of the singular value of $\mathbf{P}$ (MinPinv)
Due to
$\lambda_{i}\left(\mathbf{B^{M}DB^{M}}\right)=\sigma^{2}_{i}\left(\mathbf{B^{M}D}\right)=\sigma^{2}_{i}\left(\mathbf{\Sigma_{\mathcal{F}}O^{M}D}\right)=\sigma^{2}_{i}\left(\mathbf{P}\right)$,
in the case of uncorrelated noise, this is equivalent to minimizing
$\sum^{\left|F\right|}_{i}\sigma^{2}_{i}$,

$\mathbf{D}^{\mathrm{opt}}=\arg\min_{\mathbf{D}}~{}\sum^{\left|\mathcal{F}\right|}_{i=1}\sigma^{2}_{i}(\mathbf{P}).$
(30)

  4. 4.

Maximization of the minimum singular value of
$\mathbf{D}\mathbf{O}^{-\mathrm{M}}_{|\mathcal{F}|}$ (MaxSigMin)
The fourth strategy is designed to leverage the 2-norm to minimize NMSE, which
is equivalent to minimizing the largest singular value of $\mathbf{P}$. In
other words, we aim to maximize the smallest singular value
$\sigma_{\mathrm{min}}$ of
$\mathbf{D}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$ for every qualified
$\mathbf{D}$,

$\mathbf{D}^{\mathrm{opt}}=\arg\max_{\mathbf{D}}~{}\sigma_{\mathrm{min}}(\mathbf{DO}^{-\mathbf{M}}_{|\mathcal{F}|}).$
(31)

  5. 5.

Maximization of the singular value of
$\mathbf{\mathbf{D}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}}$ (MaxSig)
The third strategy aims to select columns of the matrix
$\mathbf{O}^{-\mathrm{M}}_{|\mathcal{F}|}$ to maximize its singular values.
Although this strategy is not directly related to NMSE optimization, it is
easy to implement and demonstrates good performance in practice,

$\mathbf{D}^{\mathrm{opt}}=\arg\max_{\mathbf{D}}~{}\sum^{\left|\mathcal{F}\right|}_{i=1}\sigma^{2}_{i}(\mathbf{\mathbf{D}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}}).$
(32)

To address the sampling strategies outlined above, we introduce five greedy
algorithms aimed at resolving this selection challenge to ascertain the
optimal sampling operator. In tackling similar matrix approximation
optimization issues, greedy algorithms have demonstrated efficacy in closely
approximating the global optimum [39]. They provide an efficient and
relatively straightforward approach to resolving the sampling set selection
dilemma. By iteratively selecting locally optimal solutions, greedy algorithms
approximate the global optimum, thereby achieving a balance between
computational complexity and the quality of the solution. Since strategies 1)
and 5), as well as 3) and 4), are essentially equivalent, we present Algorithm
1, which summarizes the three sampling strategies 2), 4), and 5). Furthermore,
Table 1 encapsulates the A, D, E, T-optimal designs [38] and objective
functions associated with the aforementioned sampling strategies.

Algorithm 1 Optimal sampling operator based on MaxSig, MaxSigMin, or MaxVol

Input

$\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$ the first $|\mathcal{F}|$ columns
of $\mathbf{O}^{-\mathbf{M}}$

$|\mathcal{S}|$ the number of samples

$N$ the number of rows in $\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}$

Output

$\mathcal{S}$ sampling set

Function

$\begin{array}[]{l}\mathbf{for~{}}k=1:|\mathcal{S}|\\\ \quad
s=\arg\max_{j\in\\{1:N\\}-\mathcal{S}}\sum^{\left|\mathcal{F}\right|}_{i=1}\sigma^{2}_{i}\left(\left(\mathbf{O^{-M}_{|\mathcal{F}|}}\right)_{\mathcal{S}\cup\left\\{j\right\\}}\right)\\\
\text{or}\
s=\arg\max_{j\in\\{1:N\\}-\mathcal{S}}\sigma_{\mathrm{min}}\left(\left(\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\right)_{\mathcal{S}+\\{j\\}}\right)\\\
\text{or}\
s=\arg\max_{j\in\\{1:N\\}-\mathcal{S}}\prod^{|\mathcal{F}|}_{i=1}\lambda_{i}\left(\left(\left(\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\right)^{\ast}_{S\cup\left\\{j\right\\}}\left(\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\right)_{S\cup\left\\{j\right\\}}\right)\right)\\\
\quad\mathcal{S}=\mathcal{S}+\\{s\\}\\\ \mathbf{end}\\\
\mathbf{return}~{}\mathcal{S}\end{array}$

Table 1: Optimal design and objective function of sampling methods Optimal Designs | Objectives | Objective Functions | Methods
---|---|---|---
A-optimal | $\min\mathrm{tr}\left[\left(\mathbf{B^{M}DB^{M}}\right)^{-1}\right]$ | $\min\ \left|\left|\mathbf{\Sigma_{\mathcal{F}}O^{M}D}\right|\right|_{F}$ | MinFro
D-optimal | $\min\mathrm{det}\left[\left(\mathbf{B^{M}DB^{M}}\right)^{-1}\right]$ | $\max\prod^{|\mathcal{F}|}_{i=1}\lambda_{i}\left(\mathbf{O}^{\mathbf{M}}_{|\mathcal{F}|}\mathbf{D}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}\right)$ | MaxVol
E-optimal | $\min\left|\left|\left(\mathbf{DB^{M}}\right)^{{\dagger}}\right|\right|_{2}$ | $\min\sum^{\left|\mathcal{F}\right|}_{i=1}\sigma^{2}_{i}\left(\mathbf{\Sigma_{\mathcal{F}}O^{M}D}\right)$ | MinPinv
E-optimal | $\min\left|\left|\left(\mathbf{DB^{M}}\right)^{{\dagger}}\right|\right|_{2}$ | $\max\ \sigma_{\min}\left(\mathbf{DO^{-M}_{|\mathcal{F}|}}\right)$ | MaxSigMin
T-optimal | $\max\mathrm{tr}\left[\left(\mathbf{B^{M}DB^{M}}\right)\right]$ | $\max\sum^{\left|\mathcal{F}\right|}_{i=1}\sigma^{2}_{i}\left(\mathbf{\mathbf{D}\mathbf{O}^{-\mathbf{M}}_{|\mathcal{F}|}}\right)$ | MaxSig

We compared our proposed five sampling strategies with random sampling under
the GLCT sampling framework. As shown in Fig. 3, the NMSE is reported, defined
as the error between the original and recovered signals divided by the
original signal for each node. (a) is based on a cycle graph with N = 32
nodes, and the signal is a $\mathbf{M}$-bandlimited signal with parameters
$\alpha=0.8,\beta=32,\xi=0.5k+1,k=1,...,N$. The frequency bandwidth is
$|\mathcal{F}|=4$. (b) is based on a swiss roll graph with N = 256 nodes, and
the signal is a $\mathbf{M}$-bandlimited signal with same parameters. The
frequency bandwidth is $|\mathcal{F}|=32$. We compare the NMSE of six sampling
methods as the number of sampled nodes increases for two different graph
signals. We observed that, with a small number of sampling nodes, all methods
except random sampling could effectively recover the signal.

(a) Ring graph.

(b) Swiss roll graph.

Figure 3: NMSE vs. number of samples for different sampling strategies.

##  5 Applications and Numerical Experiments

In this section, we conduct simulation experiments and test GLCT sampling and
recovery in a semi-supervised classification application [34] and clustering
of bus test cases [40]. Semi-supervised classification seeks to categorize a
vast dataset with the assistance of a limited number of labeled instances. In
the context of GSP, the smoothness of labels across the graph can be observed
through low-pass characteristics in the graph spectral domain. Simultaneously,
in the case of the bus test case with 118 nodes, the generator dynamically
produces smooth graph signals, justifying the validity of the
$\mathbf{M}$-bandlimited assumption. We utilize four frameworks (Laplacian-
based GFT sampling, GFT sampling, GFRFT sampling, and GLCT sampling) for GSP,
comparing the NMSE of graph signal recovery.

###  5.1 Simulations

We conducted simulations using the graph signal processing toolbox (GSPBox)
[41] for two examples. The first example involves a David sensor network with
$N=500$ nodes. We set the bandwidth $|\mathcal{F}|=60$, and the
$\mathbf{M}$-bandlimited parameters to
$\alpha=0.7,\beta=N/2,\xi=0.5k+0.5,k=1,...,N$. The graph signal is the real
part of
$\mathbf{x}=\mathrm{real}(\mathbf{f}_{1}+0.5\mathbf{f}_{2}+\mathbf{f}_{3})$
obtained through $\mathbf{M}$-bandlimited operation, where $\mathbf{f}_{i}$
represents the $i$-th column of $\mathbf{O^{-M}}$. The second example involves
the vehicular traffic data in Rome. We set the bandwidth $|\mathcal{F}|=50$,
and the $\mathbf{M}$-bandlimited parameters to
$\alpha=0.9,\beta=N,\xi=0.5k+0.7,k=1,...,N$. The graph signal is the real
vehicular traffic data in Rome [42] obtained through $\mathbf{M}$-bandlimited
operation. The original graph signals are depicted in Fig. 4.

(a) Example 1: $500$-node David sensor network graph.

(b) Example 2: Roman vehicular traffic data graph.

Figure 4: Original graph signal.

For both examples, we performed sampling and recovery of the graph signals
using four sampling frameworks based on five optimal sampling strategies. All
signal recovery formulas utilize Eq. (24) from Theorem 2. The accuracy of the
recovery was measured using the NMSE. Fig. 5 illustrates the locations of
sampled nodes. The number of sampled nodes, $|\mathcal{S}|$, is 60 for (a),
(b), (c), and (d), and 50 for (e), (f), (g), and (h). Fig. 6 displays the
signals recovered through Laplacian-based GFT sampling in (a) and (e), GFT
sampling in (b) and (f), GFRFT sampling in (c) and (g), and GLCT sampling in
(d) and (h). The NMSE results for the David sensor network and Roman vehicular
traffic data, across four sampling frameworks employing five optimal
strategies, are detailed in Tables 2 and 3, respectively. Note that the
parameters selected for our GFRFT and GLCT are consistent with those specified
in the aforementioned two examples. Since both A-optimal and T-optimal designs
target optimizing the trace of an operator, they yield identical results.
Meanwhile, MinPinv and MaxSigMin, essentially embodying E-optimal design,
produce similar outcomes. It is evident that Laplacian-based GFT, the GFT, and
the GFRFT sampling fail to perfectly recover the graph signals, as these
signals, while $\mathbf{M}$-bandlimited, are not strictly bandlimited in the
Laplacian-based GFT, the GFT and the GFRFT domains. Only the GLCT sampling
achieves perfect signal recovery.

(a) David sensor network samples from Laplacian-based GFT sampling.

(b) David sensor network samples from GFT sampling.

(c) David sensor network samples from GFRFT sampling.

(d) David sensor network samples from GLCT sampling.

(e) Roman vehicular traffic data samples from Laplacian-based GFT sampling.

(f) Roman vehicular traffic data samples from GFT sampling.

(g) Roman vehicular traffic data samples from GFRFT sampling.

(h) Roman vehicular traffic data samples from GLCT sampling.

Figure 5: The location of the sampled signals.

(a) David sensor network recovery from Laplacian-based GFT sampling.

(b) David sensor network recovery from GFT sampling.

(c) David sensor network recovery from GFRFT sampling.

(d) David sensor network recovery from GLCT sampling.

(e) Roman vehicular traffic data recovery from Laplacian-based GFT sampling.

(f) Roman vehicular traffic data recovery from GFT sampling.

(g) Roman vehicular traffic data recovery from GFRFT sampling.

(h) Roman vehicular traffic data recovery from GLCT sampling.

Figure 6: Recovered graph signals. Table 2: NMSE vs. different sampling frames for David sensor network based on sampling strategies corresponding to optimal designs Sampling methods | Laplacian-based GFT | GFT | GFRFT | GLCT
---|---|---|---|---
MinFro | $6.3526\times 10^{11}$ | $4.6075\times 10^{13}$ | $4.1884$ | $4.5432\times 10^{-22}$
MaxVol | $2.6311\times 10^{15}$ | $6.0047\times 10^{6}$ | $3.2970$ | $3.7001\times 10^{-24}$
MinPinv | $0.9631$ | $1.6822$ | $1.5754$ | $5.4518\times 10^{-25}$
MaxSigMin | $0.9635$ | $1.7311$ | $1.7956$ | $1.8202\times 10^{-25}$
MaxSig | $6.3526\times 10^{11}$ | $4.6075\times 10^{13}$ | $4.1884$ | $4.5432\times 10^{-22}$
Table 3: NMSE vs. different sampling frames for Roman vehicular traffic data based on sampling strategies corresponding to optimal designs Sampling methods | Laplacian-based GFT | GFT | GFRFT | GLCT
---|---|---|---|---
MinFro | $1.2721\times 10^{14}$ | $5.1085\times 10^{11}$ | $1.3539\times 10^{3}$ | $8.7596\times 10^{-23}$
MaxVol | $2.9563\times 10^{11}$ | $4.1516\times 10^{13}$ | $15.4921$ | $4.0667\times 10^{-24}$
MinPinv | $1.0901$ | $1.0038$ | $1.0031$ | $1.3541\times 10^{-25}$
MaxSigMin | $1.0888$ | $1.0038$ | $1.0086$ | $1.7686\times 10^{-25}$
MaxSig | $1.2721\times 10^{14}$ | $5.1085\times 10^{11}$ | $1.3539\times 10^{3}$ | $8.7596\times 10^{-23}$

###  5.2 Classifying Online Blogs

In this section, we classify a real-world political blog dataset comprising
$N=1224$ nodes, categorizing them into conservative or liberal factions [43].
Nodes within the graph denote blogs, whereas edges symbolize hyperlinks
interlinking the blogs. Label signals are defined such that 0 represents
conservative and 1 denotes liberal, with these label signals spanning the
entire frequency band. By employing sampling to reduce the number of vertices,
we facilitate easier storage and transmission. Subsequently, through signal
recovery, we reconstruct the complete set of blogs and evaluate the
performance of our sampling method based on classification accuracy. In this
case, we employ the MaxSigMin sampling strategy and adjust the parameters of
GLCT to determine the optimal performance and, as noted in [34], we know that
the best results are achieved when the parameter $\alpha$ approaches 1, albeit
not precisely equal to 1. Therefore, while adjusting other parameters, we
approximately recover the labeled signal using the lowest $|\mathcal{F}|$
frequency components by solving the following optimization problem

$\mathbf{\hat{x}^{\mathrm{opt}}_{|\mathcal{F}|}}=\arg\min_{\mathbf{\hat{x}}_{|\mathcal{F}|}\in\mathbb{R}^{|\mathcal{F}|}}~{}\left\|\mathrm{thres}\left(\mathbf{DO^{-M}_{|\mathcal{F}|}\hat{x}_{|\mathcal{F}|}}\right)-\mathbf{x}_{\mathcal{S}}\right\|^{2}_{2},$
(33)

where, the threshold function $\mathrm{thres}(\cdot)$ assigns a value of 1 to
all values greater than 0.5 and 0 otherwise. Subsequently, we set another
threshold of 0.5 [11] to assign labels to the recovered signal

$\mathbf{x}^{\mathrm{opt}}_{\mathcal{R}}=\mathrm{thres}\left(\mathbf{O^{-M}_{|\mathcal{F}|}}\mathbf{\hat{x}^{\mathrm{opt}}_{|\mathcal{F}|}}\right).$
(34)

Since the GFT and the GFRFT are special cases of the GLCT, by adjusting the
parameters of the GLCT, all methods can be encompassed. In our previous work
[28], the parameter $\xi$ only influenced the real and imaginary components of
the spectrum, while the absolute value of the overall signal remained
unchanged. Thus, we set $\xi=0.5k+1,k=1,...,N$ in this study. We varied the
parameter $\alpha$ from $0.99$ to $1$. For the parameter $\beta=m\times
2^{n},n\in\mathbb{R}$, as the accuracy remains constant when $m$ is fixed, we
considered variations in $m$, setting $m=1,3,5,N$. Fig. 7 illustrates the
spectral representations with different bandwidths $|\mathcal{F}|$ of the GLCT
with parameters $\xi=0.5k+1,\beta=2,\alpha=0.995$: (a) shows
$|\mathcal{F}|=2$, (b) $|\mathcal{F}|=612$, and (c) $|\mathcal{F}|=1224$. Fig.
8(a) illustrates the change in classification accuracy for nodes with a
frequency bandwidth $|\mathcal{F}|=2$ and a sampling size
$|\mathcal{S}|=\mathrm{M}=5$ as we vary the parameter $\alpha$. The accuracy
reaches its maximum value of 94.6078% when $\alpha=0.995$ and $\beta=2^{n}$.
In Fig. 8(b), we present a comparison of classification accuracies with
increasing sampling size $\mathrm{M}$ for GFT, GFRFT [34], and GLCT sampling.
For the GFRFT and GLCT methods, the parameters selected correspond to the
optimal values identified in Fig. 8(a). It is evident that GLCT outperforms
GFRFT and GFT across varying sampling sizes $\mathrm{M}$.

(a) The spectrum of the GLCT with parameters
$\xi=0.5k+1,\beta=2,\alpha=0.995$, and $|\mathcal{F}|=2$.

(b) The spectrum of the GLCT with parameters
$\xi=0.5k+1,\beta=2,\alpha=0.995$, and $|\mathcal{F}|=612$.

(c) The spectrum of the GLCT with parameters
$\xi=0.5k+1,\beta=2,\alpha=0.995$, and $|\mathcal{F}|=1224$.

Figure 7: Bandlimited spectrum of online political blogs.

(a) Classification accuracy as a function of fractional order $\alpha$. The
top-right graph illustrates the detailed classification accuracy of GFRFT and
GLCT with $\beta=2^{n}$. The bottom-right graph represents the detailed
classification accuracy with $\beta=3\cdot 2^{n},5\cdot 2^{n},N\cdot 2^{n}$.

(b) Classification accuracy as a function of sample size $\mathrm{M}$ for GFT,
GFRFT, and GLCT sampling. the blue line denotes GFT, the red line denotes the
GFRFT of $\alpha=0.996$ [34], and the yellow line denotes the GLCT of
$\xi=0.5k+1,\beta=2^{n},\alpha=0.995$.

Figure 8: Classification of online political blogs.

###  5.3 Clustering of Bus Test Case

The second part of our study focuses on sampling within an electrical grid,
specifically utilizing the IEEE 118 bus test case, which represents a section
of the Midwestern US power system as of December 1962. This network is
composed of 118 vertices (i.e., buses) interconnected by edges (i.e.,
transmission lines) as depicted in Fig. 9(a) [44]. The coloring of nodes is
indicative of the components of eigenvectors, which are aligned with the
eigenvalues of the GLCT operator matrix in an ascending order. These
eigenvalues are derived under the GLCT parameters configured as
$\xi=0,\beta=2,\alpha=1$. The graph signal is represented by
$\mathbf{x}=\mathbf{O}^{-\mathbf{M}}_{N}$, where the notation signifies the
$N$th column of the matrix. As suggested in [40], the dynamics of generators
yield smooth graph signals, making the assumption of
$\mathbf{M}$-bandlimitedness plausible, albeit in an approximate sense, as
demonstrated by the spectral plot of a bandlimited graph signal in Fig. 9(b).
Network parameters are available at
http://www.cise.ufl.edu/research/sparse/matrices/, with their layout derived
through a graph drawing algorithm described in [45].

(a) The graph signal characterized by vertex signals
$\mathbf{x}=\mathbf{O}^{-\mathbf{M}}_{N}$, with the parameters set as
$\xi=0,\beta=2,\alpha=1$.

(b) The spectral representation of the graph signal $\mathbf{x}$ after
$\mathbf{M}$-bandlimited.

Figure 9: Graph signal representation of IEEE 118 bus test cases.

In our demonstration, we generate random low-pass signals with a graph
spectral bandwidth $|\mathcal{F}|=8$ and collect a sample size
$|\mathcal{S}|=8$. We then count the number of triangles in the graph to
select optimal sampling parameters using different GLCT settings, employing
the MaxSigMin sampling strategy, and perform $K$-means clustering [46],
choosing $K=8$, with results illustrated in Fig. 10. And the signal recovery
formula is Eq. (24). In the context of the IEEE 118 bus test case, the pivotal
utility of cluster analysis lies in discerning the functional zones and
critical nodes within the power system. This discernment is instrumental in
guiding the optimization of grid scheduling and bolstering the system’s
stability and efficiency. Leveraging the characteristics of these categorized
groups enables the refinement of power generation scheduling strategies, aimed
at minimizing the overall generation costs while ensuring the system’s stable
operation. Moreover, cluster analysis serves a crucial role in identifying
nodes that may significantly impact the system’s stability. This
identification forms a foundational basis for the development of fault
prevention measures and the formulation of emergency response plans [47, 48,
46].

(a) $K$-means clustering of IEEE 118 bus test cases.

Figure 10: $K$-means clustering of IEEE 118 bus test cases.

Finally, we evaluate the clustering results using

$\text{Silhouette Score}=\frac{b(i)-a(i)}{\max\\{a(i),b(i)\\}},$

where, $a(i)$ represents the average distance of sample $i$ to other samples
within the same cluster (intra-cluster compactness), and $b(i)$ represents the
average distance of sample $i$ to the nearest samples in other clusters
(inter-cluster separation). A result closer to 1 indicates that samples are
correctly assigned to clusters, with good separation between clusters.

(a) Silhouette score as a function of fractional order $\alpha$, the blue line
denotes GFRFT, the red line denotes the GLCT of $\xi=0,\beta=2$, and the
yellow line denotes the GLCT of $\xi=0.5k+1,\beta=2$.

(b) Silhouette score as a function of sample size $\mathrm{M}$, the blue line
denotes GFRFT of $\alpha=1$ (GFT), the red line denotes the GFRFT of
$\alpha=0.9$, and the yellow line denotes the GLCT of
$\xi=0,\beta=2,\alpha=1$.

Figure 11: Clustering of IEEE 118 Bus Test Case.

We conducted experiments by employing the method of controlling variables,
adjusting the GLCT parameters and the number of sampling notes. Fig. 11(a)
displays the silhouette score of fixed sampling $\mathrm{M}=|\mathcal{F}|=8$
nodes with varying parameter $\alpha$. The silhouette score reaches its
maximum value of $0.8753$ when $\xi=0,\beta=2$ and $\alpha=1$. In (b), we
compare the silhouette score of the GFT, GFRFT and GLCT with optimal
parameters as the sampling size $\mathrm{M}$ varies. Overall, the GLCT
exhibits superior performance.

##  6 Conclusion

In this paper, we propose the uncertainty principle and sampling theory in the
GLCT domain. We demonstrate that the uncertainty principle in the GLCT domain
has a broader scope compared to the GFT domain and shares interesting
connections with sampling theory. We show that $\mathbf{M}$-bandlimited graph
signals in the GLCT domain can achieve perfect recovery through the GLCT with
$\mathbf{M}$ matrix coefficients. We employ experimentally designed optimal
sampling strategies to ensure perfect recovery while maximizing robustness
against noise. Various proposed sampling strategies are compared.
Subsequently, we conduct simulation experiments and test GLCT sampling and
recovery in a semi-supervised classification application and clustering of bus
test cases. When comparing its performance with GFT and GFRFT sampling, we
observe that GLCT sampling achieves superior classification accuracy under the
optimal $\mathbf{M}$ matrix.

##  Appendix A Proof of Lemma 1

The proof for this statement is available in [12]. By utilizing Eq. (14) and
Eq. (15), we derive

$\mathbf{B^{M}DB^{M}x=B^{M}Dx=B^{M}x=x}.$

Consequently, $\lambda_{\max}\left(\mathbf{B^{M}DB^{M}}\right)=1$. On the
other hand, if $\mathbf{B^{M}DB^{M}x=x}$ is known, then
$\mathbf{B^{M}DB^{M}x=B^{M}x}$, as $(\mathbf{B^{M}})^{2}=\mathbf{B^{M}}$. This
further implies $\mathbf{B^{M}x=x}$. Regarding vertex localization, utilizing
the Rayleigh-Ritz theorem, we obtain

$\max_{\mathbf{x}}\mathbf{\frac{x^{\ast}Dx}{x^{\ast}x}}=\max_{\mathbf{x}}\mathbf{\frac{x^{\ast}B^{-M}DB^{M}x}{x^{\ast}x}}=1.$

Therefore, $\mathbf{Dx=x}$, indicating perfect localization in both the vertex
and spectral domains.

##  Appendix B Proof of Lemma 2

The proof of this bound requires the use of curves [22]

$\gamma\left(\zeta\right):=\left(\left(\zeta\lambda_{\max}\right)^{\frac{1}{2}}+\left(\left(1-\zeta\right)\left(1-\lambda_{\max}\right)\right)^{\frac{1}{2}}\right)^{2},\
\zeta\in\left[\lambda_{\max},1\right].$ (35)

For a normalized vector signal $\mathbf{x}$, we consider two normalized
vectors

$\mathbf{y}=\frac{\mathbf{Dx}}{\left|\left|\mathbf{Dx}\right|\right|_{2}},\
\text{and}\ \
\mathbf{z}=\frac{\mathbf{B}^{\mathbf{M}}\mathbf{x}}{\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}.$

The usual definition of inner product
$\left<\mathbf{y,z}\right>=\mathbf{y^{\ast}z}$, angular distance is the
measure of vectors on the unit sphere, we can define the angle distance
between two vectors as

$\angle\left(\mathbf{y,z}\right)=\arccos\mathbb{R}\left<\mathbf{y,z}\right>,$
(36)

where, $\angle$ denotes the angular distance,
$\mathbb{R}\left<\mathbf{y,z}\right>$ denotes the real part. In particular,
the sum of the angular distances between vectors $\mathbf{y}$ and
$\mathbf{x}$, and $\mathbf{z}$ and $\mathbf{x}$ is always larger than the
angular distance between $\mathbf{y}$ and $\mathbf{z}$, i.e.

$\arccos\mathbb{R}\left<\mathbf{y,x}\right>+\arccos\mathbb{R}\left<\mathbf{z,x}\right>\geq\arccos\mathbb{R}\left<\mathbf{y,z}\right>.$
(37)

For $\mathbb{R}\left<\mathbf{y,z}\right>$ , the upper bound is given by the
Cauchy-Schwarz inequality

$\displaystyle\mathbb{R}\left<\mathbf{y,z}\right>\leq$
$\displaystyle\left|\mathbb{R}\left<\mathbf{y,z}\right>\right|\leq\left|\left<\mathbf{y,z}\right>\right|$
$\displaystyle=$
$\displaystyle\frac{\left|\left<\mathbf{Dx},\mathbf{B}^{\mathbf{M}}\mathbf{x}\right>\right|}{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}=\frac{\left|\left<\mathbf{D}^{2}\mathbf{x},\mathbf{B}^{\mathbf{M}}\mathbf{x}\right>\right|}{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}$
$\displaystyle=$
$\displaystyle\frac{\left|\left<\mathbf{Dx,DB}^{\mathbf{M}}\mathbf{x}\right>\right|}{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}\leq\frac{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{DB}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}$
$\displaystyle=$
$\displaystyle\frac{\left|\left|\mathbf{Dx}\right|\right|_{2}\sqrt{\left<\mathbf{D}^{2}\mathbf{B}^{\mathbf{M}}\mathbf{x,B^{M}x}\right>}}{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}=\frac{\left|\left|\mathbf{Dx}\right|\right|_{2}\sqrt{\left<\mathbf{DB}^{\mathbf{M}}\mathbf{x},(\mathbf{B^{M}})^{2}\mathbf{x}\right>}}{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}$
$\displaystyle=$
$\displaystyle\frac{\left|\left|\mathbf{Dx}\right|\right|_{2}\sqrt{\left<\mathbf{B^{M}DB^{M}x,B^{M}x}\right>}}{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}\leq\frac{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}\left|\left|\mathbf{B^{M}DB^{M}x}\right|\right|_{2}}{\left|\left|\mathbf{Dx}\right|\right|_{2}\left|\left|\mathbf{B}^{\mathbf{M}}\mathbf{x}\right|\right|_{2}}$
$\displaystyle\leq$
$\displaystyle\sqrt{\lambda_{\max}\left(\mathbf{B^{M}DB^{M}}\right)}.$

Since we assume that $\lambda_{\max}\leq\zeta^{2}\eta^{2}$ and the right-hand
expression in (37) is less than 1, we can get

$\mathbb{R}\left<\mathbf{y,z}\right>\leq\sqrt{\lambda_{\max}\left(\mathbf{B^{M}DB^{M}}\right)}\leq
1,\ \text{and}\ \mathbb{R}\left<\mathbf{y,x}\right>=\zeta,\
\mathbb{R}\left<\mathbf{z,x}\right>=\eta^{\mathbf{M}}.$

So that formula (19) is proved. Next, we prove formula (20). Because the
function $\arccos(x)$ decreases monotonically at $(0,1]$, and according to
formula (19), we get

$\arccos\eta^{\mathbf{M}}\geq\arccos\sqrt{\lambda_{\max}}-\arccos\zeta,$

substituting both sides into the cosine function respectively

$\eta^{\mathbf{M}}\leq\cos\left(\arccos\sqrt{\lambda_{\max}}-\arccos\zeta\right).$

Applying the trigonometric identity $\cos(a-b)=\cos a\cos b+\sin a\sin b$, we
finally get the inequality

$\eta^{\mathbf{M}}\leq\sqrt{\lambda_{\max}}\zeta+\sqrt{1-\lambda_{\max}}\sqrt{1-\zeta^{2}}.$

## Declaration of competing interest

The authors declare that they have no known competing financial interests or
personal relationships that could have appeared to influence the work reported
in this paper.

## Acknowledgments

This work were supported by the National Natural Science Foundation of China
[No. 62171041], the BIT Research and Innovation Promoting Project
[No.2023YCXY053], and Natural Science Foundation of Beijing Municipality [No.
4242011]. The authors would like to thank Linyu Peng, Associate Professor in
the Department of Mechanical Engineering, Faculty of Science and Technology,
Keio University, for his careful guidance.

## References

  * [1] D. I. Shuman, S. K. Narang, P. Frossard, A. Ortega, P. Vandergheynst, The emerging field of signal processing on graphs: Extending high-dimensional data analysis to networks and other irregular domains, IEEE Signal Process. Mag. 30 (3) (2013) 83-98.
  * [2] A. Sandryhaila, J. M. F. Moura, Big data processing with signal processing on graphs: Representation and processing of massive data sets with irregular structure, IEEE Signal Process. Mag. 31 (5) (2014) 80-90.
  * [3] D. K. Hammond, P. Vandergheynst, R. Gribonval, Wavelets on graphs via spectral graph theory, Appl. Comput. Harmon. Anal. 30 (2) (2011) 129–150.
  * [4] A. Sandryhaila, J. M. F. Moura, Discrete signal processing on graphs, IEEE Trans. Signal Process. 61 (7) (2013) 1644-1656.
  * [5] A. Ortega, P. Frossard, J. Kovačević, J. M. F. Moura, P. Vandergheynst, Graph signal processing: Overview, challenges, and applications, Proc. IEEE, 106 (5) (2018) 808-828.
  * [6] G. Leus, A. G. Marques, J. M. F. Moura, A. Ortega, D. I. Shuman, Graph signal processing: History, development, impact, and outlook, IEEE Signal Process. Mag. 40 (4) (2023) 49-60.
  * [7] A. Sandryhaila, J. M. F. Moura, Discrete signal processing on graphs: Frequency analysis, IEEE Trans. Signal Process. 62 (12) (2014) 3042-3054.
  * [8] D. I. Shuman, B. Ricaud, P. Vandergheynst, Vertex-frequency analysis on graphs, Appl. Comput. Harmon. Anal. 40 (2) (2016) 260-291.
  * [9] M. W. Morency, G. Leus, Graphon filters: Graph signal processing in the limit, IEEE Trans. Signal Process. 69 (2021) 1740-1754.
  * [10] Z. Yang, X. Zheng, Z. Yu, X. Li, Graph filter design by ring-decomposition for 2-connected graphs, Signal Process. 201 (2022) 108725.
  * [11] S. Chen, R. Varma, A. Sandryhaila, J. Kovačević, Discrete signal processing on graphs: Sampling theory, IEEE Trans. Signal Process. 63 (24) (2015) 6510-6523.
  * [12] M. Tsitsvero, S. Barbarossa, P. D. Lorenzo, Signals on graphs: Uncertainty principle and sampling, IEEE Trans. Signal Process. 64 (18) (2016) 4845-4860.
  * [13] A. Anis, A. Gadde, A. Ortega, Efficient sampling set selection for bandlimited graph signals using graph spectral proxies, IEEE Trans. Signal Process. 64 (14) (2016) 3775-3789.
  * [14] A. Sakiyama, Y. Tanaka, T. Tanaka, A. Ortega, Eigendecomposition-free sampling set selection for graph signals, IEEE Trans. Signal Process. 67 (10) (2019) 2679-2692.
  * [15] D. Y. Wei, Z. Yan, Generalized sampling of graph signals with the prior information based on graph fractional Fourier transform, Signal Process. 214 (2024) 109263.
  * [16] J. Shi, J. M. F. Moura, Graph signal processing: Dualizing GSP sampling in the vertex and spectral domains, IEEE Trans. Signal Process. 70 (2022) 2883-2898.
  * [17] I. Zach, T. G. Dvorkind, R. Talmon, Graph signal interpolation and extrapolation over manifold of Gaussian mixture, Signal Process. 216 (2024) 109308.
  * [18] S. Chen, A. Sandryhaila, J. M. F. Moura, J. Kovačević, Signal recovery on graphs: Variation minimization, IEEE Trans. Signal Process. 63 (17) (2015) 4609-4624.
  * [19] R. Torkamani, A. Amini, H. Zayyani, M. Korki, Graph signal recovery using variational Bayes in Fourier pairs with Cramér–Rao bounds, Signal Process. 219 (2024) 109394.
  * [20] F. Ji, Pratibha, W. P. Tay, Unlimited dynamic range signal recovery for folded graph signals, Signal Process. 198 (2022) 108574.
  * [21] A. Agaskar, Y. M. Lu, A spectral graph uncertainty principle, IEEE Trans. Inf. Theory, 59 (7) (2013) 4338–4356.
  * [22] W. Erb, Shapes of uncertainty in spectral graph theory, IEEE Trans. Inf. Theory, 67 (2) (2021) 1291-1307.
  * [23] L. L. Magoarou, R. Gribonval, N. Tremblay, Approximate fast graph Fourier transforms via multilayer sparse approximations, IEEE Trans. Signal inf. Process. Netw. 4 (2) (2018) 407-420.
  * [24] Y. Q. Wang, B. Z. Li, Q. Y. Cheng, The fractional Fourier transform on graphs, in Proc. Asia-Pacific Signal Inf. Process. Assoc. Annu. Summit Conf. (APSIPA ASC), 2017, pp. 105–110.
  * [25] Y. Zhang, B. Z. Li, The fractional Fourier transform on graphs: Modulation and convolution, in Proc. 8th IEEE Int. Conf. Signal lmage Process. (ICSIP), 2023, pp. 737-741.
  * [26] F. J. Yan, B. Z. Li, Windowed fractional Fourier transform on graphs: Properties and fast algorithm, Digit. Signal Process. 118 (2021) 103210.
  * [27] F. J. Yan, B. Z. Li, Spectral graph fractional Fourier transform for directed graphs and its application, Signal Process. 210 (2023) 109099.
  * [28] Y. Zhang, B. Z. Li, Discrete linear canonical transform on graphs, Digit. Signal Process. 135 (2023) 103934.
  * [29] H. M. Ozaktas, Z. Zalevsky, M. A. Kutay, The fractional Fourier transform with applications in optics and signal processing, Wiley Sons, Chichester, New York, 2001.
  * [30] G. B. Folland, A. Sitaram, The uncertainty principle: A mathematical survey, J. Fourier Anal. Appl. 3.3, pp. 207–238, 1997.
  * [31] F. Ji, W. P. Tay, A Hilbert space theory of generalized graph signal processing, IEEE Trans. Signal Process. 67 (24) (2019) 6188-6203.
  * [32] D. Slepian, H. O. Pollak, Prolate spheroidal wave functions, Fourier analysis and uncertainty. I, Bell Syst. Tech. J. 40 (1) (1961) 43-63.
  * [33] H. J. Landau, H. O. Pollak, Prolate spheroidal wave functions, Fourier analysis and uncertainty. II, Bell Syst. Tech. J. 40 (1) (1961) 65-84.
  * [34] Y. Q. Wang, B. Z. Li, The fractional Fourier transform on graphs: Sampling and recovery, in Proc. 14th IEEE Int. Conf. Signal Process. (ICSP), 2018, pp. 1103-1108.
  * [35] S. C. Pei, Y. C. Lai, Discrete linear canonical transforms based on dilated Hermite functions, J. Opt. Soc. Am. A, 28 (8) (2011) 1695-1708.
  * [36] X. G. Xia, Discrete chirp-Fourier transform and its application to chirp rate estimation, IEEE Trans. Signal Process. 48 (11) (2000) 3122-3133.
  * [37] M. Tsitsvero, S. Barbarossa, P. D. Lorenzo, Uncertainty principle and sampling of signals defined on graphs, in Proc. 49th Asilomar Conf. Signals, Syst. Comput. 2015, pp. 1813-1818.
  * [38] B. J. Winer, D. R. Brown, K. M. Michels, Statistical principles in experimental design, New York: McGraw-Hill, 1971.
  * [39] H. Avron, C. Boutsidis, Faster subset selection for matrices and applications, SIAM J. Matrix Anal. Appl. 34 (4) (2013) 1464–1499.
  * [40] F. Pasqualetti, S. Zampieri, F. Bullo, Controllability metrics, limitations and algorithms for complex networks, IEEE Trans. Control Netw. Syst. 1 (1) (2014) 40–52.
  * [41] N. Perraudin, J. Paratte, D. Shuman, V. Kalofolias, P. Vandergheynst, D. K. Hammond, GSPBOX: a toolbox for signal processing on graphs, 2016, arXiv preprint, arXiv: 1408.5781.
  * [42] P. D. Lorenzo, S. Barbarossa, P. Banelli, Sampling and recovery of graph signals, Coop. Graph Signal Process. P. M. Djuric and C. Richard Eds. Elsevier, pp: 261-282, 2018.
  * [43] L. A. Adamic, N. Glance, The political blogosphere and the 2004 U.S. election: divided they blog, Proc. LinkKDD, pp. 35-43, 2005.
  * [44] K. Sun, D. Z. Zheng, Q. Lu, A simulation study of OBDD-based proper splitting strategies for power systems under consideration of transient stability, IEEE Trans. Power Syst. 20 (1) (2005) 389-399.
  * [45] Y. Hu, Efficient, high-quality force-directed graph drawing, Math. J. 10 (1) (2005) 37–71.
  * [46] U. V. Luxburg, A tutorial on spectral clustering, Statist. Comput. 17 (4) (2007) 395–416.
  * [47] D. S. Yang, Y. H. Sun, B. W. Zhou, X. T. Gao, H. G. Zhang, Critical nodes identification of complex power systems based on electric cactus structure, IEEE Syst. J. 14 (3) (2020) 4477-4488.
  * [48] P. Chanfreut, J. M. Maestre, E. F. Camacho, A survey on clustering methods for distributed and networked control systems, Annu. Rev. Control, 52 (2021) 75-90.