# STAR-RIS-Enabled Simultaneous Indoor and Outdoor 3D Localization:
Theoretical Analysis and Algorithmic Design ††thanks: J. He, A. Fakhreddine,
and G. C. Alexandropoulos are with the Technology Innovation Institute, 9639
Masdar City, Abu Dhabi, United Arab Emirates. G. C. Alexandropoulos is also
with the Department of Informatics and Telecommunications, National and
Kapodistrian University of Athens, 15784 Athens, Greece. (e-mails:
{jiguang.he<EMAIL_ADDRESS>alexandg@di.uoa.gr).

Jiguang He, Aymen Fakhreddine, and George C. Alexandropoulos

###### Abstract

Recent research and development interests deal with metasurfaces for wireless
systems beyond their consideration as intelligent tunable reflectors. Among
the latest proposals is the simultaneously transmitting (a.k.a. refracting)
and reflecting reconfigurable intelligent surface (STAR-RIS) which intends to
enable bidirectional indoor-to-outdoor, and vice versa, communications thanks
to its additional refraction capability. This double functionality provides
increased flexibility in concurrently satisfying the quality-of-service
requirements of users located at both sides of the metasurfaces, for example,
the achievable data rate and localization accuracy. In this paper, we focus on
STAR-RIS-empowered simultaneous indoor and outdoor three-dimensional (3D)
localization, and study the fundamental performance limits via Fisher
information analyses and Cramér Rao lower bounds (CRLBs). We also devise an
efficient localization algorithm based on an off-grid compressive sensing (CS)
technique relying on atomic norm minimization (ANM). The impact of the
training overhead, the power splitting at the STAR-RIS, the power allocation
between the users, the STAR-RIS size, the imperfections of the STAR-RIS-to-BS
channel, as well as the role of the multi-path components on the positioning
performance are assessed via extensive computer simulations. It is
theoretically showcased that high-accuracy, up to centimeter level, 3D
localization can be simultaneously achieved for indoor and outdoor users,
which is also accomplished via the proposed ANM-based estimation algorithm.

##  I Introduction

Besides contributing to communications for improved energy efficiency (EE) and
spectrum efficiency (SE) [1, 2, 3, 4, 5], reconfigurable intelligent surfaces
(RISs) also play a critical role in radio localization as well as environment
mapping, termed as simultaneous localization and mapping (SLAM), in current
and upcoming future cellular networks [6, 7, 8, 9, 10, 11]. In these radio
localization literature, the RIS behaves in various manners, e.g., a
programmable reflector in [6, 7, 8], a cost-efficient receiver in [9], and a
simultaneous reflector and refractor in [10]. In principle, localization
performance can be significantly boosted by deploying one or multiple RISs
thanks to the subsequent reasons: i) The number of reference nodes, including
base stations (BSs), can be further increased with the introduction of RISs;
Namely, the RIS can be considered as an additional anchor upon its deployment;
Its exact location can be shared with the surrounding BSs. ii) RIS creates a
virtual line-of-sight (LoS) link in the millimeter wave (e.g., 5G frequency
range 2 (FR2)) network when direct LoS link is temporally unavailable; This
happens frequently for millimeter wave communications, known as blockage; iii)
Provided that large-sized RISs are exploited, one can obtain high resolution
on angular parameters, e.g., angles of departure (AoDs) or angles of arrival
(AoAs), associated with RISs for the purpose of user localization; iv)
Tremendous RIS beamforming gain leads to enhanced received signal strength,
which in turn boosts the localization performance.

Among different types of RIS in [6, 7, 8, 9, 12, 10, 13], simultaneously
transmitting (a.k.a. refracting) and reflecting RIS (STAR-RIS) stands out as
it provides full-dimensional coverage (i.e., $360^{\circ}$ coverage). The
application of STAR-RIS for multiple-input multiple-output (MIMO)
communications can be referred to [14] for a general overview, [15] for
channel estimation, and [16] for non-orthogonal multiple access (NOMA)
transmissions. The STAR-RIS inherently offers two operation functionalities,
i.e., reflection and refraction, controlled by two separate series of phase
shifters. Such an extraordinary property can also be leveraged for radio
localization. For instance, an outdoor BS can simultaneously localize indoor
and outdoor mobile stations (MSs) with the aid of the STAR-RIS [10]. In this
example, the STAR-RIS serves as a tunable reflector for the outdoor MS and
meanwhile a tunable refractor for the indoor MS. With the introduced
flexibility on power/energy splitting and duplex mode between the two
functionalities, the quality-of-service (QoS) requirements in terms of
localization accuracy can be met concurrently for both indoor and outdoor MSs.

To the best of the authors’ knowledge, this is the first paper introducing
STAR-RIS for simultaneous indoor and outdoor three-dimensional (3D)
localization and analyzing the structure’s theoretical performance limits
[10]. However, practical localization algorithms are left undeveloped and
several practical issues are left uninvestigated. Thus, in this paper, we
continue to focus on the STAR-RIS-enabled simultaneous indoor and outdoor 3D
localization system, which comprises one indoor MS and one outdoor MS. The
localization of the two users is performed at the BS by considering the
received sounding reference signals (SRSs) transmitted over the uplink from
the two MSs simultaneously, i.e., in a NOMA manner. We summarize the
fundamental performance limits captured by Fisher information analyses and
Cramér Rao lower bounds (CRLBs), develop effective localization algorithms
based on co-channel interference mitigation and off-grid compressive sensing
(CS) technique, named atomic norm minimization (ANM), and examine the impact
brought by the practical issues, i.e., training overhead, power splitting at
STAR-RIS and power allocation between the two users, sup-optimal/optimal STAR-
RIS design, imperfectness of STAR-RIS-to-BS channel, and multi-path components
(MPCs), on the 3D localization performance of the two MSs.

The rest of the paper is organized as follows. Section II introduces the
system model, including channel and signal models. Section III summarizes the
CRLB analyses on the positioning errors and optimizes the STAR-RIS design. In
Section IV, we provide the practical localization algorithm based on ANM,
followed by numerical study and evaluation of different practical factors in
Section V. Finally, we provide the concluding remarks and point out several
future research directions in Section VI.

Notations: A bold lowercase letter ${\mathbf{a}}$ denotes a vector, and a bold
capital letter ${\mathbf{A}}$ denotes a matrix. $(\cdot)^{\mathsf{T}}$,
$(\cdot)^{*}$, and $(\cdot)^{\mathsf{H}}$ denote the matrix or vector
transpose, conjugate, and Hermitian transpose, respectively. $(\cdot)^{-1}$
denotes inverse of a matrix, $\mathrm{tr}(\cdot)$ denotes the trace operator,
$\mathrm{diag}({\mathbf{a}})$ denotes a square diagonal matrix with the
entries of ${\mathbf{a}}$ on its diagonal, ${\mathbf{A}}\otimes{\mathbf{B}}$
and ${\mathbf{A}}\diamond{\mathbf{B}}$ denote the Kronecker and Khatri-Rao
products of ${\mathbf{A}}$ and ${\mathbf{B}}$, respectively,
$\mathbb{E}[\cdot]$ and $\mathrm{var}(\cdot)$ are the expectation and variance
operators, $\mathbf{1}$ is the all-one vector, $\mathbf{0}$ denotes the all-
zero vector or matrix, ${\mathbf{I}}_{M}$ denotes the $M\times M$ identity
matrix, $j=\sqrt{-1}$, $\|\cdot\|_{\mathrm{F}}$ denotes the Frobenius norm of
a matrix, and $\|\cdot\|_{2}$ denotes the Euclidean norm of a vector.
$[{\mathbf{a}}]_{i}$, $[{\mathbf{A}}]_{ij}$, and $[{\mathbf{A}}]_{i:j,i:j}$
denote the $i$th element of ${\mathbf{a}}$, the $(i,j)$th element of
${\mathbf{A}}$, and the submatrix of ${\mathbf{A}}$ formed by rows
$i,i+1,\ldots,j$ and columns $i,i+1,\ldots,j$. Finally, $|\cdot|$ returns the
absolute value of a complex number.

##  II System Model

Figure 1: Simultaneous indoor and outdoor 3D localization empowered by the
deployment of the STAR-RIS, where an outdoor BS localizes one outdoor MS and
one indoor MS concurrently based on the received pilot signals over the
uplink. The two MSs are connected to the outdoor BS via the STAR-RIS with all
the links marked with solid arrows. Alternatively, an indoor BS can also
localize one indoor MS and one outdoor MS with the help of one STAR-RIS, where
all links are marked with dotted arrows.

The STAR-RIS-aided 3D localization system, comprising one multi-antenna BS,
one multi-element STAR-RIS, one single-antenna indoor MS, and one single-
antenna outdoor MS, is depicted in the unshaded area in Fig. 1. By post-
processing the pilot signals received over the uplink, the outdoor BS is
capable of localizing the two MSs simultaneously, termed as simultaneous
indoor and outdoor localization. Note that the coverage of the proposed system
depends on the BS capabilities and the area of influence of the STAR-RIS [17].
Alternatively, an indoor BS can also localize one indoor MS and one outdoor MS
with the assistance of one STAR-RIS, as depicted in the shaded region in Fig.
1. In this paper, we focus on the former scenario and leave the latter for
future investigations.

In the studied 3D localization system, the BS and STAR-RIS are assumed to be
equipped with $M$ antennas and $N$ passive scattering elements, respectively.
Without loss of generality, we further assume that both the BS and STAR-RIS
employ the uniform planar array (UPA) structure parallel to the $x\text{-}z$
plane. They can also be placed parallel to the $y\text{-}z$ plane. In this
sense, the corresponding array response vectors discussed in Section II-A need
to be modified accordingly. Recall that STAR-RIS has two operation
functionalities, i.e., reflection and refraction, which can be realized
simultaneously via two separate series of phase shifters. Therefore, the
operation of STAR-RIS can be represented by two independent control matrices,
one for controlling reflection and the other for controlling refraction.

###  II-A Channel Model

The system adopts millimeter wave frequency band for its operation thanks to
the availability of substantial spectrum. Thus, we consider the Saleh-
Valenzuela parametric channel model to construct all the four individual
channels, marked by the solid arrows in Fig. 1. The direct line-of-sight (LoS)
channel between the outdoor MS and the $M$-antenna outdoor BS is denoted as
${\mathbf{h}}_{1}\in\mathbb{C}^{M\times 1}$ and is mathematically expressed as
follows:

${\mathbf{h}}_{1}=\frac{e^{-j2\pi
d_{1}/\lambda}}{\sqrt{\rho_{1}}}\boldsymbol{\alpha}_{x}(\theta_{1},\phi_{1})\otimes\boldsymbol{\alpha}_{z}(\phi_{1}),$
(1)

where $d_{1}$ (in meters) and $\rho_{1}$ (for the sake of simplicity, we
assume that $\rho_{1}=d_{1}^{2}$) are the distance and path loss between the
outdoor MS and the outdoor BS, respectively, $\lambda$ is the wavelength of
the carrier frequency, and $\theta_{1}$ and $\phi_{1}$ are the azimuth and
elevation AoAs associated with ${\mathbf{h}}_{1}$, respectively. In the
literature, we can also consider the free-space path loss, modeled as
$\rho_{1}=d_{1}^{2}f_{c}^{2}/10^{8.755}$, where $f_{c}$ (in KHz) is the
carrier frequency, defined as $f_{c}=\frac{c}{\lambda}$ with $c$ being the
speed of light. In addition, the standard 3GPP urban micro (UMi) path loss
model can be considered, according to which holds:
$\rho_{1}=10^{2.27}d_{1}^{3.67}f_{c}^{2.6}$, where $f_{c}$ needs to be
included in GHz [18].

As the BS’s antenna array is parallel to $x\text{-}y$ plane, the array
response vectors $\boldsymbol{\alpha}_{x}(\theta_{1},\phi_{1})$ and
$\boldsymbol{\alpha}_{z}(\phi_{1})$ can be written as [19]:

$\displaystyle\boldsymbol{\alpha}_{x}(\theta_{1},\phi_{1})=$
$\displaystyle\Big{[}e^{-j\frac{2\pi
d_{x}}{\lambda}(\frac{M_{x}-1}{2})\cos(\theta_{1})\sin(\phi_{1})},$
$\displaystyle\ldots,e^{j\frac{2\pi
d_{x}}{\lambda}(\frac{M_{x}-1}{2})\cos(\theta_{1})\sin(\phi_{1})}\Big{]}^{\mathsf{T}},$
(2) $\displaystyle\boldsymbol{\alpha}_{z}(\phi_{1})=$
$\displaystyle\Big{[}e^{-j\frac{2\pi
d_{z}}{\lambda}(\frac{M_{z}-1}{2})\cos(\phi_{1})},$
$\displaystyle\ldots,e^{j\frac{2\pi
d_{z}}{\lambda}(\frac{M_{z}-1}{2})\cos(\phi_{1})}\Big{]}^{\mathsf{T}},$ (3)

where $M=M_{x}M_{z}$ with $M_{x}$ and $M_{z}$ being the number of horizontal
and vertical BS antennas, respectively, and $d_{x}$ and $d_{z}$ denote their
inter-element spacing in the horizontal and vertical axes, which are set as
half-wavelength without loss of generality. Similarly, the other two channels
linking the two MSs and the STAR-RIS, e.g.,
${\mathbf{h}}_{2}\in\mathbb{C}^{N\times 1}$ and
${\mathbf{h}}_{3}\in\mathbb{C}^{N\times 1}$, can be presented in the same
manner, as:

${\mathbf{h}}_{i}=\frac{e^{-j2\pi
d_{i}/\lambda}}{\sqrt{\rho_{i}}}\boldsymbol{\alpha}_{x}(\theta_{i},\phi_{i})\otimes\boldsymbol{\alpha}_{z}(\phi_{i}),$
(4)

for $i=2$ and $3$, where $N=N_{x}N_{z}$ is the number of STAR-RIS elements
with $N_{x}$ and $N_{z}$ denoting the numbers in the horizontal and vertical
axes, respectively. Note that the array response vectors
$\boldsymbol{\alpha}_{x}(\cdot)$ and $\boldsymbol{\alpha}_{z}(\cdot)$ in (4)
possess the same format compared to those in (1) but may differ in dimension
if $M_{x}\neq N_{x}$ and $M_{y}\neq N_{y}$. Finally, $\rho_{2}$ and $\rho_{3}$
follow the same assumption as that made for $\rho_{1}$.

The channel between the STAR-RIS and BS, i.e.,
${\mathbf{H}}_{4}\in\mathbb{C}^{M\times N}$, is expressed as follows:

${\mathbf{H}}_{4}=\frac{e^{-j2\pi
d_{4}/\lambda}}{\sqrt{\rho_{4}}}\boldsymbol{\alpha}_{x}(\theta_{4},\phi_{4})\otimes\boldsymbol{\alpha}_{z}(\phi_{4})(\boldsymbol{\alpha}_{x}(\theta_{4},\phi_{4})\otimes\boldsymbol{\alpha}_{z}(\phi_{4}))^{\mathsf{H}},$
(5)

provided that the BS and the STAR-RIS are deployed in parallel without any
biased orientation in terms of their antenna (element) arrays. As seen from
Fig. 1, there does not exist a direct LoS path between the outdoor BS and the
indoor MS due to the blockage incurred by the wall in between them. The only
path connecting them is the refraction route via the STAR-RIS. Unlike the
indoor MS, there exist one direct LoS path and one reflection path via the
STAR-RIS between the outdoor MS and the outdoor BS.

In the channel model, we consider only LoS paths for all the individual
channels; the extension for the multipath scenario are only examined in the
numerical study in Section V by adding random errors. We ignore the possible
orientations between the two UPAs (one for BS and the other for STAR-RIS),
since they can be known in priori and compensated when implementing practical
estimation algorithms for the angular parameters, i.e., azimuth and elevation
AoAs.

###  II-B Geometric Relationship

Given a pair of nodes, the geometric relationship is built between their
Cartesian coordinates and the latent channel parameters, e.g., $d_{1}$,
$\theta_{1}$, and $\phi_{1}$ in (1). The Cartesian coordinates of the BS and
STAR-RIS as well as the outdoor and indoor MSs are
${\mathbf{p}}_{\text{B}}=(x_{\text{B}},y_{\text{B}},z_{\text{B}})^{\mathsf{T}}$,
${\mathbf{p}}_{\text{R}}=(x_{\text{R}},y_{\text{R}},z_{\text{R}})^{\mathsf{T}}$,
${\mathbf{p}}_{\text{U},1}=(x_{\text{U},1},y_{\text{U},1},z_{\text{U},1})^{\mathsf{T}}$,
and
${\mathbf{p}}_{\text{U},2}=(x_{\text{U},2},y_{\text{U},2},z_{\text{U},2})^{\mathsf{T}}$,
respectively. The relationship between the distances and a pair of Cartesian
coordinates are listed below:

$\displaystyle d_{1}$
$\displaystyle=\|{\mathbf{p}}_{\text{B}}-{\mathbf{p}}_{\text{U},1}\|_{2},$ (6)
$\displaystyle d_{i}$
$\displaystyle=\|{\mathbf{p}}_{\text{R}}-{\mathbf{p}}_{\text{U},i-1}\|_{2},\;\text{for}\;i=2,3,$
(7) $\displaystyle d_{4}$
$\displaystyle=\|{\mathbf{p}}_{\text{B}}-{\mathbf{p}}_{\text{R}}\|_{2}.$ (8)

By introducing the three-element vector
$\boldsymbol{\xi}_{i}\triangleq[\cos(\theta_{i})\cos(\phi_{i}),\\\
\sin(\theta_{i})\cos(\phi_{i}),\sin(\phi_{i})]^{\mathsf{T}}$ for $i=1,2,3,4$,
the geometric relationship between the angular parameters and the Cartesian
coordinates of the nodes can be expressed as

$\displaystyle{\mathbf{p}}_{\text{R}}$
$\displaystyle={\mathbf{p}}_{\text{B}}+d_{4}\boldsymbol{\xi}_{4},$ (9)
$\displaystyle{\mathbf{p}}_{\text{U},1}$
$\displaystyle={\mathbf{p}}_{\text{B}}+d_{1}\boldsymbol{\xi}_{1}={\mathbf{p}}_{\text{R}}+d_{2}\boldsymbol{\xi}_{2},$
(10) $\displaystyle{\mathbf{p}}_{\text{U},2}$
$\displaystyle={\mathbf{p}}_{\text{R}}+d_{3}\boldsymbol{\xi}_{3}.$ (11)

The geometric relationship plays an important role in localization. According
to (10) and (11), the BS can calculate the coordinate of the MSs based on the
estimates of channel parameters ($d_{i},\theta_{i},\phi_{i}$, for $i=1,2,3$),
and the pre-known coordinates of the anchors (${\mathbf{p}}_{\text{B}}$ and/or
${\mathbf{p}}_{\text{R}}$).

###  II-C Signal Model

It is known that the STAR-RIS has two operation functionalities, i.e.,
reflection and refraction, which are realized by two separate series of phase
shifters. We introduce two phase control matrices, i.e.,
$\boldsymbol{\Omega}_{1}\in\mathbb{C}^{N\times N}$ for controlling refraction
and $\boldsymbol{\Omega}_{2}\in\mathbb{C}^{N\times N}$ for controlling
reflection, which are diagonal matrices with each diagonal element satisfying
the unit-modulus constraints, e.g.,
$|[\boldsymbol{\Omega}_{1}]_{jj}|=|[\boldsymbol{\Omega}_{2}]_{jj}|=1$,
$\forall j=1,2,\ldots,N$. However, non-ideal reflection and refraction bring
attenuation loss, results in reduced modulus, i.e.,
$|[\boldsymbol{\Omega}_{1}]_{jj}|<1$ and $|[\boldsymbol{\Omega}_{2}]_{jj}|<1$
[20]. We consider the 3D localization via the uplink transmission, where the
two users send their SRSs towards the BS in a NOMA manner. The received signal
during the $k$th time slot, for $k=1,2,\ldots,K$, can be mathematically
expressed as

${\mathbf{y}}_{k}={\mathbf{h}}_{1}x_{1,k}+\epsilon_{2}{\mathbf{H}}_{4}\boldsymbol{\Omega}_{2,k}{\mathbf{h}}_{2}x_{1,k}+\epsilon_{1}{\mathbf{H}}_{4}\boldsymbol{\Omega}_{1,k}{\mathbf{h}}_{3}x_{2,k}+{\mathbf{n}}_{k},$
(12)

where $x_{1,k}$ is the SRS from the outdoor MS, $x_{2,k}$ is the SRS from the
indoor MS, and coefficients $\epsilon_{1}$ (for refraction) and $\epsilon_{2}$
(for reflection) are used to control the power splitting for the two different
operational modes of the STAR-RIS, normalized as
$\epsilon_{1}^{2}+\epsilon_{2}^{2}=1$. The received signal at the BS is
further corrupted by the white Gaussian noise ${\mathbf{n}}_{k}$, and each
element of ${\mathbf{n}}_{k}$ follows complex Gaussian distribution
$\mathcal{CN}(0,\sigma^{2})$ with zero mean and $\sigma^{2}$ variance. During
the $k$th time slot, the refraction matrix $\boldsymbol{\Omega}_{1,k}$ and the
reflection matrix $\boldsymbol{\Omega}_{2,k}$ are considered at the STAR-RIS.
In order to ensure good estimates of the channel parameters and the locations
of the MSs, $\boldsymbol{\Omega}_{1,k}$ and $\boldsymbol{\Omega}_{2,k}$ vary
from one time slot to another, i.e.,
$\boldsymbol{\Omega}_{1,1}\neq\boldsymbol{\Omega}_{1,2}\neq\ldots\neq\boldsymbol{\Omega}_{1,K}$
and
$\boldsymbol{\Omega}_{2,1}\neq\boldsymbol{\Omega}_{2,2}\neq\ldots\neq\boldsymbol{\Omega}_{2,K}$.
The design of the refractive/reflective beam sweeping will be optimized in
Section III-C and verified through our numerical results in Section V.

The received signal vector ${\mathbf{y}}_{k}$ in (12) can be further expressed
as

$\displaystyle{\mathbf{y}}_{k}=$
$\displaystyle{\mathbf{h}}_{1}x_{1,k}+\epsilon_{2}{\mathbf{H}}_{4}\mathrm{diag}({\mathbf{h}}_{2})\boldsymbol{\omega}_{2,k}x_{1,k}$
$\displaystyle+\epsilon_{1}{\mathbf{H}}_{4}\mathrm{diag}({\mathbf{h}}_{3})\boldsymbol{\omega}_{1,k}x_{2,k}+{\mathbf{n}}_{k},$
(13)

where $\boldsymbol{\Omega}_{1,k}=\mathrm{diag}(\boldsymbol{\omega}_{1,k})$ and
$\boldsymbol{\Omega}_{2,k}=\mathrm{diag}(\boldsymbol{\omega}_{2,k})$, $\forall
k$. By stacking all ${\mathbf{y}}_{k}$’s column by column, we get the
expression:

$\displaystyle{\mathbf{Y}}=$
$\displaystyle\eta_{1}\sqrt{P}{\mathbf{h}}_{1}\mathbf{1}^{\mathsf{T}}+\eta_{1}\sqrt{P}\epsilon_{2}{\mathbf{H}}_{4}\mathrm{diag}({\mathbf{h}}_{2})\bar{\boldsymbol{\Omega}}_{2}$
$\displaystyle+\eta_{2}\sqrt{P}\epsilon_{1}{\mathbf{H}}_{4}\mathrm{diag}({\mathbf{h}}_{3})\bar{\boldsymbol{\Omega}}_{1}+{\mathbf{N}},$
(14)

where $\mathbf{1}$ denotes the $K$-element all-one vector,
${\mathbf{Y}}=[{\mathbf{y}}_{1},\ldots,{\mathbf{y}}_{K}]$,
${\mathbf{N}}=[{\mathbf{n}}_{1},\ldots,{\mathbf{n}}_{K}]$,
$\bar{\boldsymbol{\Omega}}_{1}=[\boldsymbol{\omega}_{1,1},\ldots,\boldsymbol{\omega}_{1,K}]$,
and $\bar{\boldsymbol{\Omega}}_{2}=[\boldsymbol{\omega}_{2,1},\ldots,\\\
\boldsymbol{\omega}_{2,K}]$ with
$|[\bar{\boldsymbol{\Omega}}_{1}]_{mn}|=|[\bar{\boldsymbol{\Omega}}_{2}]_{mn}|=1,\forall
m,n$. Without loss of generality, we assume that the sum transmit power
constraint is applied for each time slot, i.e.,
$|x_{1,k}|^{2}+|x_{2,k}|^{2}=P$, $\forall k$, and introduce coefficients
$\eta_{1}$ and $\eta_{2}$ for characterizing the power allocation between the
two MSs, satisfying $|x_{1,k}|^{2}=\eta_{1}^{2}P$,
$|x_{2,k}|^{2}=\eta_{2}^{2}P$, and $\eta_{1}^{2}+\eta_{2}^{2}=1$. Based on the
received signals across $K$ time slots, the BS estimates the Cartesian
coordinates of both the indoor and outdoor users, enabling 3D localization.

Applying vectorization to ${\mathbf{Y}}$ in (II-C), we get the following
expression:

$\displaystyle{\mathbf{y}}=$
$\displaystyle\eta_{1}\sqrt{P}(\mathbf{1}\otimes{\mathbf{I}}_{M}){\mathbf{h}}_{1}+\eta_{1}\sqrt{P}\epsilon_{2}(\bar{\boldsymbol{\Omega}}_{2}^{\mathsf{T}}\otimes{\mathbf{I}}_{M})({\mathbf{I}}_{N}\diamond{\mathbf{H}}_{4}){\mathbf{h}}_{2}$
$\displaystyle+\eta_{2}\sqrt{P}\epsilon_{1}(\bar{\boldsymbol{\Omega}}_{1}^{\mathsf{T}}\otimes{\mathbf{I}}_{M})({\mathbf{I}}_{N}\diamond{\mathbf{H}}_{4}){\mathbf{h}}_{3}+{\mathbf{n}},$
(15)

where ${\mathbf{y}}=\mathrm{vec}({\mathbf{Y}})$ and
${\mathbf{n}}=\mathrm{vec}({\mathbf{N}})\sim\mathcal{CN}(\boldsymbol{0},\sigma^{2}{\mathbf{I}}_{KM})$.
The expression in (II-C) can be re-written as:

${\mathbf{y}}=\sqrt{P}\eta_{1}{\mathbf{A}}_{1}{\mathbf{h}}_{1}+\sqrt{P}\eta_{1}\epsilon_{2}{\mathbf{A}}_{2}{\mathbf{h}}_{2}+\sqrt{P}\eta_{2}\epsilon_{1}{\mathbf{A}}_{3}{\mathbf{h}}_{3}+{\mathbf{n}},$
(16)

by introducing the following three new notations:

$\displaystyle{\mathbf{A}}_{1}$
$\displaystyle=(\mathbf{1}\otimes{\mathbf{I}}_{M})\in\mathbb{C}^{KM\times M},$
(17) $\displaystyle{\mathbf{A}}_{2}$
$\displaystyle=(\bar{\boldsymbol{\Omega}}_{2}^{\mathsf{T}}\otimes{\mathbf{I}}_{M})({\mathbf{I}}_{N}\diamond{\mathbf{H}}_{4})\in\mathbb{C}^{KM\times
N},$ (18) $\displaystyle{\mathbf{A}}_{3}$
$\displaystyle=(\bar{\boldsymbol{\Omega}}_{1}^{\mathsf{T}}\otimes{\mathbf{I}}_{M})({\mathbf{I}}_{N}\diamond{\mathbf{H}}_{4})\in\mathbb{C}^{KM\times
N}.$ (19)

As we can see from (17) to (19), ${\mathbf{A}}_{1}$ is independent of the
STAR-RIS design, while ${\mathbf{A}}_{2}$ and ${\mathbf{A}}_{3}$ are functions
of $\bar{\boldsymbol{\Omega}}_{2}$ and $\bar{\boldsymbol{\Omega}}_{1}$,
respectively.

Upon the STAR-RIS deployment, we assume that the BS knows the exact/precise
location of the STAR-RIS. Thus, we assume that the BS has exact information on
${\mathbf{H}}_{4}$ in terms of the parameters $\theta_{4}$, $\phi_{4}$, and
$d_{4}$. Thus, ${\mathbf{A}}_{1}$, ${\mathbf{A}}_{2}$, and ${\mathbf{A}}_{3}$
in (16) are known measurement matrices to the BS (the BS also knows the
refractive/reflective phase configurations due to its interaction with the
STAR-RIS controller) in the theoretical performance limit analyses in Section
III and localization algorithm development in Section IV. However, this
assumption will be relaxed and its effect will be examined in Section V since
perfect information on ${\mathbf{H}}_{4}$ is usually infeasible in practice.

##  III Cramér Rao Lower Bound Analyses

In this section, we summarize the CRLBs on the estimation of the intermediate
channel parameters from [10], followed by the 3D Cartesian coordinates’
estimation. This two-step approach can be commonly seen in the literature [7,
8, 21]. We also present the refraction/reflection optimization of the STAR-RIS
for the 3D localization objective, where the case $K\geq 2N+1$ is considered
and its optimal solution is found. 111In general, a reasonable training
overhead is required in order to achieve satisfactory localization performance
for both users. Thus, in this work, we only focus on the scenario with $K$
slightly larger than $2N+1$, which fits well with the aforementioned
statement. As said, we can find the optimal STAR-RIS design for such a case,
detailed in Section III-C.

###  III-A Estimation of Channel Parameters

The unknown channel parameters to be estimated are those included in
${\mathbf{h}}_{1}$, ${\mathbf{h}}_{2}$, and ${\mathbf{h}}_{3}$, i.e., the
nine-tuple
$\boldsymbol{\nu}\triangleq[\theta_{1},\phi_{1},d_{1},\theta_{2},\phi_{2},d_{2},\theta_{3},\\\
\phi_{3},d_{3}]^{\mathsf{T}}$. Since the additive noise follows complex
Gaussian distributed, by introducing
$\boldsymbol{\mu}(\boldsymbol{\nu})\triangleq\eta_{1}{\mathbf{A}}_{1}{\mathbf{h}}_{1}+\eta_{1}\epsilon_{2}{\mathbf{A}}_{2}{\mathbf{h}}_{2}+\eta_{2}\epsilon_{1}{\mathbf{A}}_{3}{\mathbf{h}}_{3}$
from (16), the Fisher information matrix for $\boldsymbol{\nu}$ is obtained
as:

$[{\mathbf{J}}(\boldsymbol{\nu})]_{i,j}=\frac{P}{\sigma^{2}}\Re\Big{\\{}\frac{\partial\boldsymbol{\mu}^{\mathsf{H}}}{\partial\nu_{i}}\frac{\partial\boldsymbol{\mu}}{\partial\nu_{j}}\Big{\\}}.$
(20)

The information on the partial derivatives in (20) related to parameters in
${\mathbf{h}}_{1}$, ${\mathbf{h}}_{2}$, and ${\mathbf{h}}_{3}$ can be referred
to [10] for more details.

For any unbiased estimator (denoted by $\hat{\boldsymbol{\nu}}({\mathbf{y}})$)
for the channel parameters, we can calculate the CRLB on the error covariance
matrix as follows:

$\mathbb{E}\\{(\boldsymbol{\nu}-\hat{\boldsymbol{\nu}}({\mathbf{y}}))(\boldsymbol{\nu}-\hat{\boldsymbol{\nu}}({\mathbf{y}}))^{\mathsf{H}}\\}\succeq{\mathbf{J}}^{-1}(\boldsymbol{\nu}),$
(21)

where the notation ${\mathbf{A}}\succeq{\mathbf{B}}$ for square matrices
${\mathbf{A}}$ and ${\mathbf{B}}$ means
${\mathbf{a}}^{\mathsf{H}}{\mathbf{A}}{\mathbf{a}}\geq{\mathbf{a}}^{\mathsf{H}}{\mathbf{B}}{\mathbf{a}}$
for any valid vector ${\mathbf{a}}$. The expression (21) indicates that the
estimation error variance for each individual channel parameter in
$\boldsymbol{\nu}$ is lower bounded by the corresponding diagonal element in
${\mathbf{J}}^{-1}(\boldsymbol{\nu})$, which indicates the best performance
any unbiased estimator can reach in theory.

###  III-B Estimation of 3D Cartesian Coordinates

Our ultimate goal is to estimate the 3D Cartesian coordinates of the indoor
and outdoor MSs. Therefore, after estimating the channel parameters, we need
to map them to 3D Cartesian coordinates, e.g.,
$\boldsymbol{\kappa}=[x_{\text{U},1},y_{\text{U},1},z_{\text{U},1},x_{\text{U},2},y_{\text{U},2},z_{\text{U},2}]^{\mathsf{T}}$,
based on the geometrical relationship among the BS, the STAR-RIS, and the two
MSs, discussed in Section II-B. For the CRLB evaluation of
$\boldsymbol{\kappa}$, we resort to the Jacobian matrix ${\mathbf{T}}$, which
links the connection between the channel parameters $\boldsymbol{\nu}$ and the
3D Cartesian coordinates of the two MSs $\boldsymbol{\kappa}$. Each $(i,j)$th
element of ${\mathbf{T}}$ is expressed as:

$[{\mathbf{T}}]_{ij}=\frac{\partial[\boldsymbol{\nu}]_{j}}{\partial[\boldsymbol{\kappa}]_{i}}.$
(22)

Again, we omit the details on the calculation of each individual derivative in
(22), which can be found in [10].

In addition, it can be easily seen that only the channel parameters
$\\{\theta_{1},\phi_{1},d_{1},\theta_{2},\phi_{2},d_{2}\\}$ are related to the
coordinates $(x_{\text{U},1},y_{\text{U},1},z_{\text{U},1})$ of the outdoor
MS, and only the parameters $\\{\theta_{3},\phi_{3},d_{3}\\}$ are related to
the coordinates $(x_{\text{U},2},y_{\text{U},2},z_{\text{U},2})$ of the indoor
MS, as concluded from (10) and (11). Therefore, the Jacobian matrix
${\mathbf{T}}$ has the following form:

${\mathbf{T}}=\begin{bmatrix}{\mathbf{T}}_{1}&\mathbf{0}\\\
\mathbf{0}&{\mathbf{T}}_{2}\end{bmatrix},$ (23)

where the submatrix ${\mathbf{T}}_{1}\in\mathbb{R}^{3\times 6}$ consists of
the partial derivatives related to the outdoor MS, and the submatrix
${\mathbf{T}}_{2}\in\mathbb{R}^{3\times 3}$ consists of the partial
derivatives related to the indoor MS.

The Fisher information of $\boldsymbol{\kappa}$ can be then expressed as [8]

${\mathbf{J}}(\boldsymbol{\kappa})={\mathbf{T}}{\mathbf{J}}(\boldsymbol{\nu}){\mathbf{T}}^{\mathsf{T}}.$
(24)

Similar to (21), we have the inequality for the CLRB:

$\mathbb{E}\\{(\boldsymbol{\kappa}-\hat{\boldsymbol{\kappa}}({\mathbf{y}}))(\boldsymbol{\kappa}-\hat{\boldsymbol{\kappa}}({\mathbf{y}}))^{\mathsf{T}}\\}\succeq{\mathbf{J}}^{-1}(\boldsymbol{\kappa}).$
(25)

The performance lower bounds on the root mean square error (RMSE) of the
position estimation of the outdoor and indoor MSs are:

$\displaystyle\text{RMSE}_{\text{U},1}=\sqrt{\mathrm{var}(\hat{{\mathbf{x}}}_{\text{U},1})}$
$\displaystyle\geq\sqrt{\mathrm{tr}\\{[{\mathbf{J}}^{-1}(\boldsymbol{\kappa})]_{1:3,1:3}\\}},$
(26)
$\displaystyle\text{RMSE}_{\text{U},2}=\sqrt{\mathrm{var}(\hat{{\mathbf{x}}}_{\text{U},2})}$
$\displaystyle\geq\sqrt{\mathrm{tr}\\{[{\mathbf{J}}^{-1}(\boldsymbol{\kappa})]_{4:6,4:6}\\}}$
(27)

where $\hat{{\mathbf{x}}}_{\text{U},1}$ and $\hat{{\mathbf{x}}}_{\text{U},2}$
are the unbiased estimates of ${\mathbf{x}}_{\text{U},1}$ and
${\mathbf{x}}_{\text{U},2}$, respectively.

###  III-C Localization-Optimal Design for the STAR-RIS

In this subsection, we consider the optimization of STAR-RIS during the pilot
transmission, aiming at maximizing the overall 3D localization performance of
the two MSs. The scenario with $K\geq 2N+1$ is considered, where we aim at
optimizing the STAR-RIS based on the inverse matrix of the Fisher information
matrix. By introducing
${\mathbf{G}}_{1}\triangleq[\frac{\partial\boldsymbol{\mu}}{\partial\theta_{1}},\frac{\partial\boldsymbol{\mu}}{\partial\phi_{1}},\frac{\partial\boldsymbol{\mu}}{\partial
d_{1}},\frac{\partial\boldsymbol{\mu}}{\partial\theta_{2}},\frac{\partial\boldsymbol{\mu}}{\partial\phi_{2}},\frac{\partial\boldsymbol{\mu}}{\partial
d_{2}}]$,
$\hat{{\mathbf{G}}}_{1}\triangleq{\mathbf{G}}_{1}{\mathbf{T}}_{1}^{\mathsf{H}}$,
${\mathbf{G}}_{2}\triangleq[\frac{\partial\boldsymbol{\mu}}{\partial\theta_{3}},\frac{\partial\boldsymbol{\mu}}{\partial\phi_{3}},\frac{\partial\boldsymbol{\mu}}{\partial
d_{3}}]$, and
$\hat{{\mathbf{G}}}_{2}\triangleq{\mathbf{G}}_{2}{\mathbf{T}}_{2}^{\mathsf{H}}$,
the expression of ${\mathbf{J}}^{-1}(\boldsymbol{\kappa})$ in (21) can be
expressed as [22, 23]

${\mathbf{J}}^{-1}\\!(\boldsymbol{\kappa})\\!=\\!\frac{\sigma^{2}}{P}\begin{bmatrix}(\hat{{\mathbf{G}}}_{1}^{\mathsf{H}}({\mathbf{I}}\\!-\\!{\mathbf{P}}_{\hat{{\mathbf{G}}}_{2}})\hat{{\mathbf{G}}}_{1})^{-1}&*\\\
*&(\hat{{\mathbf{G}}}_{2}^{\mathsf{H}}({\mathbf{I}}\\!-\\!{\mathbf{P}}_{\hat{{\mathbf{G}}}_{1}})\hat{{\mathbf{G}}}_{2})^{-1}\end{bmatrix}\\!\\!,$
(28)

where
${\mathbf{P}}_{\hat{{\mathbf{G}}}_{i}}=\hat{{\mathbf{G}}}_{i}(\hat{{\mathbf{G}}}_{i}^{\mathsf{H}}\hat{{\mathbf{G}}}_{i})^{-1}\hat{{\mathbf{G}}}_{i}^{\mathsf{H}}$
is the orthogonal projection onto the column space of $\hat{{\mathbf{G}}}_{i}$
for $i=1$ and $2$.

In order to simplify the diagonal terms in (28), we rewrite ${\mathbf{G}}_{1}$
and ${\mathbf{G}}_{2}$ as
${\mathbf{G}}_{1}=[\eta_{1}{\mathbf{A}}_{1}{\mathbf{Q}}_{1},\;\eta_{1}\epsilon_{2}{\mathbf{A}}_{2}{\mathbf{Q}}_{2}]$
and ${\mathbf{G}}_{2}=\eta_{2}\epsilon_{1}{\mathbf{A}}_{3}{\mathbf{Q}}_{3}$,
where ${\mathbf{Q}}_{1}\in\mathbb{C}^{M\times 3}$,
${\mathbf{Q}}_{2}\in\mathbb{C}^{N\times 3}$, and
${\mathbf{Q}}_{3}\in\mathbb{C}^{N\times 3}$ contain the remaining terms of
$\\{\frac{\partial\boldsymbol{\mu}}{\partial\theta_{1}},\frac{\partial\boldsymbol{\mu}}{\partial\phi_{1}},\frac{\partial\boldsymbol{\mu}}{\partial
d_{1}}\\},\\{\frac{\partial\boldsymbol{\mu}}{\partial\theta_{2}},\frac{\partial\boldsymbol{\mu}}{\partial\phi_{2}},\frac{\partial\boldsymbol{\mu}}{\partial
d_{2}}\\},\\{\frac{\partial\boldsymbol{\mu}}{\partial\theta_{3}},\frac{\partial\boldsymbol{\mu}}{\partial\phi_{3}},\\\
\frac{\partial\boldsymbol{\mu}}{\partial d_{3}}\\}$ (excluding
$\eta_{1}{\mathbf{A}}_{1}$, $\eta_{1}\epsilon_{2}{\mathbf{A}}_{2}$, and
$\eta_{2}\epsilon_{1}{\mathbf{A}}_{3}$), respectively. Our goal is to optimize
$\bar{\boldsymbol{\Omega}}_{1}$ in ${\mathbf{A}}_{3}$ and
$\bar{\boldsymbol{\Omega}}_{2}$ in ${\mathbf{A}}_{2}$ in order to obtain the
best theoretical localization performance. By further dividing
${\mathbf{T}}_{1}$ into two submatrices, i.e., as
${\mathbf{T}}_{1}=[\tilde{{\mathbf{T}}}_{1}\in\mathbb{R}^{3\times
3},\bar{{\mathbf{T}}}_{1}\in\mathbb{R}^{3\times 3}]$, we can derive the
expressions
$\hat{{\mathbf{G}}}_{1}=\eta_{1}{\mathbf{A}}_{1}{\mathbf{Q}}_{1}\tilde{{\mathbf{T}}}_{1}^{\mathsf{H}}+\eta_{1}\epsilon_{2}{\mathbf{A}}_{2}{\mathbf{Q}}_{2}\bar{{\mathbf{T}}}_{1}^{\mathsf{H}}$
and
$\hat{{\mathbf{G}}}_{2}=\eta_{2}\epsilon_{1}{\mathbf{A}}_{3}{\mathbf{Q}}_{3}{\mathbf{T}}_{2}^{\mathsf{H}}$.
For the sake of tractability for the STAR-RIS optimization (i.e.,
$\bar{\boldsymbol{\Omega}}_{1}$ and $\bar{\boldsymbol{\Omega}}_{2}$), we
maximize the principal angle (within $[0,\pi/2]$) between the subspaces of
$\hat{{\mathbf{G}}}_{1}$ and $\hat{{\mathbf{G}}}_{2}$ by following [22], which
is equivalent to minimizing
$\|\hat{{\mathbf{G}}}_{2}^{\mathrm{H}}\hat{{\mathbf{G}}}_{1}\|_{\mathrm{F}}$,
i.e.,

$\big{\|}\eta_{1}\eta_{2}\epsilon_{1}{\mathbf{T}}_{2}{\mathbf{Q}}_{3}^{\mathrm{H}}{\mathbf{A}}_{3}^{\mathrm{H}}{\mathbf{A}}_{1}{\mathbf{Q}}_{1}\tilde{{\mathbf{T}}}_{1}^{\mathsf{H}}+\eta_{1}\epsilon_{2}\eta_{2}\epsilon_{1}{\mathbf{T}}_{2}{\mathbf{Q}}_{3}^{\mathrm{H}}{\mathbf{A}}_{3}^{\mathrm{H}}{\mathbf{A}}_{2}{\mathbf{Q}}_{2}\bar{{\mathbf{T}}}_{1}^{\mathsf{H}}\big{\|}_{\mathrm{F}}.$
(29)

We focus on ${\mathbf{A}}_{3}^{\mathrm{H}}{\mathbf{A}}_{1}$ in the first term
and ${\mathbf{A}}_{3}^{\mathrm{H}}{\mathbf{A}}_{2}$ in the second term of
(29), which are detailed as

$\displaystyle{\mathbf{A}}_{3}^{\mathrm{H}}{\mathbf{A}}_{1}$
$\displaystyle=({\mathbf{I}}_{N}\diamond{\mathbf{H}}_{4}^{\mathrm{H}})(\bar{\boldsymbol{\Omega}}_{1}^{*}\mathbf{1}\otimes{\mathbf{I}}_{M}),$
(30) $\displaystyle{\mathbf{A}}_{3}^{\mathrm{H}}{\mathbf{A}}_{2}$
$\displaystyle=({\mathbf{I}}_{N}\diamond{\mathbf{H}}_{4}^{\mathrm{H}})(\bar{\boldsymbol{\Omega}}_{1}^{*}\bar{\boldsymbol{\Omega}}_{2}^{\mathrm{T}}\otimes{\mathbf{I}}_{M})({\mathbf{I}}_{N}\diamond{\mathbf{H}}_{4}),$
(31)

by using the following properties:
$({\mathbf{B}}{\mathbf{C}})^{\mathrm{H}}={\mathbf{C}}^{\mathrm{H}}{\mathbf{B}}^{\mathrm{H}}$,
$({\mathbf{B}}\times{\mathbf{C}})({\mathbf{D}}\times{\mathbf{E}})=({\mathbf{B}}{\mathbf{D}}\otimes{\mathbf{C}}{\mathbf{E}})$,
$({\mathbf{B}}\otimes{\mathbf{C}})^{\mathrm{H}}={\mathbf{B}}^{\mathrm{H}}\otimes{\mathbf{C}}^{\mathrm{H}}$,
and
$({\mathbf{B}}\diamond{\mathbf{C}})^{\mathrm{H}}={\mathbf{B}}^{\mathrm{H}}\diamond{\mathbf{C}}^{\mathrm{H}}$
for matrices ${\mathbf{B}}$, ${\mathbf{C}}$, ${\mathbf{D}}$, and
${\mathbf{E}}$ meeting the dimension requirements.

When training overhead $K\geq 2N+1$, we can easily find the optimal solution
for $\bar{\boldsymbol{\Omega}}_{1}$ and $\bar{\boldsymbol{\Omega}}_{2}$ that
minimize the Frobenius norm in (29), e.g.,
$\bar{\boldsymbol{\Omega}}_{1}=[{\mathbf{W}}]_{2:N+1,:}$ and
$\bar{\boldsymbol{\Omega}}_{2}=[{\mathbf{W}}]_{N+2:2N+1,:}$ with
${\mathbf{W}}$ being the transpose of a $K\times K$ discrete Fourier transform
(DFT) matrix or Hadamard matrix (only for certain integer $K$ values). By
following this, we have
$\bar{\boldsymbol{\Omega}}_{1}^{*}\mathbf{1}=\mathbf{0}$ in (30) and
$\bar{\boldsymbol{\Omega}}_{1}^{*}\bar{\boldsymbol{\Omega}}_{2}^{\mathrm{T}}=\mathbf{0}$
in (31), which in turn yields
$\|\hat{{\mathbf{G}}}_{2}^{\mathrm{H}}\hat{{\mathbf{G}}}_{1}\|_{\mathrm{F}}=0$
(i.e., we get the largest principal angle between the subspaces of
$\hat{{\mathbf{G}}}_{1}$ and $\hat{{\mathbf{G}}}_{2}$, which is $\pi/2$ [24]).

##  IV Localization Algorithm Development

By aligning with the theoretical analyses studied in Section III, we develop a
two-step localization algorithm to estimate the location of the two MSs. In
the first step, based on (16) one estimates the channel parameters in
${\mathbf{h}}_{1}$, ${\mathbf{h}}_{2}$, and ${\mathbf{h}}_{3}$. In the second
step, by following the geometric relationship in Section II-B, we map the
estimates of channel parameters to the coordinates of the MSs. Special
attention needs to be paid to the first step for the co-channel interference
mitigation/cancellation, i.e., null the other two when estimating one channel.
Also, for the channel parameter estimation, we resort to off-grid CS ANM for
recovering the sparsity-one channel vectors and root multiple signal
classification (root-MUSIC) for extracting the corresponding angular channel
parameters, followed by distance estimates using the least squares (LS)
principle [25]. For the location mapping of the outdoor MS, we adopt a
weighted sum criterion based on the two distance estimates associated with
${\mathbf{h}}_{1}$ and ${\mathbf{h}}_{2}$.

###  IV-A Co-Channel Interference Cancellation

As an essential process for channel parameter estimation, we need to null the
interference for each MS in (16). By multiplying the interference nulling
matrix ${\mathbf{U}}_{i}$ to ${\mathbf{y}}$, for $i=1,2,3$, we get the
following

${\mathbf{U}}_{i}{\mathbf{y}}=\gamma_{i}{\mathbf{U}}_{i}{\mathbf{A}}_{i}{\mathbf{h}}_{i}+{\mathbf{U}}_{i}{\mathbf{n}},$
(32)

where $\gamma_{1}=\sqrt{P}\eta_{1}$,
$\gamma_{2}=\sqrt{P}\eta_{1}\epsilon_{2}$, and
$\gamma_{3}=\sqrt{P}\eta_{2}\epsilon_{1}$. We rely on (32) in the next
subsection for channel parameter estimation. The co-channel interference can
be completely cancelled in (32) in certain situations depending on the
relationship among the vector spaces spanned by the column vectors of
${\mathbf{A}}_{1}$, ${\mathbf{A}}_{2}$, and ${\mathbf{A}}_{3}$, respectively
[26]. Otherwise, we merge the residual interference to the additive noise term
${\mathbf{U}}_{i}{\mathbf{n}}$. One way to design ${\mathbf{U}}_{1}$ is to
perform singular value decomposition (SVD) on the combined matrix
$[{\mathbf{A}}_{2},{\mathbf{A}}_{3}]$ and take the left singular vectors
corresponding to zero singular values [26]. A similar approach can be applied
for constructing ${\mathbf{U}}_{2}$ and ${\mathbf{U}}_{3}$.

###  IV-B Channel Parameter Estimation

All the channel vectors in (16) are in the form of Kronecker product of two
array response vectors. It can be seen as sparsity-one signal (i.e., a linear
combination of one atom having the same structure as the channel component in
the form of Kronecker product of two array response vectors), which is an
extreme case under the framework of CS. Therefore, we follow the off-grid CS
techniques, i.e., ANM, for recovering the sparsity-one channel vector. This
technique has been used in the literature [27, 25] for channel estimation and
localization purposes. We first introduce the atomic set [28, 29, 27], as

$\mathcal{A}\triangleq\\{\boldsymbol{\alpha}_{x}(x_{1},x_{2})\otimes\boldsymbol{\alpha}_{z}(x_{2}),x_{1}\in[0,\pi],x_{2}\in[-\pi/2,\pi/2]\\},$
(33)

where each atom possesses the same structure with the linear term in
${\mathbf{h}}_{i}$, for $i=1,2,3$. For any vector ${\mathbf{h}}_{i}$ of the
form
${\mathbf{h}}_{i}=\sum_{l}\eta_{l}\boldsymbol{\alpha}_{x}(x_{1,l},x_{2,l})\otimes\boldsymbol{\alpha}_{z}(x_{2,l})$
with each $\eta_{l}>0$ being a coefficient, $x_{1,l}\in[0,\pi]$, and
$x_{2,l}\in[-\pi/2,\pi/2]$, its atomic norm with respect to the atomic set
$\mathcal{A}$ is written as follows:

$\displaystyle\|{\mathbf{h}}_{i}\|_{\mathcal{A}}=$
$\displaystyle\mathrm{inf}_{\mathcal{B}}\Big{\\{}\frac{1}{2T_{i}}\mathrm{Tr}(\mathrm{Toep}(\mathcal{U}_{2}))+\frac{t}{2}\Big{\\}},$
$\displaystyle\text{s.t.}\;\begin{bmatrix}\mathrm{Toep}(\mathcal{U}_{2})&{\mathbf{h}}_{i}\\\
{\mathbf{h}}_{i}^{\mathsf{H}}&t_{i}\end{bmatrix}\succeq\mathbf{0},$ (34)

where set $\mathcal{B}\triangleq\\{\mathcal{U}_{2}\in\mathbb{C}^{T_{i}\times
T_{i}},t_{i}\in\mathbb{R}\\}$ with $\mathcal{U}_{2}$ being a $2$-way tensor
and $\mathrm{Toep}(\mathcal{U}_{2})$ is a $2$-level block Toeplitz matrix,
which results from the Vandermonde decomposition lemma for positive
semidefinite Toeplitz matrices [29]. The value for $T_{i}$ depends on the
dimension of ${\mathbf{h}}_{i}$, i.e., $T_{i}=M$ for $i=1$, and $T_{i}=N$, for
$i=2,3$.

The ANM based channel estimation can be formulated as a regularized
optimization problem:

$\displaystyle\hat{{\mathbf{h}}}_{i}=\arg\min_{{\mathbf{h}}_{i}\in\mathbb{C}^{T_{i}},\;\mathcal{B}}$
$\displaystyle\mu_{i}\|{\mathbf{h}}_{i}\|_{\mathcal{A}}+\frac{1}{2}\|{\mathbf{U}}_{i}{\mathbf{y}}-\gamma_{i}{\mathbf{U}}_{i}{\mathbf{A}}_{i}{\mathbf{h}}_{i}\|_{2}^{2}$
$\displaystyle\text{s.t.}\;\begin{bmatrix}\mathrm{Toep}(\mathcal{U}_{2})&{\mathbf{h}}_{i}\\\
{\mathbf{h}}_{i}^{\mathsf{H}}&t_{i}\end{bmatrix}\succeq\mathbf{0},$ (35)

where $\mu_{i}\propto\sigma\sqrt{T_{i}\log(T_{i})}$ is the regularization term
of the atomic norm penalty, and $\hat{{\mathbf{h}}}_{i}$ is the estimate of
${\mathbf{h}}_{i}$. This problem can be efficiently solved using the Matlab
CVX toolbox.

Based on the $\hat{{\mathbf{h}}}$, the elevation and azimuth AoAs can be
extracted by following root-MUSIC algorithm [25] and the distance $d_{i}$ can
be estimated by following LS principle, as

$\hat{d}_{i}=\sqrt{T_{i}/{\mathbf{h}}_{i}^{\mathsf{H}}{\mathbf{h}}_{i}},$ (36)

where $\hat{d}_{i}$ is the estimate of $d_{i}$.

###  IV-C Location Mapping

Based on the estimate of ${\mathbf{h}}_{i}$ in (IV-B), we can further resort
to root-MUSIC for extracting the angular parameters, denoted by
$\hat{\theta}_{i}$ and $\hat{\phi}_{i}$ [25]. Together with the estimate of
$d_{i}$ in (36), the location of the MSs can be calculated by following the
geometric relationship in Section II-B.

Specifically, for the location estimate of the outdoor MS, since both
${\mathbf{h}}_{1}$ and ${\mathbf{h}}_{2}$ contribute to it, we apply weighted
sum principle, as

${\mathbf{p}}_{\text{U},1}=w_{1}({\mathbf{p}}_{\text{B}}+\hat{d}_{1}\hat{\boldsymbol{\xi}}_{1})+(1-w_{1})({\mathbf{p}}_{\text{R}}+\hat{d}_{2}\hat{\boldsymbol{\xi}}_{2}),$
(37)

where the weight
$w_{1}=\frac{\hat{d}_{2}^{2}}{\hat{d}_{1}^{2}+\hat{d}_{2}^{2}}$ is set in a
heuristic way by following that the weight is reversely proportional to the
path loss.

##  V Numerical Results

In this section’s numerical investigation, we set the system parameters as
follows: ${\mathbf{p}}_{\text{B}}=(0,0,8)^{\mathsf{T}}$,
${\mathbf{p}}_{\text{R}}=(2,2,5)^{\mathsf{T}}$,
${\mathbf{p}}_{\text{U},1}=(5,1,2)^{\mathsf{T}}$, and
${\mathbf{p}}_{\text{U},2}=(1,5,2)^{\mathsf{T}}$. The numbers of BS antennas,
STAR-RIS elements, and SRSs from each MS are set as $M=16$, $N=36$, and
$K=100,130$. The signal-to-noise ratio (SNR) is defined as $P/\sigma^{2}$. The
parameter setup is summarized in Table I.

TABLE I: Parameter Setup. Parameter | Value | Parameter | Value
---|---|---|---
$M$ | $16$ | $N$ | $36$
$K$ | $100,130$ | ${\mathbf{p}}_{\text{B}}$ | $(0,0,8)^{\mathsf{T}}$
${\mathbf{p}}_{\text{R}}$ | $(2,2,5)^{\mathsf{T}}$ | ${\mathbf{p}}_{\text{U},1}$ | $(5,1,2)^{\mathsf{T}}$
${\mathbf{p}}_{\text{U},2}$ | $(1,5,2)^{\mathsf{T}}$ |  |

###  V-A Effect of Training Overhead

As shown in Section III-C, we find the optimal design of the STAR-RIS for the
training overhead $K\geq 2N+1$. We pick up two $K$ values meeting this
requirement and compare their impact on the localization performance. The
simulation results, including both theoretical and practical, are shown in
Fig. 2 for the training overhead $K=100,130$, where
$\epsilon_{1}=\sqrt{0.9},\eta_{1}=\sqrt{0.5}$. In the legend, “ANM” denotes
the proposed ANM based 3D localization scheme while “Theo” stands for the
theoretical performance limit analyses. From the theoretical ones
characterized by CRLBs, we know that higher training overhead can bring better
localization performance, up to centimeter level for both MSs. Moreover, the
indoor MS can achieve better performance than the outdoor MS in such a
unbalanced setup on $\epsilon_{1}$ and $\eta_{1}$ since the STAR-RIS power
splitting coefficient is large (Namely, more percent of energy is refracted
towards the BS via the STAR-RIS). The practical results from ANM are
consistent with the theoretical studies. The performance gain brought by
increasing the training overhead from $100$ to $130$ is not so obvious from
the practical results, especially in the low SNR regime. However, a constant
gain (roughly $4$ dB) is observed from the theoretical CRLB results across all
the SNR values for both users.

Figure 2: The effect of training overhead on 3D localization with
$\epsilon_{1}=\sqrt{0.9}$ and $\eta_{1}=\sqrt{0.5}$.

###  V-B Effect of Power Splitting and Allocation

The performance of 3D localization is not only affected by the training
overhead but also the two parameters, i.e., $\epsilon_{1}$ for controlling
power splitting at STAR-RIS and $\eta_{1}$ for controlling power allocation
between the two users. For the training overhead $K=100$, the simulation
results with different setups on $\epsilon_{1}$ and $\eta_{1}$ are shown in
Fig. 3. As we can see that when a balanced setup is adopted, i.e.,
$\epsilon_{1}=\eta_{1}=\sqrt{0.5}$, the performance gap between the two MSs’
position estimation is small. However, in the other setup, i.e.,
$\epsilon_{1}=\sqrt{0.9},\eta_{1}=\sqrt{0.5}$, the gap is obvious. By
carefully choosing the values for the two parameters, we can simultaneously
achieve the QoS requirements for both users.

Figure 3: 3D localization performance with different pairs of $\epsilon_{1}$
and $\eta_{1}$, where the training overhead is set as $K=100$.

###  V-C Effect of STAR-RIS Design

In this subsection, we evaluate the effect of the STAR-RIS design (under the
condition of training overhead $K=100$) by considering two different cases: i)
$\bar{\boldsymbol{\Omega}}_{1}$ and $\bar{\boldsymbol{\Omega}}_{2}$ are
designed according Section III-C; ii) the phases of
$\bar{\boldsymbol{\Omega}}_{1}$ and $\bar{\boldsymbol{\Omega}}_{2}$ are
randomly generated. The simulation results are shown in Fig. 4. As expected,
the first case outperforms the second one, since according to Section III-C it
is optimal. The performance gain between the two cases is obvious with
different setups for $\epsilon_{1}$ and $\eta_{1}$. For instance, when
$\epsilon_{1}=\sqrt{0.9}$ and $\eta_{1}=\sqrt{0.5}$, more than $5$ dB gain in
terms of SNR can be obtained for the indoor MS by following the optimal STAR-
RIS design compared to the random phase design.

Figure 4: The effect of STAR-RIS design on 3D localization with training
overhead $K=100$ and different pairs of $\epsilon_{1}$ and $\eta_{1}$.

###  V-D Effect of Imperfect STAR-RIS-to-BS Channel

In this section, we examine the effect of imperfectness of the STAR-RIS-to-BS
channel on the 3D localization performance. Different from the previous
subsections, here, we assume that the STAR-RIS-to-BS channel matrix is
available but in the form of imperfectness during the localization process. We
introduce individual random variation to each of the channel parameters in the
STAR-RIS-to-BS channel, e.g., $d_{4}$, $\theta_{4}$, $\phi_{4}$. Such
variations can be introduced by practical algorithms for estimating them. We
further assume that the variations follow uniform distribution, i.e., $\Delta
d_{4}\sim\mathcal{U}[-\hat{d},\hat{d}]$,
$\Delta\theta_{4},\Delta\phi_{4}\sim\mathcal{U}[-\hat{\varphi},\hat{\varphi}]$.
The values of $\hat{d}$ and $\hat{\varphi}$ jointly determine the level of
imperfectness of the STAR-RIS-to-BS channel. In this experiment, we evaluate
two cases: i) $\hat{d}=0.5$ in meter and $\hat{\varphi}=0.2$ in radian, ii)
$\hat{d}=1$ in meter and $\hat{\varphi}=0.4$ in radian, under training
overhead $K=100$. The simulation results are provided in Fig. 5, where
$\epsilon_{1}=\sqrt{0.9}$ and $\eta_{1}=\sqrt{0.5}$. As expected, the
localization performance degrades when $\hat{d}$ and $\hat{\varphi}$
increases, especially in the high SNR regime.

Figure 5: The effect of imperfect STAR-RIS-to-BS channel on 3D localization
with training overhead $K=100$ and different pairs of $\hat{d}$ and
$\hat{\varphi}$.

###  V-E Effect of MPCs

The proposed localization algorithm purely relies channel parameters
associated with the LoS path. Therefore, the strength of the non-line-of-sight
(NLoS) path will negatively affect the localization performance. In principle,
the stronger the strength, the worse the localization performance. In this
subsection, we pick up different setups on the average sum power of the NLoS
paths (revealed by the distance) and evaluate accordingly their negative
effect on the localization performance. We introduce two MPCs to
${\mathbf{h}}_{1}$, ${\mathbf{h}}_{2}$, and ${\mathbf{h}}_{3}$ while keeping
perfect LoS condition for ${\mathbf{H}}_{4}$. The introduced MPCs for
${\mathbf{h}}_{i}$, for $i=1,2,3$, comprise the following channel parameters:
case i) $\\{10d_{i},\theta_{i}+\pi/6,\phi_{i}+\pi/6\\}$ and
$\\{10d_{i},\theta_{i}+\pi/3,\phi_{i}+\pi/3\\}$ and case ii)
$\\{5d_{i},\theta_{i}+\pi/6,\phi_{i}+\pi/6\\}$ and
$\\{5d_{i},\theta_{i}+\pi/3,\phi_{i}+\pi/3\\}$. The simulation results are
shown in Fig. 6, where we observe that the stronger the MPCs, the worse the
localization performance. Even though there are works on leveraging NLoS paths
for enhancing localization performance [30, 31], we will leave such an
extension for our future investigation.

Figure 6: The effect of MPCs on 3D localization with different pairs of
$\hat{d}$ and $\hat{\varphi}$ with $K=100$, $\epsilon_{1}=\sqrt{0.9}$, and
$\eta_{1}=\sqrt{0.5}$.

##  VI Conclusion and Future Work

In this paper, we studied the fundamental 3D localization performance limits
of STAR-RIS-empowered millimeter wave MIMO systems for simultaneously serving
one indoor and one outdoor MSs. We presented a practical localization
algorithm based on ANM, which approaches the theoretical performance limits
characterized by our derived CRLBs. In addition, we thoroughly investigated
the effect of training overhead, energy splitting at the STAR-RIS, the power
allocation between the two MSs, the STAR-RIS design, the imperfectness of
STAR-RIS-to-BS channel, as well as the role of the presence of MPCs on the
localization performance of the two MSs, offering some useful insights for
future practical implementations. In future works, we will extend the
theoretical CRLB analyses to general multipath scenarios and exploit the
availability of MPCs for further enhancing the localization performance of the
two MSs in the design of practical localization algorithms.

## References

  * [1] C. Huang, A. Zappone, G. C. Alexandropoulos, M. Debbah, and C. Yuen, “Reconfigurable intelligent surfaces for energy efficiency in wireless communication,” _IEEE Trans. Wireless Commun._ , vol. 18, no. 8, pp. 4157–4170, Aug. 2019.
  * [2] M. Di Renzo et al., “Smart radio environments empowered by reconfigurable AI meta-surfaces: an idea whose time has come,” _EURASIP J. Wireless Comm. and Networking_ , vol. 2019, p. 129, 2019.
  * [3] G. C. Alexandropoulos, G. Lerosey, M. Debbah, and M. Fink, “Reconfigurable intelligent surfaces and metamaterials: The potential of wave propagation control for 6G wireless communications,” _IEEE ComSoc TCCN Newslett._ , vol. 6, no. 1, pp. 25–37, Jun. 2020.
  * [4] Q. Wu, S. Zhang, B. Zheng, C. You, and R. Zhang, “Intelligent reflecting surface aided wireless communications: A tutorial,” _IEEE Trans. Commun._ , vol. 69, no. 5, pp. 3313–3351, May 2021.
  * [5] E. Calvanese Strinati, G. C. Alexandropoulos, H. Wymeersch, B. Denis, V. Sciancalepore, R. D’Errico, A. Clemente, D.-T. Phan-Huy, E. D. Carvalho, and P. Popovski, “Reconfigurable, intelligent, and sustainable wireless environments for 6G smart connectivity,” _IEEE Commun. Mag._ , vol. 59, no. 10, pp. 99–105, Oct. 2021.
  * [6] H. Wymeersch, J. He, B. Denis, A. Clemente, and M. Juntti, “Radio localization and mapping with reconfigurable intelligent surfaces: Challenges, opportunities, and research directions,” _IEEE Veh. Technol. Mag._ , vol. 15, no. 4, pp. 52–61, Dec. 2020.
  * [7] J. He, H. Wymeersch, L. Kong, O. Silvén, and M. Juntti, “Large intelligent surface for positioning in millimeter wave MIMO systems,” in _Proc. IEEE VTC-Spring_ , 2020, pp. 1–5.
  * [8] A. Elzanaty, A. Guerra, F. Guidi, and M.-S. Alouini, “Reconfigurable intelligent surfaces for localization: Position and orientation error bounds,” _IEEE Trans. Signal Process._ , vol. 69, Aug. 2021.
  * [9] G. C. Alexandropoulos, I. Vinieratou, and H. Wymeersch, “Localization via multiple reconfigurable intelligent surfaces equipped with single receive RF chains,” _IEEE Wireless Commun. Lett._ , vol. 11, no. 5, pp. 1072–1076, 2022.
  * [10] J. He, A. Fakhreddine, and G. C. Alexandropoulos, “Simultaneous indoor and outdoor 3D localization with STAR-RIS-assisted millimeter wave systems,” in _proc. IEEE Vehicular Technology Conference (Fall)_ , September, 2022, pp. 1–6.
  * [11] M. Jian, G. C. Alexandropoulos, E. Basar, C. Huang, R. Liu, Y. Liu, and C. Yuen, “Reconfigurable intelligent surfaces for wireless communications: Overview of hardware designs, channel models, and estimation techniques,” _Int. Conv. Netw._ , vol. 3, no. 1, pp. 1–32, Mar. 2022.
  * [12] G. C. Alexandropoulos, N. Shlezinger, and P. del Hougne, “Reconfigurable intelligent surfaces for rich scattering wireless communications: Recent experiments, challenges, and opportunities,” _IEEE Commun. Mag._ , vol. 59, no. 6, pp. 28–34, Jun. 2021.
  * [13] G. C. Alexandropoulos, N. Shlezinger, I. Alamzadeh, M. F. Imani, H. Zhang, and Y. C. Eldar, “Hybrid reconfigurable intelligent metasurfaces: Enabling simultaneous tunable reflections and sensing for 6G wireless communications,” 2021. [Online]. Available: https://arxiv.org/pdf/2104.04690
  * [14] Y. Liu, X. Mu, J. Xu, R. Schober, Y. Hao, H. V. Poor, and L. Hanzo, “STAR: Simultaneous transmission and reflection for 360∘ coverage by intelligent surfaces,” _IEEE Wireless Commun._ , vol. 28, no. 6, pp. 102–109, Dec. 2021.
  * [15] C. Wu, C. You, Y. Liu, X. Gu, and Y. Cai, “Channel estimation for star-ris-aided wireless communication,” _IEEE Commun. Lett._ , vol. 26, no. 3, pp. 652–656, 2022.
  * [16] X. Yue, J. Xie, Y. Liu, Z. Han, R. Liu, and Z. Ding, “Simultaneously transmitting and reflecting reconfigurable intelligent surface assisted noma networks,” _IEEE Trans. Wireless Commun._ , pp. 1–1, 2022.
  * [17] G. C. Alexandropoulos, M. Crozzoli, D.-T. Phan-Huy, K. D. Katsanos, H. Wymeersch, P. Popovski, P. Ratajczak, Y. Bénédic, M.-H. Hamon, S. Herraiz Gonzalez, R. D’Errico, and E. Calvanese Strinati, “Smart wireless environments enabled by riss: Deployment scenarios and two key challenges,” in _Proc. Joint European Conference on Networks and Communications & 6G Summit_, Grenoble, France, Jun. 2022, pp. 1–6.
  * [18] M. R. Akdeniz, Y. Liu, M. K. Samimi, S. Sun, S. Rangan, T. S. Rappaport, and E. Erkip, “Millimeter wave channel modeling and cellular capacity evaluation,” _IEEE J. Sel. Areas Commun._ , vol. 32, no. 6, pp. 1164–1179, Jun. 2014.
  * [19] Y. Tsai, L. Zheng, and X. Wang, “Millimeter-wave beamformed full-dimensional MIMO channel estimation based on atomic norm minimization,” _IEEE Trans. Commun._ , vol. 66, no. 12, pp. 6150–6163, Dec. 2018.
  * [20] Q. Wu and R. Zhang, “Beamforming optimization for wireless network aided by intelligent reflecting surface with discrete phase shifts,” _IEEE Trans. Commun._ , pp. 1–1, 2019.
  * [21] A. Shahmansoori, G. E. Garcia, G. Destino, G. Seco-Granados, and H. Wymeersch, “Position and orientation estimation through millimeter-wave MIMO in 5G systems,” _IEEE Trans. Wireless Commun._ , vol. 17, no. 3, pp. 1822–1835, 2018.
  * [22] L. L. Scharf and L. T. McWhorter, “Geometry of the Cramer-Rao bound,” _Signal Process._ , vol. 31, no. 3, pp. 301–311, Apr. 1993.
  * [23] P. Pakrooh, A. Pezeshki, L. L. Scharf, D. Cochran, and S. D. Howard, “Analysis of Fisher information and the Cramér–Rao bound for nonlinear parameter estimation after random compression,” _IEEE Trans. Signal Process._ , vol. 63, no. 23, pp. 6423–6428, Dec. 2015.
  * [24] G. H. Golub and C. F. Van Loan, _Matrix computations_. JHU press, 2013.
  * [25] J. He, A. Fakhreddine, C. Vanwynsberghe, H. Wymeersch, and G. C. Alexandropoulos, “3D localization with a single partially-connected receiving RIS: Positioning error analysis and algorithmic design,” _arXiv preprint arXiv:2212.02088_ , 2022.
  * [26] H. Sung, S.-H. Park, K.-J. Lee, and I. Lee, “Linear precoder designs for k-user interference channels,” _IEEE Trans. Wireless Commun._ , vol. 9, no. 1, pp. 291–301, 2010.
  * [27] J. He, H. Wymeersch, and M. Juntti, “Channel estimation for RIS-aided mmWave MIMO systems via atomic norm minimization,” _IEEE Trans. Wireless Commun._ , pp. 1–1, 2021.
  * [28] Z. Yang and L. Xie, “Exact joint sparse frequency recovery via optimization methods,” _IEEE Trans. Signal Process._ , vol. 64, no. 19, pp. 5145–5157, Oct 2016.
  * [29] G. Tang, B. N. Bhaskar, P. Shah, and B. Recht, “Compressed sensing off the grid,” _IEEE Trans. Inf. Theory_ , vol. 59, no. 11, pp. 7465–7490, 2013.
  * [30] K. Witrisal, P. Meissner, E. Leitinger, Y. Shen, C. Gustafson, F. Tufvesson, K. Haneda, D. Dardari, A. F. Molisch, A. Conti, and M. Z. Win, “High-accuracy localization for assisted living: 5G systems will turn multipath channels from foe to friend,” _IEEE Signal Process. Mag._ , vol. 33, no. 2, pp. 59–70, 2016.
  * [31] R. Mendrzik, H. Wymeersch, G. Bauch, and Z. Abu-Shaban, “Harnessing NLOS components for position and orientation estimation in 5G millimeter wave MIMO,” _IEEE Trans. Wireless Commun._ , vol. 18, no. 1, pp. 93–107, 2019\.


# A Meta-MARL approach for data collection and trajectory optimizations for
UAVs assisted networks

Mouhamed Naby Ndiaye  School of Computer Science, Mohammed VI Polytechnic
University (UM6P), Benguerir, Morocco,
emails<EMAIL_ADDRESS>El Houcine
Bergou  School of Computer Science, Mohammed VI Polytechnic University (UM6P),
Benguerir, Morocco,
emails<EMAIL_ADDRESS>Hajar El
Hammouti  School of Computer Science, Mohammed VI Polytechnic University
(UM6P), Benguerir, Morocco,
emails<EMAIL_ADDRESS>
###### Abstract

Unmanned aerial vehicle (UAV) is seen as a promising technology to accomplish
a plethora of tasks in many fields. In telecommunications in particular, UAVs
are often used as base stations or as data collectors. In this work we
consider a deployment of a group of UAVs to collect data as fresh as possible
generated by IoT devices on the ground. We formulate our problem as a multi-
objective optimization problem where the main objectives are to maximize the
total data collected and minimize the Age-of-Update of the system. In order to
efficiently solve our optimization problem, we propose a cooperative Multi
Agent Reinforcement Learning (MARL) algorithm based on the popular on-policy
RL algorithm Policy Proximal Optimization (PPO). We also propose a meta-
training procedure that allows UAVs to generalize their knowledge and quickly
adapt to a new environment. Simulation results show that our proposed
algorithm performs very well compared to the traditional off-policy MARL
algorithm and, with our meta-procedure, UAVs are able to adapt quickly to
unknown environments.

##  1 Introduction

Nowadays, Unmanned Aerials Vehicles (UAVs), due to their flexibility and
manoeuvrability, offer a plethora of possibilities in the field of
telecommunication [bajracharya20226g]. They can be used as base stations or as
data collectors. Their agility allows them to fly close to Users Equipment
(UEs), collect the data they generate and relay it to the Base Station (BS).
However, despite their agility, UAVs need their trajectory to be efficiently
designed to carry out their mission. Moreover, drones, like all machines, have
no intelligence and are therefore not autonomous. Furthermore, as the data
generated by UEs are time sensitive and may lose their value after a certain
period of time, it is crucial to design the trajectory of UAVs in such a way
as to collect the freshest possible information. Therefore, in this work, we
present a practical framework for the deployment of autonomous UAVs capable of
adapting to previously unknown environments working cooperatively to achieve
several goals, which are maximising the total data collected and minimising
the age of the collected data information captured by a metric called Age-of-
Updtae (AoU).

###  1.1 Related Work

Many problems related to UAV-assisted networks have been studied in the
existing literature [Wang2020EnergyEfficientDC, Yuan2021NovelOT,
Wang2022JointOO, Jiang2021PeakAO, Cao2021PeakAO]. On the one hand, work has
focused on the design of trajectories for efficient data collection and on the
other hand on the freshness of the information collected, often quantified by
a metric called Age-of-Information (AoI). In the paper
In[Wang2020EnergyEfficientDC], The authors propose that the UAV trajectory and
the device transmission schedule be jointly optimized, reducing the maximum
energy consumption of all devices while ensuring the reliability of the
collected data. The authors of this paper defined their problem as a non-
convex mixed integer optimisation problem and solved it effectively using a
differential evolution approach.

However, many of these works are limited to a single UAV and none of them
design UAV trajectories in cooperative mode. In this perspective, a handful of
works have been done on the cooperation of UAVs to accomplish their mission,
usually based on multi-agent reinforcement learning (MARML)
[Cui2020MultiAgentRL, Seid2021MultiAgentDF]. The authors
in[Cui2020MultiAgentRL] formulate their problem as a stochastic game to
maximise expected rewards and propose a MARL-based solution to solve it.
Similarly in[Seid2021MultiAgentDF] formulate their problem as a natural
extension of the Markov decision process (MDP) concerning stochastic game, to
minimize the long-term computation cost in terms of energy and delay, and
propose a multi-agent deep reinforcement learning (MADRL)-based approach to
ensure the quality of service (QoS) requirements of IoT devices while
minimizing the overall network computation cost.

Nevertheless, most existing MARL-based solutions as those proposed, agents can
search for strategies independently based on their own actions and states, as
shown in the MARL-based solutions proposed in[Hu2020CooperativeIO,
Lin2021DecentralizedPD]. However, as agents act in a decentralised manner and
only optimise their individual goals, this type of solution will not allow
drones to maximise the total data collected nor minimise the overall AoU of
the system. In contrast, the solutions proposed in [Hu2021DistributedMM,
Hu2022DistributedAD] allow agents to share neither their states, nor their
rewards, nor their strategies, while achieving the global goals. However in
the real life, particularly in tasks which need collaboration, agents need to
communicate for achieving their mission. These types of solutions that allow
agent to share knowledge when they learn and act independently, are generally
based on the Centralised Learning and Decentralised Execution
(CTDE)[Wong2021MultiagentDR] paradigm , where agents learn centrally and act
independently.

Furthermore, in general, within the MARL paradigm, Policy Gradient (PG)
methods are generally considered less efficient in terms of sampling than
Value Decomposition (VD) methods [Hu2021DistributedMM], which are not policy-
based, due to their nature. However, some recent empirical studies
[Witt2020IsIL, Papoudakis2021BenchmarkingMD] show that with appropriate hyper-
parametric parameters and input representation, PG method can achieve
surprisingly high performance compared to non-policy based VD methods. These
can learn a blending function that blends local Q nets into a global Q
function. To satisfy the Individual-Global-Max (IGM) principle, which ensures
that greedily choosing the optimal action locally for each agent, the optimal
joint action can be computed, this mixture function is usually applied. In
contrast, PG methods directly apply the policy gradient to learn a centralized
value function for each agent and an individual policy. For an accurate
estimate of the global value, the value function can take as input the
concatenation of all local observations, or the global state.

Moreover, the traditional solutions proposed, based on RL, do not allow agents
to adapt quickly and efficiently to unknown environments, as they are often
too adapted to the tasks they have been trained for. Indeed, in order for
these algorithms to adapt to learning tasks, hyper-parameters and
initialization are adapted manually. Thus, when the agents are confronted with
a new environment, in other words an unknown task, these algorithms converge
very slowly or may even fail to converge. In contrast, in the work
[Hu2021DistributedMM, Zhou2022MultiagentFM], the proposed solutions allow
agents to easily adapt to unknown environments through model-agnostic meta-
learning. By meta-training agents on multiple tasks and not on a single task
as in traditional RL solutions, the solutions proposed in
[Hu2021DistributedMM, Zhou2022MultiagentFM] allow rapid adaptation of drones
to unknown environments.

Finally, none of these works consider a multi-objective optimization problem.
In most of these works, agents learn to achieve a single goal. In this work,
we propose a framework in which agents learn to achieve several goals at the
same time in a collaborative way. Our proposed MARL algorithm use the actor
critic architecture, where the critic is central (centralized value function)
and the actor is decentralized (individual policy for each UAVs). Furthermore,
we propose a meta-training procedure that allows UAVs to adapt quickly to
unknown environments.

###  1.2 Contribution

The major contribution of this paper is a new distributed framework for
designing the trajectories of a swarm of cooperative UAVs in highly dynamic
and especially unpredictable environments. Our major contributions can be
listed in the points below:

  * •

We consider a practical UAV-assisted network system in which a team of UAVs is
deployed to act as a relay between low-capacity user equipment and remote BS
in order to collect as much fresh data generated by these user equipment on
the ground as possible, despite strict mission time constraints and limited
information on their environments. The UAVs can collect data generated by UEs
while adapting their trajectory to minimize the system’s AoU. The problem is
formulated as a multi-objective optimization problem, which turns out to be
non-convex. To solve this formulated non-convex multi-objective optimization
problem, we propose a distributed reinforcement learning algorithm based on
the popular single-agent reinforcement learning algorithm called Proximal
Policy Optimization (PPO). The proposed algorithm, MAPPO-TO-AoU, is built on
the Central Training and Decentralized Execution (CTDE) paradigm based on the
actor-critic architecture, where a critic learns a central value function and
an actor finds the local optimal policy. This proposed MAPPO-TO-AoU algorithm
allows UAVs to find a local optimal solution to the studied multi-objective
problem using only their own local states information.

  * •

In order to find an optimal policy and value function initialization with
appropriate estimation over a multitude of environments, we propose to meta-
train our proposed MAPPO-TO-AoU algorithm using model-agnostic meta-learning
(MAML). This meta-training procedure improves the convergence speed of the
MAPPO-TO-AoU algorithm for unseen environments. The advantage of this meta-
training procedure is that it is less complex than other meta-training
procedures [ritter2018been] because it does not require additional neural
networks. Thus, using this meta-training approach, drones are able to cope
with various environments (handle various tasks) instead of just one, due to
the use of meta-training initialization. The simulation results show that the
proposed MAPPO-TO-AoU algorithm requires fewer iterations to converge and
obtain higher rewards compared to the existing traditional off-policy MARL
solutions used as benchmark. In addition with MAPPO-TO-AoU, UAVs collect more
data, quickly minimize the AoU of the system and serve more UEs. The results
also show that our meta-trained MAPPO-TO-AoU improves the MAPPO-TO-AoU
algorithm by converging faster and achieving even higher rewards and data
collected.

###  1.3 Organization

The remainder of the paper is organized as follows. First, the studied system
model is described in Section 2. The mathematical formulation of the problem
is given in Section 3. In Section LABEL:Algo, we describe the joint
association probability and 3D locations to solve the underlying optimization
efficiently. Next, in Section LABEL:Simu, we show the performance of our
algorithm for various scenarios. Finally, concluding remarks are provided in
Section 6.

##  2 System Model

We consider a system concerning of set $\mathcal{I}$ of $I$ UEs and a set
$\mathcal{U}$ of $U$ UAVs. Since the BS is far from UEs, a team of UAVs which
acted as data collector are deployed. Each UAVs must collect data and reach
it’s final position which are located at $(x_{u,F},y_{u,F})$. For avoiding
collision each UAV fly as the same altitude given by $H_{u}$. Thus, at the
time step $t$, the location of the UAV $u$ is given by
$(x_{u}[t],y_{u}[t],z_{u}[t])$. The location of UE $i$ is given by,
$(x_{i},y_{i},0)$. Therefore the distance between UAV $u$ and user $i$, at
time step $t$ is given by

$d_{i,u}[t]=\sqrt{(x_{u}[t]-x_{i})^{2}+(y_{u}[t]-y_{i})^{2}+(H_{u}[t])^{2}}.$

We assume that the initial position of UAVs is the Center and Docker Station
(CDS) located at $(x_{u,CDS},y_{u,CDS})$. We also assume that the speed of
UAVs are constant and denoted by $V$. UAVs must find a global and local
strategy to collect the maximum amount of data possible while ensuring the
freshness of the data collected. Let $S_{u}$ is the set of visited location by
the $u$-th UAV, the time needed to fly from one location to another is given
by:

$T^{C}_{u}=\sum_{t=1}^{\left|S_{u}\right|}\frac{\sqrt{(x_{u}[t]-x_{u}[t-1])^{2}+(y_{u}[t]-y_{u}[t-1])^{2}}}{V}$

###  2.1 Communication Model

Each UAV communicates with user via the air-to-ground channel. To model the
air-to-ground channel, we assume a Rician fading distribution
$\widehat{\Xi}_{i,u}[t]$

$\widehat{\Xi}_{i,u}[t]=\left(\sqrt{\frac{\Phi}{\Phi+1}}\bar{\Xi}_{i,u}[t]+\sqrt{\frac{1}{\Phi+1}}\widetilde{\Xi}_{i,u}[t]\right),$

where $\Phi$ is the Rician factor, $\bar{\Xi}_{i,u}[t]$ is the line-of-sight
(LoS) component with $\left|\bar{\Xi}_{i,u}[t]\right|=1$, and
$\widetilde{\Xi}_{i,u}[t]$ the random non-line-of-sight (NLoS) component.

Each user $i$ transmits its updates with a power $P_{i}[t]$ during time
interval $t$. Thus, the received power at UAV $u$ is
$\left|\widehat{\Xi}_{i,u}[t]\right|^{2}\beta_{0}\left(d_{i,u}[t]\right)^{-2}P_{i}[t]$,
where $\beta_{0}$ is the average channel power gain at a reference distance
$d_{0}=1\mathrm{~{}m}$.

The signal-to-noise ratio (SNR) of IoT device $i$ with respect to UAV $u$ is
given by

$\Theta_{i,u}[t]=P_{i}[t]\left|\widehat{\Xi}_{i,u}[t]\right|^{2}\beta_{0}\left(d_{i,u}[t]\right)^{-2}/\sigma^{2},$

where $\sigma^{2}$ is the thermal noise power. Therefore, the rate of IoT
device $i$ when it is associated with UAV $u$ during time slot $k$ can be
written

$\overline{R}_{i,u}[t]=B_{i,u}[t]\log_{2}\left(1+\Theta_{i,u}[t]\right),$

where $B_{i,u}[t]$ is the allocated bandwidth between device $i$ and UAV $u$
during time slot $t$.

###  2.2 Age-of-Updates Metric

The objective of this work is to maximize the total data collected, while
maximizing the freshness of the data and ensure that all UAVs complete their
mission. In particular, our aim is to minimize the AoU.

Let $a_{i,u}[t]$ be a binary variable equal to $1$ where the UAV $u$ is
associate with UE $i$ at time slot $t$.

To model the updates’ freshness, we use the same definition as in
[yang2020age]. Particularly, if device $i$ is associated with UAV $u$ during
time interval $t$, its AoU evolves as follows,

$T_{i,u}[t]=\left(T_{i,u}[t-1]+1\right)\left(1-a_{i,u}[t-1]\right),$

where $T_{i,u}[0]=0$. $T_{i,u}[t]$ is the age of the updates received by drone
$u$ and collected from device $i$ during time interval $t$. Specifically, when
the updates of device $i$ are not collected during time interval $t$ (i.e.,
$a_{i,u}[k]=0$), the AoU is increased by one unit of time. Inversely, when the
updates are transmitted, the AoU is reinitialized to zero.

##  3 Problem Formulation

In this work, our aim is to optimize the joint trajectory optimization, data
collection and AoU. In this kind of problem, the objective must be designed to
be dependent on both the total sum-rate of ground UE associated with UAVs and
the AoU of collected data.

Thus the problem can be seen as a multi objective optimization problem. The
objective of each UAVs at each time slot is to maximize the following:

$\mathfrak{r}_{u}[t]=\sum_{i}^{I}\lambda
a_{i,u}[t]\overline{R}_{i,u}[t]-(1-\lambda)T_{i,u}[t]$

Thus the objective of the UAVs team is given by:

$\sum_{u}^{U}\sum_{t}^{T}r_{u}[t]$ (1)

The optimization problem formulated as game is as follow:

$\displaystyle\begin{split}\begin{aligned}
\max_{\mathbf{a},\mathbf{x},\mathbf{y}}\sum_{t}^{T}\sum_{u}^{U}\sum_{i}^{I}\Bigg{(}\lambda
a_{i,u}[t]\overline{R}_{i,u}[t]+(1-\lambda)T_{i,u}[t]\Bigg{)}\end{aligned}\end{split}$
(2a) $\displaystyle\begin{split}\text{s . t.}\quad R_{i,u}[t]>=R_{min},\forall
i\in\mathcal{I},u\in U,\forall t\in T\end{split}$ (2b)
$\displaystyle\begin{split}\sum_{u\in U}a_{i,u}[t]\leq 1,\;\forall
i\in\mathcal{I},\forall k\in\mathcal{K}\end{split}$ (2c)
$\displaystyle\begin{split}a_{i,u}[t]\in{0,1},\forall i\in\mathcal{I},u\in
U,\forall t\in T\end{split}$ (2d) $\displaystyle\begin{split}T_{u}^{c}\leq
T,\forall u\in U\end{split}$ (2e)
$\displaystyle\begin{split}(x_{u}[T],y_{u}[T])=(x_{u,F},y_{u,F}),\forall u\in
U\end{split}$ (2f) $\displaystyle\begin{split}x_{u\text{min}}\leq x_{u}[t]\leq
x_{u\text{max}},\forall u\in U,\forall t\in T\end{split}$ (2g)
$\displaystyle\begin{split}y_{u\text{min}}\leq y_{u}[t]\leq
y_{u\text{max}},\forall u\in U,\forall t\in T\end{split}$ (2h)

where $\mathbf{a}$ are the UAVs-UE association vectors. In this formulation,
the constraint (2b) ensures that the uplink rate is above a predefined
threshold $R_{min}$, so that data is transmitted on time. The constraint (2c)
ensures that no more than one drone serves the same IoT device at the same
time. The constraint (2d) ensures that each UAV complete their mission within
a predefined time. Moreover each UAV need to reach their final position at the
end of the predefined time. Finally, the constraints (2g) and (2h) ensure to
limit the movement of the UAVs.

##  4 Proposed Approach

###  4.1 Muti-Agent Proximal Policy Optimization Algorithm for trajectory and
AoU optimization (MAPPO-TO-AoU)

In this section we describe the proposed algorithm to solve the multi-
objective optimisation and the meta procedure used two generalize the
knowledge of UAVs for quickly adapt in unseen environment.

The reward gleaned by each agent in a multi-agent system is not only related
to its own actions but also to those of other agents. The selection of the
optimal policy of the other agents will be affected by the change of the
policy of an agent and the estimation of the value function will be imprecise,
so it will be difficult to ensure the convergence of the algorithm. In order
to solve this problem, MAPPO[yu2021surprising] uses the idea of CTDE as shown
in figure 1 and extend the Proximal Policy Optimization Algorithm
(PPO)[schulman2017proximal] for the multiagent system. The Critic network
learns a central value function. It can observe global information, including
information from other agents and environmental information. In the case of
the Actor network, each agent calculates the policy only from its own local
observations.

Figure 1: Example of CTDE

The proposed MAPPO-TO-AoU algorithms consist of:

  * •

Agents: $\mathcal{U}$ is the set of agent represented by the UAVs swarm. Each
UAV in the swarm is an agent.

  * •

States: All UAVs have a state constituted by their location. Thus the set of
state of all UAVs is $\mathcal{S}=\\{\boldsymbol{S}[t],t=1\ldots T\\}$. The
state of the UAV $u$ at the time step $t$ is
$s_{u}[t]=\\{x_{u}[t],y_{u}[t]\\}$ and the global state of UAVs at timestep
$t$ is $\boldsymbol{S}[t]=[{s}_{u}[t],u=1,\ldots,U]$.

  * •

Actions: The action of each UAV is the UE that it seeks to collect data or
return to CDS after completing their mission. Each UAVs can move along $x$ or
$y$ axis. At time step $t$ the action taken by UAV $u$ is given by
$a_{u}[t]\in\left\\{0,1,2,3,4\right\\}$, where action $0$ is to stay as its
position and action $1,2,3,4$ design respectively going to UP, DOWN, RIGHT,
and LEFT. The joint action take by UAVs at timestep $t$ is,
$\boldsymbol{a}[t]=[a_{u}[t],u=1,\ldots,U]$

  * •

Rewards: The team global reward $r(\boldsymbol{a},\mathcal{S})$ represent the
total AoU of all UAVs as defined in . Hence the reward of each

  * •

Policy function: We define a policy function
$\boldsymbol{\pi}_{\boldsymbol{\theta}}(a_{u},s_{u})$ represented by the actor
network and parameterised by $\boldsymbol{\theta}$, it is used to generate the
UAV strategy. The agent share the same network but each UAV $u$ takes action
based on its own local state $s_{u}$, and optimize its accumulated discounted
accumulated reward
$J(\theta)=\mathbb{E}_{a[t],s[t]}\left[\sum_{t}\gamma^{t}R\left(s[t],a[t]\right)\right]$.

  * •

Value function: The centralized value function
$V_{\boldsymbol{\phi}}(\boldsymbol{S})$ is represented by the critic network
and parameterised by $\boldsymbol{\phi}$. These network is used to estimate
the achievable future reward by UAVs at every state $\boldsymbol{S}$. The goal
of the critic network is to find the global optimal strategy that maximizes
the expected reward.

UAVs interact with their environment by selecting actions based on the
selected strategy generated by $\boldsymbol{\pi}_{\boldsymbol{\theta}}$. In
each step the UAVs select actions and receive the team global reward and store
its experiences.

The implementation of MAPPO-TO-AoU closely resembles the structure of PPO in
single-agent settings by learning a policy
$\boldsymbol{\pi}_{\boldsymbol{\theta}}$ and a value function
$V_{\boldsymbol{\phi}}(\boldsymbol{S})$. These two functions are represented
as two separate networks, respectively an actor network and a critic network.

$V_{\boldsymbol{\phi}}(\boldsymbol{S})$ is used for variance reduction and is
only utilized during training; hence, it can take as input extra global
information not present in the agent’s local state, allowing PPO in multi-
agent domains to follow the CTDE structure. If these networks are homogenous,
they can be shared by all agents, but each agent can also have a unique
collection of actor and critic networks. For the purpose of notational
simplicity, we assume that all UAVs share the critic and actor networks.

In, MAPPO-TO-AoU algorithm, The actor network is trained to maximize

$\mathcal{L}(\boldsymbol{\theta})=\frac{1}{UB}\sum_{b=1}^{B}\sum_{u=1}^{U}\left[\min\left(r_{\boldsymbol{\theta},u}^{b}A_{u}^{b},\operatorname{clip}\left(r_{\boldsymbol{\theta},u}^{b},1-\epsilon,1+\epsilon\right)A_{u}^{b}\right)+\sigma
E\left[\pi_{\theta}\left(s_{u}^{b}\right)\right)\right]$

where
$r_{\boldsymbol{\theta},u}^{b}=\frac{\pi_{\boldsymbol{\theta}}\left(a_{u}^{b}\mid
s_{u}^{b}\right)}{\pi_{\boldsymbol{\theta}_{old}}\left(a_{u}^{b}\mid
s_{u}^{b}\right)}$. The advantage $A_{u}^{b}$ is computed using the GAE
method, $E$ is the policy entropy, and $\sigma$ is the entropy coefficient
hyperparameter.

The critic network is trained to minimize the loss function

$\mathcal{L}(\boldsymbol{\phi})=\frac{1}{BU}\sum_{b=1}^{B}\sum_{u=1}^{U}\max\left[\left(V_{\boldsymbol{\phi}}\left(s_{u}^{b}\right)-\hat{r}^{b}\right)^{2},\left(\operatorname{clip}\left(V_{\boldsymbol{\phi}}\left(s_{u}^{b}\right),V_{\boldsymbol{\phi}_{\text{old
}}}\left(s_{u}^{b}\right)-\varepsilon,V_{\boldsymbol{\phi}_{\text{old
}}}\left(s_{u}^{b}\right)+\varepsilon\right)-\hat{r}^{b}\right)^{2}\right]$

where $\hat{r}^{b}$ is the discounted reward-to-go.

The proposed MAPPO-TO-AoU solution is summarized in algorithm 2.

Algorithm 1 MAPPO-TO-AoU

1:Input:  User positions, $T$ and $R_{\text{min}}$.

2:Init:  Initialize actor network parameter $\boldsymbol{\phi}$, and initial
critic network parameter $\boldsymbol{\theta}$.

3:for training epoch e=1:E do

4: Initialize data Buffer D

5: for b=1:B do

6: for t=1:T do

7: UAVs Carry experience by executing joint action $\boldsymbol{a}[t]$ and
observe reward $r[t]$ and the next state $\boldsymbol{s}[t+1]$ using initial
policy $\boldsymbol{\pi}_{\boldsymbol{\theta}}$

8: Store experience $\boldsymbol{e}[t]$ in $\tau$

9: Compute Advantage using GAE on experience in $\tau$

10: Compute Discounted reward on experience in $\tau$

11: Store $\tau$ in D

12: Update $\boldsymbol{\theta}$ with data from random mini-batch sample on D

13: Update $\boldsymbol{\phi}$ with data from random mini-batch sample on D

###  4.2 Meta training

In this section we present our proposed meta training approach that meta train
our MAPPO-TO-AoU solution using model agnostic meta-learning
(MAML)[finn2017model]. The meta training procedure will allow the agent to
reach quickly the optimal team strategy in various environment. During the
meta training an environment is created using sample $\boldsymbol{p}_{m}$
drawn from $\mathcal{P}=[\boldsymbol{T},\boldsymbol{R_{min}}]$. Then UAVs
collect experiences $\boldsymbol{e}_{m}$ using their initial policy
$\boldsymbol{\pi}_{\boldsymbol{\theta},m}$ and update their policy and value
function parameterized respectively by $\boldsymbol{\theta}^{\prime}_{m}$ and
$\boldsymbol{\phi}^{\prime}_{m}$ using algorithm 2. Then in a new experiences
$\boldsymbol{e}^{\prime}_{m}$ collected using
$\boldsymbol{\pi}_{\boldsymbol{\theta}^{\prime},m}$ MAPPO-TO-AoU test the
updated policy and value function at current environment sampling
$\boldsymbol{p}_{c}$. The values of the loss functions that measure the
distance between the policy and the updated value functions that produce
optimal team strategies, relative to the policy and the optimal value
functions, constitute the feedback from these tests to the UAVs. These updates
are done as follows:

$\mathcal{L}(\boldsymbol{\theta}^{\prime}_{m})=\frac{1}{UB}\sum_{b=1}^{B}\sum_{u=1}^{U}\left[\min\left(r_{\boldsymbol{\theta}^{\prime}_{m},u}^{b}A_{u}^{b},\operatorname{clip}\left(r_{\boldsymbol{\theta}^{\prime}_{m},u}^{b},1-\epsilon,1+\epsilon\right)A_{u}^{b}\right)+\sigma
E\left[\pi_{\boldsymbol{\theta}^{\prime}_{m}}\left(s_{u}^{b}\right)\right)\right]$
(3)

$\mathcal{L}(\boldsymbol{\phi}^{\prime}_{m})=\frac{1}{BU}\sum_{b=1}^{B}\sum_{u=1}^{U}\max\left[\left(V_{\boldsymbol{\phi}^{\prime}_{m}}\left(s_{u}^{b}\right)-\hat{r}^{b}\right)^{2},\left(\operatorname{clip}\left(V_{\boldsymbol{\phi}^{\prime}_{m}}\left(s_{u}^{b}\right),V_{\boldsymbol{\phi}^{\prime}_{m,\text{old
}}}\left(s_{u}^{b}\right)-\varepsilon,V_{\boldsymbol{\phi}^{\prime}_{m,\text{old
}}}\left(s_{u}^{b}\right)+\varepsilon\right)-\hat{r}^{b}\right)^{2}\right]$
(4)

The initial critic and actor parameter are updated as follow:

$\mathcal{L}(\boldsymbol{\phi})=\sum_{\boldsymbol{p}_{m}\sim\mathcal{P}}\mathcal{L}({\boldsymbol{\phi}^{\prime}_{m}})$
(5)

$\mathcal{L}(\boldsymbol{\theta})=\sum_{\boldsymbol{p}_{m}\sim\mathcal{P}}\mathcal{L}^{\mathrm{CLIP}}(\boldsymbol{\theta}^{\prime}_{m})$
(6)

The proposed Meta-MAPPO-TO-AoU solution is summarized in algorithm

Algorithm 2 META-MAPPO-TO-AoU

1:Input:  User positions, $T$ and $R_{\text{min}}$.

2:Init:  Initialize actor network parameter $\boldsymbol{\phi}$, and initial
critic network parameter $\boldsymbol{\theta}$.

3:for training epoch e=1:E do

4: for environment sample m=1:M do

5: sample $\boldsymbol{p}_{m}\sim\mathcal{P}$.

6: Initialize data Buffer D.

7: for b=1:B do

8: for t=1:T do

9: UAVs carry experience using initial policy
$\boldsymbol{\pi}_{\boldsymbol{\theta}}$.

10: Store experience $\boldsymbol{e}_{m}[t]$ in $\tau$.

11: Compute Advantage using GAE on experience in $\tau$.

12: Compute reward-to-go on experience in $\tau$.

13: Store $\tau$ in D.

14: Update value function using (4) with data from random mini-batch sample on
D.

15: Update policy function using (3) with data from random mini-batch sample
on D.

16: for t=1:T do

17: UAVs carry experience using initial policy
$\boldsymbol{\pi}^{\prime}_{\boldsymbol{\theta}_{m}}$.

18: Store experience $\boldsymbol{e}^{\prime}_{m}[t]$ in $\tau$.

19: Store $\tau$ in D.

20: Compute loss using $\mathcal{L}({\boldsymbol{\phi}^{\prime}_{m}})$ and
$\mathcal{L}^{\mathrm{CLIP}}(\boldsymbol{\theta}^{\prime}_{m})$ using (3) and
(4) with data from random mini-batch sample on D.

21: Update $\boldsymbol{\theta}$ with data from random mini-batch sample on D
using (6).

22: Update $\boldsymbol{\phi}$ with data from random mini-batch sample on D
using (5).

Figure 2: Mean reward over iterations

Figure 3: Mean data collected over iterations

Figure 4: Mean AoU over iterations

Figure 5: $\lambda_{1}=0.1$, $\lambda_{2}=0.4$ and $\lambda_{3}=0.5$

Figure 6: Mean reward over iterations

Figure 7: Mean data collected over iterations

Figure 8: Mean AoU over iterations

Figure 9: $\lambda_{1}=0.1$, $\lambda_{2}=0.6$ and $\lambda_{3}=0.3$

##  5 Simulation Results

$I$ | $20$ | $U$ | $3$
---|---|---|---
$R_{min}$ | $100000$ | $\lambda_{3}$ | $0.4$
$x$ | $5$ | $\beta_{0}$ | $1$
$y$ | $5$ | $B_{i,u}$ | $[1500,1700]$
$H$ | $[80,100]$ | $\alpha$ | $2$
$\lambda_{1}$ | $0.3$ | $\sigma^{2}$ | $1e^{-15}$
$\lambda_{2}$ | $0.3$ | $T$ | $100$
$P_{i}$ | $[0,1]$ | $\widehat{\Xi}_{i,u}$ | $[0,1]$
Table 1: Experiment setup

In order to evaluate the performance of our algorithm, we perform a simulation
with the parameters described in the table 1. We performed simulations on two
sets of lambda as shown in figure 5 and 9, where we plot the reward gleaned by
the UAVs over iterations, as well as the total data collected and the AoU of
the system. The results show that, in both situations, our algorithm converges
faster and achieves higher rewards (figure 5 and 9), also the results show
that the UAVs can find quickly a global strategy that minimises the AoU
(figure 5 and 9) and maximises the total data collected by the UAVs (figure 5
and 9), compared to the algorithms used as a benchmark. This confirms the fact
that on-policy methods can perform very well when solving tasks where
cooperation is required. Furthermore, our meta-training procedure shows that
the agents adapt very quickly to a new environment and manage to collect even
more data while minimising the AoU of the system. With these two sets of
lambda, we notice that the AoU and the total data collected sometimes follow
the same evolution. The AoU formulation allows us to say that the AoU is
minimal when the number of users served is maximal and furthermore it can be
said that serving more users is synonymous with collecting more data. But as
in our case the IoT devices generate different data over time, the total data
collected can increase or decrease without the AoU increasing or decreasing,
we can see this in figures 5 and 5 where the AoU remains constant over time,
but the total data collected fluctuates. This shows that drones can find an
overall strategy to serve all users of the system but are not always sure to
collect the same amount of data.

For meta learning, we sample $100$ tasks for
$\boldsymbol{R_{min}}=[100000,120000]$ and $\boldsymbol{T}=[100,200]$. The
results presented in figures LABEL:fig:reward_meta,
LABEL:fig:data_collected_meta and LABEL:fig:AoU_meta show that our meta-
training procedure converges quickly and achieves the best rewards among all
algorithms.

##  6 Conclusion

In this work, we develop a practical framework in which drones can act
cooperatively to achieve multiple objectives. We formulate our problem into a
multi-objective optimization problem that is shown to be non-convex. We then
propose a solution based on MARML to efficiently solve our formulated problem.
In order to enable our UAVs to perform well in unknown environments, we
propose a meta-procedure that allows meta-training of UAVs. Simulation results
show that our proposed solution converges quickly and solves the task better
than the reference solution. As we used a solution based on the actor-critic
framework, the extension of this work is to consider the optimisation of the
3D position of the UAVs rather than the 2D position.