# RAPF: Efficient path planning for lunar microrovers

Thomas Manteaux  Signal Processing Systems, Faculty of EEMCS, Delft University
of Technology, The Netherlands  Advanced Quantum Architecture Lab (AQUA),
Ecole Polytechnique Fédérale de Lausanne (EPFL), Switzerland  David Rodríguez-
Martínez  Advanced Quantum Architecture Lab (AQUA), Ecole Polytechnique
Fédérale de Lausanne (EPFL), Switzerland  Raj Thilak Rajan  Signal Processing
Systems, Faculty of EEMCS, Delft University of Technology, The Netherlands

###### Abstract

Efficient path planning is key for safe autonomous navigation over complex and
unknown terrains. Lunar Zebro (LZ), a project of the Delft University of
Technology, aims to deploy a compact rover, no larger than an A4 sheet of
paper and weighing not more than $3$ kilograms. In this work, we introduce a
Robust Artificial Potential Field (RAPF) algorithm, a new path-planning
algorithm for reliable local navigation solution for lunar microrovers. RAPF
leverages and improves state of the art Artificial Potential Field (APF)-based
methods by incorporating the position of the robot in the generation of
bacteria points and considering local minima as regions to avoid. We perform
both simulations and on field experiments to validate the performance of RAPF,
which outperforms state-of-the-art APF-based algorithms by over 15% in
reachability within a similar or shorter planning time. The improvements
resulted in a 200% higher success rate and 50% lower computing time compared
to the conventional APF algorithm. Near-optimal paths are computed in real-
time with limited available processing power. The bacterial approach of the
RAPF algorithm proves faster to execute and smaller to store than path
planning algorithms used in existing planetary rovers, showcasing its
potential for reliable lunar exploration with computationally constrained and
energy constrained robotic systems.

###### Index Terms:

Path planning, Artificial Potential Field, lunar rover, microrovers

##  I Introduction

Despite significant advancements in space technology and robotics, navigating
efficiently across unknown and unstructured environments remains an ongoing
pursuit. This paper explores the importance of efficient path planning for
small lunar rovers ($<$100 kg) and microrovers ($<$10 kg), which have gained
significance with the advent of cheaper and commercial exploration missions
over the last decade. Several state-of-the-art lunar rover path planning
solutions have been developed, drawing from diverse fields such as robotics,
artificial intelligence, and space exploration. For instance, Yu et al. [1]
propose a heuristic-based approach that considers terrain complexity and
energy efficiency for global path planning on the lunar surface. Another
notable work by Peng et al. [2] utilizes deep neural networks to to perform
data-driven path planning on the basis of heuristic search. These advancements
showcase the ongoing efforts to enhance the autonomy and efficiency of lunar
exploration missions through innovative path planning techniques. Although
recent trends, such as deep learning, show promise, they remain
computationally intensive, as compared to more resource efficient algorithms
e.g., Rule-based approaches, including Dijkstra-based [3], Artificial
Potential Field (APF)-based [13], and sampling-based [5] algorithms. Small
lunar rovers, operating in complex environments with communication delays,
restricted computational, memory and power capabilities, and limited sensor
suites demand robust, real-time path planning algorithms. These algorithms
must cope with the irregularities of the lunar surface—characterized by uneven
and unconsolidated terrains filled with craters and rocks of different
sizes—while addressing measurement noise from sensors on-board the rover. Due
to their small size and reduced mobility, microrovers face challenges in
negotiating rough or uneven terrain.

Figure 1: a) Field test layout. b) Potential map computed by the RAPF
algorithm. Dark blue represents low-potential regions while red represents
high-potential regions. The orange line shows the computed path, black line
defines the ground truth, and the yellow cross shows where a local minimum was
found. c) Lunar Zebro prototype used for testing

Small size constrains sensor range, limiting the perception of surrounding
environments for micro-rovers preventing comprehensive perception and posing
challenges for effective path planning in navigation. The power constraint of
lunar rovers demands low-complexity and efficient path planning, striking a
balance between computational cost and performance. Path planning emerges,
therefore, as a critical aspect impacting mission efficiency and success.

###  I-A Contributions

To robustly and efficiently compute the path to be followed by a lunar
microrover, we propose a new lightweight solution, RAPF, based on adapting and
expanding an existing APF-based algorithm. Our approach is adapted to run
onboard the LZ walking rover on a Raspberry Pi 4B. We demonstrate the
performance of our algorithm both in simulations and during field tests across
a lunar-like environment. Results show the greater accuracy and lower memory
allocation requirements of our approach, displaying comparable performances
even when compared to optimal terrestrial algorithms such as A*.

###  I-B Related work

Different approaches exist for computing the optimal path to follow between
points of interest.

Combinatorial methods—which include the Dijkstra, A* [3], and D* [4]
algorithms—explore discrete paths through a graph representation of the
environment. Graph nodes represent specific configurations of the robot. Edges
represent feasible transitions between configurations. Only transitions
between neighboring graph nodes are allowed. Graph search algorithms explore
the graph by visiting nodes and searching for the optimal path based on
certain criteria such as distance or cost. Despite ensuring optimal paths and
completeness, combinatorial methods demand significant memory and
computational capabilities, struggling in environments where the map is
frequently recomputed.

Sampling-based algorithms divide a continuous space into sample points forming
a roadmap or tree, with Rapidly Exploring Random Tree (RRT) [5] and
Probabilistic Roadmaps (PRM) [6] as prominent examples of this approach.
Sampling involves a deterministic or stochastic selection of points within the
configuration space, representing potential robot configurations. These
sampled points stand as potential configurations that the robot can take. The
objective is to comprehensively cover the configuration space to unveil
feasible paths. The connection of these points forms a graph structure,
balancing optimality and completeness. However, these algorithms can be
sensitive to the initial configuration [7] and may necessitate time-intensive
pre-processing.

Biologically inspired methods such as spiking neural models [8] use artificial
neural networks to mimic the brain’s structure, enabling robots to adaptively
plan paths based on sensory inputs. Swarm intelligence, inspired by collective
behavior, uses algorithms like Particle Swarm Optimization [9] and Ant Colony
Optimization [10]. Genetic algorithms [11], derived from natural selection
principles [12], iteratively refine potential solutions. As they replicate
problem-solving strategies found in nature, biologically inspired methods
offer a large adaptability making them well-suited in complex and dynamic
environments. They provide robustness against unpredicted obstacles and
changing surroundings. However, they require a large number of parameters to
be fine-tuned and their implementation is usually complex. They therefore have
a large computational cost. Biologically inspired methods do not always
outperform traditional path planning methods – combinatorial, sampling-based
and APF-based methods – in static environments and their complexity might
hinder their implementation on embedded systems working in real-time.

Lastly, potential field techniques [13] guide robots using artificial forces.
In the Artificial Potential Field (APF) method, an attractive potential pulls
the robot towards the target, while repulsive potentials keep it away from
obstacles. The APF approach ensures collision-free paths in real-time with low
computational complexity. However, they face local minima issues, potentially
trapping the robot in sub-optimal positions. Local minimas occur due to
various reasons such as symmetry between the different obstacles, closely
spaced obstacles, or in highly cluttered environments. Tuning the balance
between repulsive and attractive potentials is crucial for successful path
planning and to avoid local minima. The next section presents the original APF
algorithm and some of its more promising variations. These form the basis of
the approach proposed herein.

##  II Artificial Potential Field (APF) algortihms

The classical APF was originally proposed in [13], in which a robot navigates
in an artifical force field, with attractive forces being emitted by the
target and repulsive forces by the obstacles. These forces are directly
derived from quadratic potential functions, which are parameterized by the
robots distance to a target or an obstacle. The robot, therefore, aims at
minimizing its potential cost function at any given point along its’ path. A
local minima event occurs, when no adjacent point exists in the neighborhood
of the rover with a lower potential as compared to the current position of the
robot. The attractive potential acting on a robot at the global position
$\textbf{x}=[x,y]^{T}$ in a 2-dimensional scenario, from a target located at
$\textbf{x}_{t}=[x_{t},y_{t}]^{T}$ is defined as

$J_{a}(\textbf{x})=\frac{1}{2}k_{a}d_{t}^{2}(\textbf{x}),$ (1)

where $k_{a}$ is the attractive gain,
$d_{t}(\textbf{x})=||\textbf{x}-\textbf{x}_{t}||_{2}$ is the Euclidean
distance between the robot and the target. The repulsive potential acting on a
robot located at x from the nth detected obstacle at location
$\textbf{x}_{o}^{n}=[x_{o}^{n},y_{o}^{n}]^{T}$ is defined as

$J_{o}^{n}(\textbf{x})=\begin{cases}\frac{1}{2}k_{r}(\frac{1}{d_{n}(\textbf{x})}-\frac{1}{\rho_{0}})^{2}&\text{if
$d_{n}(\textbf{x})\leq\rho_{0}$},\\\ 0&\text{otherwise},\end{cases}$ (2)

where $k_{r}$ is the repulsive gain, $\rho_{0}$ the maximum distance within
which the potential field of the obstacle exerts its influence, and
$d_{n}(\textbf{x})=||\textbf{x}-\textbf{x}_{o}^{n}||_{2}$ is the Euclidean
distance between the robot and the nth detected obstacle. Let $N_{obs}$ be the
number of obstacles detected by the robot at position x, then the total
repulsive potential acting on the robot is

$J_{r}(\textbf{x})=\sum_{n=1}^{N_{obs}}J_{o}^{n}(\textbf{x})$ (3)

and the total potential function is the sum of the attractive and the
repulsive potentials as stated below

$J(\textbf{x})=J_{a}(\textbf{x})+J_{r}(\textbf{x}).$ (4)

The robot’s direction of motion is determined by calculating the gradient
relative to the potential of its current position.

TABLE I: Overview of state-of-the-art APF-based algorithms Algorithm |  APF |  RVF |  CRBAPF* |  RAPF
---|---|---|---|---
Key concept |  Attractive/Repulsive potential |  Rotated vector field |  Bacteria points + RWT |  Bacteria points + artificial obstacles
Cost function type |  Quadratic |  Quadratic |  Gaussian |  Gaussian
Handling Local minima |  - |  Rotation of the field |  Random walk technique |  Artificial obstacles
Trajectory smoothness |  Abrupt changes |  Smooth trajectories |  Oscillations for low $N_{B}$ |  Smooth trajectories

###  II-A RVF

An alternative potential field method for escaping local minima is proposed in
[14] and [15], called the Rotated Vector Field (RVF), which reorients the
force field directly derived from the classical APF potential functions as
described below. The resulting attractive force is defined as the gradient
function of the attractive potential

$\textbf{F}_{a}(\textbf{x})=-\nabla J_{a}(\textbf{x}),$ (5)

where $J_{a}(\textbf{x})$ is defined in (1). Similar to the attractive force,
the resulting repulsive force is directly derived from the potential function
as

$\textbf{F}_{r}^{n}(\textbf{x})=-\nabla J_{o}^{n}(\textbf{x}),$ (6)

where $J_{o}^{n}(\textbf{x})$ is defined in (2). Similar to the total force
field defined for APF in (4), the total force is the sum of the attractive and
repulsive forces written as

$\textbf{F}(\textbf{x})=-\nabla J(\textbf{x})=-(\nabla
J_{a}(\textbf{x})+\nabla J_{r}(\textbf{x})).$ (7)

The main benefit of this approach is that instead of being repelled when
approaching an obstacle, the robot is forced to move around it. This is
achieved by rotating the vector field around the obstacles in one direction or
another. This creates a vortex field at the location of the nth detected
obstacle given by

$F_{n}(\textbf{x})=\pm\begin{pmatrix}\frac{\partial
J_{t,n}(\textbf{x})}{\partial y}\\\ \frac{-\partial
J_{t,n}(\textbf{x})}{\partial x}\\\ \end{pmatrix}.$ (8)

The intensity of the vector field remains the same, only the direction
changes. The vortex field is bounded at a distance $\rho_{0}$ from the
obstacles, the maximum distance within which the potential field of the
obstacle exerts its influence, such that the entire vector field does not take
on this rotary shape.

###  II-B CRBAPF*

In [16] and [17] a set of enhanced APF algorithms called Bacteria APFs (BAPFs)
are proposed. These algorithms involve what is called bacteria points, i.e.,
virtual points around the robot representing a potential solution to the
optimization problem of the APF algorithm. A new APF algorithm based on the
BAPF concept was introduced in [18], named Changing Radii Bacteria APF
(CRBAPF*), which actively modifies the potential functions and enables a
random walk technique (RWT) to escape local minima. The potential functions
are computed for the bacteria points by taking into account their distances
from the target and any obstacles detected in the environment. The bacteria
with the lowest cost is chosen as the next position for the robot, and the
latter moves towards that selected point. The cost functions from the target
and the obstacles are expressed as Gaussian functions instead of quadratic
functions in the classical APF method. Using Gaussian over quadratic cost
functions presents several advantages such as a more gradual increase in
repulsion as the robot approaches an obstacle and a smoother transition from
attractive to repulsive forces. This solution results in smoother
trajectories, a more natural behavior, and a better handling of narrow
passages. The attractive potential acting on a robot located at
$\textbf{x}=[x,y]^{T}$ from a target located at
$\textbf{x}_{t}=[x_{t},y_{t}]^{T}$ is defined as

$J_{a}(\textbf{x})=-\alpha_{a}e^{-\mu_{a}d_{t}(\textbf{x})},$ (9)

with $\alpha_{a}$ and $\mu_{a}$ constant parameters established through
empirical means and they can be seen as the height and the width of the
potential source. $d_{t}(\textbf{x})$ is the squared Euclidean distance
between the robot and the target. The repulsive potential acting on a robot
located at x from the nth detected obstacle at location
$\textbf{x}_{o}^{n}=[x_{o}^{n},y_{o}^{n}]^{T}$ is defined as

$J_{o}^{n}(\textbf{x})=\begin{cases}0&\text{if
$\sqrt{d_{n}(\textbf{x})}>\rho_{u}$},\\\
\alpha_{o}e^{-\mu_{o}d_{n}(\textbf{x})}&\text{if
$\rho_{l}<\sqrt{d_{n}(\textbf{x})}<\rho_{u}$},\\\ \infty&\text{if
$\sqrt{d_{n}(\textbf{x})}<\rho_{l}$},\end{cases}$ (10)

with $\rho_{l}$ the lower radius and $\rho_{u}$ the upper radius setting the
boundaries of the potential source. $\alpha_{o}$ and $\mu_{o}$ are constants
tuned by testing and dependent on the average density of obstacles, and
$d_{n}(\textbf{x})$ is the squared Euclidean distance between the robot and
the nth detected obstacle. For a number $N_{obs}$ of obstacles detected by the
robot at a position x, the total repulsive potential acting on the robot is

$J_{r}(\textbf{x})=\sum_{n=1}^{N_{obs}}J_{o}^{n}(\textbf{x}),$ (11)

and the total potential function is the sum of the attractive and the
repulsive potentials as stated below.

$J_{r}(\textbf{x})=J_{a}(\textbf{x})+J_{r}(\textbf{x}).$ (12)

Bacteria potential

For a robot located at $\textbf{x}=[x,y]^{T}$, the position of the jth the
bacteria is $\textbf{x}_{b,j}=[x_{b,j},y_{b,j}]^{T}$, defined as

$\begin{pmatrix}x_{b,j}\\\ y_{b,j}\\\
\end{pmatrix}=\begin{cases}x+\rho_{b}cos(\frac{2\pi n_{j}}{N_{B}}),\\\
y+\rho_{b}sin(\frac{2\pi n_{j}}{N_{B}}),\end{cases}$ (13)

where $N_{B}$ is the number of bacteria, $n_{j}=1,2,...,N_{B}$ and $\rho_{b}$
the radius around the robot where bacteria are generated. Similarly to the
rover, the total potential function at the location of the jth bacteria is
expressed as

$J_{b,j}(\textbf{x}_{b,j})=J_{a}(\textbf{x}_{b,j})+J_{r}(\textbf{x}_{b,j}).$
(14)

As the robot is looking for the lowest potential towards the target, the next
selected bacteria point j shall meet the criterion

$J_{b,j}(\textbf{x}_{b,j})-J_{r}(\textbf{x})<0.$ (15)

If several bacteria points meet (15), the one with the lowest Euclidean
distance toward the target is selected.

Random walk technique
To prevent the robot from getting stuck in a local minimum, a random walk
technique is employed. The robot randomly selects a bacteria point around it,
moves to it, and repeats the process for a pre-defined amount of time. The
random selection of the bacteria point follows a uniform distribution. The
only condition to be met follows from (10).

###  II-C Robust and Lightweight APF

Our proposed approach in this work, the Robust and Lightweight APF (RAPF)
method, differs from the CRBAPF* algorithm in two different ways. First, RAPF
uses a different local minima escape approach instead of the random walk
technique. In RAPF, when local minimas are encountered, then they are labeled
as artificial obstacles. These artificial obstacles generate a repulsive
potential (Eq. 10) that influences the total potential function (Eq. 12). The
existing planned path is cleared and a new path is computed from the current
position, which includes the local minima as an obstacle. Another improvement
introduced into the RAPF algorithm consists of taking into account the
position of the robot with respect to the target during the generation of the
bacteria points. The $N_{B}$ bacteria points are generated such that there is
always a bacteria point along the robot-target vector. This removes the
oscillation effect from the CRBAPF* algorithm when a small number of bacteria
points are used. This feature enables the RAPF algorithm to use a reduced
number of bacteria points ($N_{B}$ = 8), resulting in a more lightweight
method. Algorithm 1 describes the working principle of the RAPF algorithm.

At the onset crucial parameters are set-up like observed obstacles, the number
of bacterial agents, movement step size, and maximum duration. It initializes
key variables such as the current position, a timer for elapsed time tracking,
and a boolean flag for path determination. Within each iteration, the
algorithm assesses the robot’s position by calculating the cost function (J)
considering obstacle influence. It computes the cost function for each
bacteria point ($Idx_{B}$) and selects a new robot position based on these
evaluations, appending it to the path. If a local minimum is detected, the
algorithm resets the path to the starting position and updates observed
obstacles. It also checks if the maximum time has been exceeded or if the
robot is within the target margin. If any condition is met, the algorithm
terminates and returns the path; otherwise, it iterates until finding a
satisfactory path or reaching the time limit.

The key properties of the different APF algorithms including RAPF are listed
in Table I.

Algorithm 1 RAPF algorithm

1:$\mathbf{x},\mathbf{x}_{t},N_{obs},N_{B},StepSize,MaxTime,\\{\mathbf{x}_{i},r_{i}\\}_{i\in
N_{obs}},$ $GoalMargin$

2:$PathToTarget$

3:$\mathcal{N}_{obs}=\\{1,2,\ldots,N_{obs}\\}$

4:$\mathcal{N}_{B}=\\{1,2,\ldots,N_{B}\\}$

5:$Path=\mathbf{x}$

6:$Timer=CurrentTime$

7:$FindingPath=True$

8:while $FindingPath$ do

9: $J_{a,robot}\leftarrow get\\_Ja(Path,\mathbf{x}_{t})$

10: $J_{r,robot}\leftarrow get\\_Jr(Path,\\{\mathbf{x}_{i},r_{i}\\}_{i\in
N_{obs}})$

11: $J_{robot}\leftarrow J_{a,robot}+J_{r,robot}$

12: for $Idx_{b}\in\mathcal{N}_{B}$ do

13: $\mathbf{x_{b}}\leftarrow get\\_bacteria(Path,StepSize,N_{B},Idx_{b})$

14: $J_{a,bacteria}\leftarrow get\\_Ja(\mathbf{x_{b}},\mathbf{x}_{t})$

15: $J_{r,bacteria}\leftarrow
get\\_Jr(\mathbf{x_{b}},\\{\mathbf{x}_{i},r_{i}\\}_{i\in N_{obs})}$

16: $J_{bacteria}\leftarrow J_{a,bacteria}+J_{r,bacteria}$

17: end for

18: $new_{x}\leftarrow
pick\\_bacteria(Path,J_{robot},J_{bacteria},\mathcal{N}_{obs})$

19: $Path\leftarrow Add(new_{x},Path)$

20: if $local\\_minimum(Path)$ then

21: $Path=\mathbf{x}$

22: $\mathcal{N}_{obs}\leftarrow
Add(\mathbf{x}_{loc\\_min},\\{\mathbf{x}_{i},r_{i}\\}_{i\in N_{obs}})$

23: end if

24: if $Timer-CurrentTime>MaxTime$ then

25: $FindingPath=False$

26: end if

27: if $||\mathbf{x}_{t}-\mathbf{x}_{j}||_{2}<GoalMargin$ then

28: return $Path$

29: end if

30:end while

31:return $None$

##  III Simulations

###  III-A Virtual environment

To test the performance of our algorithms in a more realistic, we perform
simulations in a virtual environment consisting of rocks and craters on Lunar
surface. We achieve this by studying the size-frequency and cumulative number
distributions of these geographical features, which are derived from past
lunar mission data. Based on data from [19], the cumulative fractional area
covered by rocks with a diameter greater than D versus rock diameter (D) is
given by

$F_{r}(D)=k_{r}e^{-q_{r}(k_{r})*D},$ (16)

where $k_{r}$ is the rock abundance [%] and $q_{r}$ [$m^{-1}$] is a
coefficient determined empirically. From (16) the cumulative number of rocks
per square meter versus rock diameter ($D_{r}$) is then given by

$N_{r}(D)=\frac{4q_{r}(k_{r})k_{r}}{\pi}\left(\frac{e^{-q_{r}(k_{r})D}}{D}-q_{r}(k_{r})Ei(-q_{r}(k_{r})D)\right).$
(17)

The rock abundance is expected to be smaller than 10% as observed in previous
NASA and ESA missions [20, 21]. Extending values from previous missions
(Apollo, Chang’e, and Luna) to the LZ landing site we obtain $k_{r}$ = 0.02
and $q_{r}$ = 1.6 $m^{-1}$. The height-to-diameter ratio allows us to
determine that rocks with a diameter greater than $D_{r,crit}=6.5$cm are
considered obstacles for LZ, i.e. a height greater than 35 mm. The crater
abundance model uses the same (16) and (17). A Monte-Carlo simulation is run
to create three different scenarios in which to compare the previously
introduced algorithms. The map is 30x30$m^{2}$ with a distribution of the
obstacles within $(x_{min},y_{min})=(5,5)$ and $(x_{max},y_{max})=(25,25)$.
The rocks cover 1.8% and the craters 15% of this 20x20$m^{2}$ area. The
starting position is located at $(x_{s},y_{s})$ = (2,2) and the target
position is a circle centered at $(x_{g},y_{g})$ = (28,28), with a radius of
0.5m. The radius of the rover is 0.20m. A* is used as the optimal reference
for this analysis. Table II presents the distribution of rocks and craters for
each scenario. Despite the large differences between the numbers, the obstacle
coverage of the surface is constant as previously defined ($A_{rock}=1.8\%$,
$A_{crater}=11\%$). 500 simulations are run on each environment for each
algorithm. The position of obstacles is randomly generated for each new
simulation.

TABLE II: Number of obstacles for each scenario Scenario | Rocks | Craters | Total obstacles
---|---|---|---
A | 42 | 38 | 80
B | 88 | 32 | 120
C | 137 | 24 | 161

###  III-B Performance metrics

The following metrics are used to evaluate the performance of the algorithms.

Reachability, $R_{s}$: It defines the ratio of successful navigation trials
$N_{s}$ over $N_{m}$ attempts. A navigation trial where the rover succeeds in
reaching the target without colliding with obstacles defines a successful
navigation trial.

$R_{s}=\frac{N_{s}}{N_{m}},$ (18)

Computational complexity, $T_{c}$ : It defines the sum of the time the
algorithm takes to compute each of the successful paths to the target until
the rover reaches it.

Path length, $M_{s}$: It defines the total distance covered by the rover to
reach its target. The average path length for $N_{s}$ successful navigation
trials is given by

$\overline{M_{s}}=\frac{1}{N_{s}}\sum_{k=1}^{N_{s}}M_{s}(k),$ (19)

with $M_{s}(k)$ the path length of the kth successful run.

Safety, SF: Considering the robot has detected N obstacles during a successful
navigation attempt k of length $M_{s}$ at locations
$r_{1},r_{2},...,r_{M_{s}}$, the safety parameter refers to the average
minimum distance between the robot and the detected obstacles. For $N_{s}$
successful navigation trials it is defined as

$\overline{SF}=\frac{1}{N_{s}}\sum_{k=1}^{N_{s}}\frac{\sum_{m=1}^{N(k)}d_{min}(r_{o}^{m})}{N(k)},$
(20)

with
$d_{min}(r_{o}^{m})=min\\{||r_{1}-r_{o}^{m}||,||r_{2}-r_{o}^{m}||,...,||r_{M_{s}}-r_{o}^{m}||\\}$
and $r_{o}^{m}$ and the location of the mth detected obstacle. Safety will be
used to compare field testing with simulation performance of the algorithms.

TABLE III: Reachability, mean planning time and mean path length for each scenario and each algorithm | Reachability [%] | Mean planning time [ms] | Mean path length [m]
---|---|---|---
Algorithm |  Scenario A |  Scenario B |  Scenario C |  Scenario A |  Scenario B |  Scenario C |  Scenario A |  Scenario B |  Scenario C
A* |  100 |  100 |  100 |  3906.6 |  5444.3 |  7699.1 |  38.6 |  38.5 |  38.7
APF |  62.2 |  51.0 |  30.6 |  1225.8 |  1669.6 |  2180.0 |  39.2 |  40.0 |  40.5
RVF |  73.6 |  63.4 |  46.6 |  1634.8 |  2237.5 |  2966.6 |  41 |  41.7 |  42.2
CRBAPF* |  83 |  83.2 |  80 |  579.6 |  888.1 |  1328.6 |  43.5 |  44.2 |  45.2
RAPF |  96.4 |  93.8 |  91.8 |  589.6 |  889.5 |  1276.5 |  39.8 |  40.3 |  41.1

###  III-C Simulation results

Simulation results are listed in Table III. As anticipated, our reference
algorithm, A*, consistently reaches the target in every scenario. The
classical APF approach exhibits poor results, reaching the target less than
one in three instances as cluttering increases. RVF demonstrates relatively
high reachability in non-cluttered environments but experiences a drastic
decrease with cluttering. In contrast, CRBAPF* maintains relatively high
reachability independent of the environment. Notably, RAPF yields excellent
results, exceeding CRBAPF* reachability by more than 15%. Planning time for
all the algorithms increases with cluttering (from Scenario A to Scenario C),
with A* having the longest computation time, as expected. Both RAPF and
CRBAPF* share the same planning time, outperforming both RVF and A*. On
average, they are 6 times faster than A* and 2.5 times faster than RVF. RVF
has a larger computation time than APF as it computes the field of forces and
the rotation around each obstacles.

Stability and outliers: CRBAPF* have a low variability but a certain amount of
outliers, meaning that some configurations may lead to an abnormal computation
time. RVF is the most stable and consistent algorithm in terms of planning
time as it has a very low variability and almost no outliers. A* presents the
same characteristics as RAPF and CRBAPF* in terms of outliers, but shows much
higher median and variability, making it less consistent than RAPF. Outliers
for CRBAPF* come from the random walk technique implemented to escape the
local minima, which our algorithm solves.

Path length: In terms of path length, the difference between algorithms
becomes less evident. A* and APF yield the shortest paths, but his result may
be less reliable for APF, given that the mean path is computed based on
successful paths. However, APF exhibits an overall low number of completed
paths, indicating that path length is only computed in favorable environments.
RVF produces paths between 5% to 10% longer than A*. This can be attributed to
RVF providing the same direction of rotation to the vector field of clusters
of obstacles, leading to detours, while A* finds a direct path between
obstacles. CRBAPF* yields a path length between 10% and 15% greater than A*,
attributable to the random walk technique implemented in CRBAPF*. RAPF
generates paths between 3% and 6% longer than those computed by A*,
significantly lower than the CRBAPF* paths. The path length for A* remains
independent of the environment, whereas RVF, RAPF, and CRBAPF* exhibit a
slight increase with cluttering. Path length distributions are highly
asymmetrical. The median length is lower for CRBAPF* than RVF, while the
opposite holds for their means. This discrepancy arises from the substantial
and widespread number of outliers generated by CRBAPF*, attributed to the tail
effect stemming from the random walk technique. This technique can produce
extremely long paths as the path generation relies solely on randomness.

Planning time: In terms of planning time, A* presents a large number of
outliers but with low scattering. In contrast, RAPF displays half as many
outliers and a considerably lower Interquartile Range (IQR) than CRBAPF*,
where the random walk technique is not used. The results suggest that RAPF is
more consistent between paths, demonstrating increased stability and
robustness. The RAPF algorithm displays approximately 50% lower navigation
time, up to 200% higher success rate, and similar path lengths compared to the
conventional APF algorithm. RAPF also outperforms RVF in success rate, path
length, and planning time, while the planning time distribution closely
mirrors that of the CRBAPF*. Notably, it also outperforms the CRBAPF* in path
length, producing a near-optimal path.

##  IV Field validation

###  IV-A Testbed description and testing approach

The rover is tested on a 3x4 m sandy terrain featuring rocks, craters, and
slopes mimicking the lunar surface (see Fig. 1). The bottom left corner
defines the origin of the 2D Cartesian reference frame. A Raspberry Pi 4B 2GB
is used as an on-board computer on which the different path planning
algorithms are implemented. Six motors drive the six legs. An IMU and an RF
Tag are carried on board. The RF Tag communicates with 4 anchors located in
each corner of the testbed for localization purposes. A threshold filter is
applied to the raw data from the RF Tag to reduce noise. A Kalman filter
merges the data from the motor encoders, IMU, and filtered RF Tag data at a
frequency of 5Hz. Obstacles are detected with an on-board camera. Objects of a
given size are used as obstacles to allow for the estimation of the distance
and the angle from the robot. Obstacles are detected in a range of 0.8 m and a
field-of-view of 62°. The software runs on Robotic Operating System (ROS)
Noetic. Localization data, which consists of the 2D coordinates of the rover
($x_{r}$, $y_{r}$) and its heading ($\theta_{r}$), is shared with the path
planning node. The distances, angles, and sizes of all detected obstacles are
transmitted upon request from the planning node.

We completed a set of tests with different starting points, targets, and
obstacle layouts. Arrangements that produce local minima are also tested to
ensure the robustness of the RAPF algorithm. Tests have been conducted with up
to 16 obstacles of varying sizes to comply with the lunar model presented
before. Different terrain topologies were also tested (flat surface,
overcoming small rocks, and slopes). Twenty tests were performed using virtual
obstacle detection. This means that instead of using the on-board camera to
detect obstacles, their positions are measured offline and provided to the
planning algorithm. The detection is simulated in software as a perfect
detection since all potential noise or uncertainty is removed. Simulated
object detection assumes the same range and field-of-view of the on-board
camera (i.e., 0.8 m, 62°). This enables us to evaluate the performance of the
path planning algorithm irrespective of potential errors induced by other
parts of the navigation pipeline. Ten tests were conducted using the on-board
camera to prove that RAPF is robust to sensor-induced noise and uncertainties.
Details of each test are provided in the following section.

###  IV-B Testing results and discussion

Robustness: In a first set of tests, we evaluated the robustness of RAPF
against variations in the environment. Ten field tests are performed for two
different obstacle layouts. These tests are performed with the virtual
obstacle detection configuration only to avoid introducing further
uncertainties into the analysis. Obstacles, starting points, and target
layouts for the two scenarios, together with the computed and walked paths,
are plotted in Fig. 2.

Figure 2: Mean computed path (red) , mean walked path (black) and standard
deviation of the mean walked path (gray) for 10 repetitions for two different
scenarios (Path1 & Path 2)

A t-test is performed to compare the ratio ($\Omega$) of the actual walked
path length over the computed path length for two different scenarios, to
determine if the RAPF algorithm’s behavior is influenced by the environment,
with a null hypothesis suggesting no significant difference in the ratio of
computed path length to actual walked path length between two map
configurations. For a sample size of 10 for each group, a calculated
t-statistic of 0.346 and a critical value of 2.110 for a two-tailed test on 17
degrees of freedom fails to reject the null hypothesis, indicating no
statistical difference in the analyzed ratio.

The main challenge brought by real-world testing is the uncertainties. In a
second set of tests, we evaluated the robustness of RAPF against sensor-
induced uncertainties. In this case, data is collected for both rover
configurations (virtual and on-board obstacle detection) and across Path 2. A
Welch’s t-test is employed with a sample size of 10 for each group, resulting
in a calculated t-statistic of 0.64 and a degree of freedom of 15. With a
p-value of 0.05 and a critical t-value of 2.131, the null hypothesis,
suggesting no statistical difference in path length between the two scenarios,
cannot be rejected as the calculated t-value is smaller than the critical
value.
Lastly, it is worth highlighting that the minimum distance maintained between
the rover and the nearest obstacle was 7cm on average when the on-board
obstacle detection was used and 9cm when the virtual detection was used. The
minimum distance kept from all obstacles was 20cm on average for the on-board
detection and 26cm for the virtual detection. The hardware detection shows the
lowest median and largest IQR, which means that when RAPF is implemented
alongside on-board obstacle detection, the robot tends to travel closer to
obstacles along its path. The on-board hazard detection leads to a smaller
safety distance due to the greater variability in the estimation of the
position of the obstacles relative to the rover.

Computational complexity: The path planning algorithm is run on the rover on-
board computer, a Raspberry Pi 4B 2GB. The total computation time from the
starting point to the target for the environments previously defined was
always less than 0.2 s. The time to recompute one path took a maximum of 0.025
s.

For APF-based path planning methods, memory usage and CPU exploitation are
intrinsically related to the size of the potential map. The smaller the map,
the lower the overall computational resource needed to compute a path. In
APFs, a map of the surroundings of the rover is stored in memory. Every time a
new waypoint is reached, all path planning data is cleared from memory and new
paths are recomputed. The RAPF algorithm is distinguished by the fact that at
bacteria points only the potential is computed. The number of planned paths
depends on the level of cluttering of the environment as a new path is
recomputed every time an obstacle is detected. For each new planned path, the
number of potential values computed can be derived from:

$N_{RAPF}=\frac{N_{B}M_{s}}{\rho_{B}},$ (21)

where $N_{B}$ is the number of bacteria, $M_{s}$ is the path length, and
$\rho_{B}$ is the distance between the rover and the bacteria point. The RAPF
execution speed highly relies, therefore, on the distance between the rover
and the targeted waypoint. Our field experiments showed that to reach a
waypoint located about 3 m from the rover, the number of computed potential
values for RAPF is 610 on average. In comparison, the estimated number of
computed potential values for the grid-based APF and RVF algorithms would be
between 1100 and 1300 for a 5x5$m^{2}$ stored map. Each obstacle is stored in
6 bytes. The size of the list of obstacles stored depends on the number of
obstacles $N_{o}$. From the lunar model derived, up to 35 obstacles can be
considered within a map of 5x5 $m^{2}$. The computed path to follow is also
made of float16 values. It consists of a list of coordinates
($X_{r}$,$Y_{r}$), each couple having a size of 4 bytes. The total size
depends on the length of the path to follow. For a 5x5 $m^{2}$ and a waypoint
located at around 3 m from the rover, RAPF approximately requires 0.5 KB of
storage. In comparison, for storing a map of the same size, Field D* needs
around 10 KB of memory [22].

Reachability: The path planning module running onboard is capable of
estimating the average time it would take to complete a traced path. In
combination with the relative location of the rover and the current drawn by
each motor, this can help the rover know if it is physically stuck or blocked.
The algorithm knows the number of walking and turning steps required to reach
a waypoint as well as the number of times a picture needs to be taken to check
the path and the time needed to perform these tasks. These data are subject to
slippage, actual resource allocation, and terrain uncertainties. Therefore a
100% margin is added to this estimated time. If the designated waypoint is not
reached within the allocated time for executing the path, an error flag is
raised and notified to the on-board computer. Data gathered by the rover while
navigating are plotted in Fig.3. Potential maps as displayed in this figure
are not stored on the rover and only computed for bacteria points; only the
list of obstacles is stored on board. The potential maps plotted here are for
environmental visualization purposes only. Regions between the obstacles and
$\rho_{0}$ normally have infinite potential, represented in these plots by a
very large finite value. On the left image, no local minima escape method is
implemented (i.e., CRBAPF path planner). The rover gets stuck at the first
local minimum encountered. On the right image, the RAPF algorithm is used,
which implements a potential function assigned to local minima as described
before. The path leading to the local minimum is cleared and a new path is
computed considering the potential generated by the local minimum. This leads
to a new path capable of reaching the target, as can be seen in the figure. In
some cases, the algorithm does not manage to plan a path that avoids the local
minimum despite the local minima technique implemented. In this case, an error
flag ’local_minimum’ is set by a timer. If no local minimum-free path is
output after 10 seconds, the path planning algorithm notifies the OBC that no
path can be found. This flag prevents the rover from continuing to search an
impossible path and avoids unnecessary resource consumption.

Figure 3: Potential map, computed path (orange), and path followed by the
rover (black). Left: CRBAPF. Right: RAPF

##  V Summary

The RAPF algorithm enables real-time computation of new paths for small and
resource-limited mobile robots. The proposed approach drastically increases
the success rate compared to state-of-the-art APF methods by effectively
addressing the local minima issue. The bacterial approach of the RAPF
algorithm, minimizes memory requirements, making it ideally suited for robotic
systems with strong computational limitations. Simulation results show a 200%
higher success rate and 50% lower computing time than the classical APF
algorithm, with only a 5% longer path length than the optimal algorithm, A*,
in a lunar-like environment. It also outperforms state-of-the-art APF-based
algorithms by more than 15% in reachability and 10% in path length for a
similar or shorter planning time. During the field validation, the RAPF
algorithm displayed robust performance, effectively adapting to uncertainties
introduced by various sensors. The RAPF algorithm performs well in various
situations. The rover walks in a straight line towards the target when no
obstacles are detected. It navigates smoothly around small discrete obstacles.
The path planning performance of the rover do not depend on the relief of the
terrain. Its consistent performance across diverse environments underscores
its capability for real-time autonomous navigation over challenging terrains.

## Acknowledgment

We express our appreciation to Lunar Zebro [23] for providing the rover
prototype and granting access to their testing facility. The project is
partially funded by the Moonshot program supported by the Delft Space
Institute.

## References

  * [1] X. Yu, Q. Huang, P. Wang and J. Guo, ”Comprehensive Global Path Planning for Lunar Rovers,” 2020 3rd International Conference on Unmanned Systems (ICUS), Harbin, China, 2020, pp. 505-510
  * [2] Peng, Q., Wang, S., Zhang, J., ”Fast Lunar Path Planning Based on Intelligent Heuristic Map,” International Conference on Autonomous Unmanned Systems 2022.
  * [3] Dijkstra, E., ”A note on two problems in connexion with graphs,” Numerische Mathematik, vol. 1, 1959
  * [4] P. Hart, N. Nilsson, B. Raphael, B., ”A Formal Basis for the Heuristic Determination of Minimum Cost Paths,” IEEE Transactions on Systems Science and Cybernetics, vol. 4, pp. 100-107, 1968
  * [5] LaValle, S., ”Rapidly-exploring random trees : a new tool for path planning,” The annual research report, 1998
  * [6] L. Kavraki, P. Svestka, J. Latombe, M. Overmars, ”Probabilistic roadmaps for path planning in high-dimensional configuration spaces,” IEEE Transactions on Robotics and Automation, vol. 12, 1996
  * [7] S. Karaman, E. Frazzoli, ”Sampling-based Algorithms for Optimal Motion Planning”, International Journal of Robotics Research, 2011
  * [8] W. McCulloch, W. Pitts, ”A logical calculus of the ideas immanent in nervous activity,” Bulletin of Mathematical Biophysics 5, 115–133, 1943
  * [9] J. Kennedy, R. Eberhart, ”Particle swarm optimization”, Proceedings of ICNN’95 - International Conference on Neural Networks, vol. 4, pp 1942-1948, 1995
  * [10] M. Dorigo, M. Birattari, T. Stutzle, ”Ant colony optimization,” IEEE Computational Intelligence Magazine, vol. 1, pp. 28-39, 2006
  * [11] S. Katoch, S. Chauhan, V. Kumar, ”A review on genetic algorithm: past, present, and future,” Multimed Tools Appl 80, 8091–8126, 2021
  * [12] Holland, John H. “Genetic Algorithms.” Scientific American, vol. 267, no. 1, 1992, pp. 66–73. JSTOR, http://www.jstor.org/stable/24939139. Accessed 16 Jan. 2024.
  * [13] Khatib, O., ”Real-time obstacle avoidance for manipulators and mobile robots,” Proceedings. 1985 IEEE International Conference on Robotics and Automation, vol. 2, pp 500-505, 1985
  * [14] L. Nasuha, L. Priambodo, G. Pratama, “Vortex artificial potential field for mobile robot path planning,” Journal of Physics: Conference Series, vol. 2406, 2022.
  * [15] Gelling, L., “A lunar surface path finding and obstacle avoidance algorithm,” M.S. thesis, Delft University of Technology, 2023.
  * [16] A. Hossain, I. Ferdous, “Autonomous robot path planning in dynamic environment using a new optimization technique inspired by bacterial foraging technique,” Robotics and Autonomous Systems, vol. 64, pp. 137–141, 2015.
  * [17] Ahmad, Y. “Matlab implementation of artificial potential field,” https://github.com/Yaaximus/artificial-potential-field-matlab
  * [18] M. Diab, M. Mohammadkarimi, R. Rajan, “Artificial potential field-based path planning for cluttered environments,” IEEE Aerospace Conference, 2023, pp. 1–8
  * [19] M. Golombek, D. Rapp, ”Size-frequency distributions of rocks on Mars and Earth analog sites: Implications for future landed missions,” Journal of Geophysical Research: Planets, vol. 102, pp. 4117-4129, 1997
  * [20] M. Golombek, A. Haldemann, N. Forsberg-Taylor, E. DiMaggio, R. Schroeder, B. Jakosky, M. Mellon, J. Matijevic, ”Rock size-frequency distributions on Mars and implications for Mars Exploration Rover landing safety and operations,” Journal of Geophysical Research: Planets, vol. 108, 2003
  * [21] Y. Li, B. Wu, ”Analysis of Rock Abundance on Lunar Surface From Orbital and Descent Images Using Automatic Rock Detection,” Journal of Geophysical Research: Planets, vol. 123, 2018
  * [22] J. Carsten, A. Rankin, D. Ferguson, A. Stentz, ”Global Path Planning on Board the Mars Exploration Rovers”, IEEE Aerospace Conference, 2007, pp. 1-11
  * [23] R.T.Rajan, A.Noroozi, P.Sachdeva, A.Menicucci, C.Verhoeven, ”Lunar Zebro - An Autonomous Moon Rover”, 75th International Astronautical Congress, 2024


&\geq \Pr\left[ \alpha + \sum_{i=1}^{T(1-\alpha)} \frac{y^v_i}{T} < \theta_v \bigwedge  \alpha + \sum_{i=1}^{T(1-\alpha)} \frac{y^u_i}{T} < \theta_u\right] \cdot \sum_{0\leq \beta\leq \alpha} \Pr[\tilde x_e = \beta ]\\
&\geq \sum_{0\leq \beta \leq \alpha} \Pr\left[ \frac{2\tau}{p} + \sum_{i=1}^{T(1-\beta)} \frac{y^v_i}{T} < \theta_v \bigwedge \frac{2\tau}{p} + \sum_{i=1}^{T(1-\beta)} \frac{y^u_i}{T} < \theta_u \right] \cdot \Pr[\tilde x_e = \beta]\\
&= \sum_{0\leq \beta \leq  \alpha} \Pr\left[ \frac{2\tau}{p} + \tilde x_{v \setminus e} < \theta_v \bigwedge \frac{2\tau}{p} + \tilde x_{u \setminus e }  < \theta_u \mid \tilde x_e = \beta\right] \cdot \Pr[\tilde x_e = \beta]\\
&= \Pr\left[\frac{2\tau}{p} + \tilde x_{v \setminus e} < \theta_v \bigwedge \frac{2\tau}{p} + \tilde x_{u \setminus e} < \theta_u \mid \tilde x_e < \frac{2\tau}{p}\right].
\end{align*}

The first inequality holds becasue $\sum_{0\leq \beta\leq \alpha}\Pr[\tilde x_e = \beta ] \leq 1$. The second ineqialty holds because $\frac{2\tau}{p}>\alpha $. Moreover, the first equality and the second last equality follows from the fact that $y_i^v$ for all $i \in[T]\setminus Q^e$ are independent and identically distributed.  Therefore, by defination of $s_e$ and $s_e'$, we obatin,

\[
\Pr[s_e \geq \theta \mid \tilde x_e = \alpha] \geq  \Pr\left[s'_e \geq \theta \mid \tilde x_e < \frac{2\tau}{p} \right]
\]

and so

\[
\E\left[s_e \mid \tilde x_e=\alpha\right] \geq \E\left[s'_e \mid \tilde x_e < \frac{2\tau}{p}\right]\geq \E[s_e'] - \eps \cdot q_e.
\]

We are now ready to prove the main theorem.

For the sake of exposition, we let $\alpha =\frac{2\tau}{p}-\frac{1}{T}  $.

\begin{align*}
\Pr[e\in M_\nc] &= \E_{\tilde x_e}[ \tilde x_e \E[s_e \mid \tilde x_e]] \cdot \Pr[e\in R]\\
&\geq \sum_{\tilde x_e \leq \alpha}\tilde x_e \cdot \E[s_e \mid \tilde x_e ] \Pr[\tilde x_e] \cdot \Pr[e\in R]\\
&\geq \sum_{\tilde x_e \leq \alpha}\tilde x_e \cdot \E\left[s_e \mid \tilde x_e =\alpha  \right] \Pr[\tilde x_e] \cdot \Pr[e\in R]&& (\text{Claim~\ref{claim:monotonicitys-e}})\\
&\geq \E_{\tilde x_e} \left[\tilde x_e \mid \tilde x_e \leq \alpha \right] \cdot \E\left[s_e \mid \tilde x_e =\alpha  \right] \cdot \Pr[e\in R]\\
% &= \E\left[s_e \mid \tilde x_e =\frac{2\u}{p} \right]\cdot \E_{\tilde x_e} [\tilde x_e ]\\
&\geq \E\left[s_e \mid \tilde x_e =\alpha  \right]\cdot (1-3\eps) \cdot q_e && (\text{Claim~\ref{claim:a2}})\\
&\geq (\E\left[s'_e\right] - \eps \cdot q_e) \cdot (1-3\eps) \cdot q_e &&  (\text{Claim~\ref{claim:monotonicitys}})\\
&\geq (1-2\eps )^2\cdot (1-3\eps) \cdot q_e - \eps \cdot q_e &&(\text{Claim~\ref{claim:s-eis'1}})\\
&\geq (1-12\eps)q_e.
\end{align*}

Lastly, we will prove the monotonicity property (<ref>) of the matching constructed by Algorithm <ref>.
Let $e\in \nc$ be an arbitrary non-crucial edge and $S \subseteq N(e)$ be an arbitrary subset of edges incident to $e$. Note that when we condition on $Q_1 \dots Q_T$ then $f_e := \sum_{i=1}^T \mathbbm{1}[e \in \optsol(Q_i)]/T$ becomes deterministic for all non-crucial edges $e \in \nc$. For the rest of the proof, we condition on the set $Q_1,\dots, Q_T$.

Recall that $\Pr[e \in M_\nc] = \Ex[\tilde x_e \cdot s_e]$, therefore it is enough to show that $s_e$ weakly increases once we condition on $S \cap R=\emptyset$, i.e. all edges in the set $S$ are inactive. Note that the probability and expectation are taken over the random set $R$.  Formally, we let $\mathcal D$ be the distribution of active edges $R\subseteq E$ and $\D'$ be the distribution of active edges $R\subseteq E\setminus S$ conditioned on $S\cap R = \emptyset$. We will show that,

\begin{equation}
\label{eq:stoc_dominance}
\E_{R \sim \mathcal D}[s_e] \leq \E_{R' \sim \mathcal D'}[s_e].
\end{equation}

We couple two probability distributions $\mathcal D$ and $\mathcal D'$ as follows: we first sample $R \sim \mathcal D$ and construct $R' = R \setminus S$. Since each edge $e\in E$ appears in the active set of edges $R\sim \D$ independently with probability $p$, it follows that $R'\sim \D'$. Therefore the coupling $(\mathcal D, \mathcal D')$ is valid with the correct marginals. We denote this joint distribution of $(\D,\D')$ as $\gamma$.

For any sample $(R, R')\sim \gamma$, from the definition of $\tilde x_{v} ( R) \geq \tilde x_v ( R')$, where $\tilde x_v (R)$ denotes the value of $\tilde x_v$ conditioned on the set of active edges being $R$. Note that $\tilde x_e ( R)$ and $\tilde x_v (R')$ are deterministic as we have already conditioned on $Q_1 ,\dots , Q_T$. This further implies that $s_e ( R )\leq s_e (R')$ where $\tilde x_v ( R')$, $s_e ( R )$ and $ s_e (R')$ are defined similar to the $\tilde x_v (R)$.  This completes the proof.

§.§ Proof of Proposition <ref>

Given an instance of stochastic weighted matching $\match$ with $p\leq \eps^4$ and $0<\eps<1/2$, for any non-crucial edge $e=(u,v)\in \nc$ with $N(e) = \{e_{1},\dots, e_{k}\}$ be the incident edges on vertices $u$ and $v$ in the graph $G$,
$$\Pr\left[ Q_{\crs}\cap N(e) \leq \frac{(1+\eps)\cdot 2}{p}\right] \geq 1-\epsilon$$

Let $e:=(u,v) \in \nc$ be an arbitrary non-crucial edge of the problem $\match$. For each edge $e_i\in N(e)$, let $Z_{i}$ be the indicator random variable which denotes whether the edge $e_i$ belongs to $Q_{\crs}$ or not, i.e., $Z_{i} = \mathbbm{1}[e_i \in Q_{\crs} \cap R]$. Note that $Z_i\sim Ber(\frac{q_e}{p})$ and independent across $i\in [k]$ and $\sum_{i\in [k]}\E[Z_i] = \sum_{i\in [k]} q_e /p \leq 2/p$. Therefore, by setting $\eps' = \frac{2(1+\eps) - \sum_{i\in [k]} q_{e_i}}{\sum_{i\in [k]} q_{e_i}}$ and applying Chernoff bound, we get

\begin{align*}
\Pr\left[\sum_{i=1}^k Z_i \geq (1+\eps')\cdot \sum_{i=1}^k \frac{q_{e_i}}{p}\right] \leq \exp\left(-(\eps')^2 \cdot  \frac{\sum_{i\in [k]}q_{e_i}}{3p}\right).
\end{align*}

Now, as $\eps ' \geq \frac{2\cdot \eps}{\sum_{i\in [k]} q_{e_i}}$, we obtain,

\begin{align*}
\Pr\left[\sum_{i=1}^k Z_i \geq (1+\eps')\cdot \sum_{i=1}^k \frac{q_{e_i}}{p}\right]&= \Pr\left[\sum_{i=1}^k Z_i \geq (1+\eps)\cdot  \frac{2}{p}\right] \\
&\leq   \exp\left(-4\eps^2 \cdot  \frac{1}{3\eps^4\sum_{i\in [k]} q_{e_i}}\right)\\
&\leq \exp\left( -\frac{2}{3\eps^2}\right) \leq \eps.
\end{align*}

We can express $\Pr[\mathcal E_e]$ as:

\begin{align*}
\Pr[\mathcal E_e] &= \Pr[\mathcal E_e\mid e\notin  Q_{\crs}]\cdot \Pr[e\notin  Q_{\crs}] +  \Pr[\mathcal E_e\mid e\in  Q_{\crs}]\cdot \Pr[e\in  Q_{\crs}] \\
&\leq \Pr[\mathcal E_e\mid e\notin  Q_{\crs}] + \Pr[e\in  Q_{\crs}] \leq \Pr[\mathcal E_e\mid e\notin  Q_{\crs}]  + \epsilon, \\
\end{align*}

where the last inequality follows due to definition of non-crucial edges, i.e. $\Pr[e\in Q_{\crs}] =q_e/p \leq \epsilon$. Combining above inequality with Claim <ref>, we obtain

\begin{equation*}
\Pr[\mathcal E_e\mid e\notin  Q_{\crs}] \geq \Pr[\mathcal E_e] -\epsilon \geq 1-2\epsilon.
\end{equation*}


#  Adversarial Predictions of Data Distributions Across Federated Internet-of-
Things Devices

Samir Rajani1, Dario Dematties12, Nathaniel Hudson34, Kyle Chard34, Nicola
Ferrier12,
Rajesh Sankaran12, and Peter Beckman12  1Northwestern Argonne Institute of
Science and Engineering, Northwestern University; Evanston, IL, United States
2Mathematics and Computer Science Division, Argonne National Laboratory;
Lemont, IL, United States  3Department of Computer Science, University of
Chicago; Chicago, IL, United States  4Data Science and Learning Division,
Argonne National Laboratory; Lemont, IL, United States

###### Abstract

Federated learning (FL) is increasingly becoming the default approach for
training machine learning models across decentralized Internet-of-Things (IoT)
devices. A key advantage of FL is that no raw data are communicated across the
network, providing an immediate layer of privacy. Despite this, recent works
have demonstrated that data reconstruction can be done with the locally
trained model updates which are communicated across the network. However, many
of these works have limitations with regard to how the gradients are computed
in backpropagation. In this work, we demonstrate that the model weights shared
in FL can expose revealing information about the local data distributions of
IoT devices. This leakage could expose sensitive information to malicious
actors in a distributed system. We further discuss results which show that
injecting noise into model weights is ineffective at preventing data leakage
without seriously harming the global model accuracy.

###### Index Terms:

Federated Learning, Data Leakage, Internet-of-Things, Data Privacy, Sensor
Networks

##  I Introduction

Federated Learning (FL) is a method for decentralized machine learning which
can seamlessly train models across many Internet-of-Things (IoT) devices. In
conventional decentralized machine learning, IoT devices send their training
data to a central location (e.g., remote cloud server) to train a deep neural
network or some other model. This approach is unattractive due to its high
communication cost [1] and its requirement that all data—including possibly
sensitive data—be shared across a network [2]. Under FL, no raw data are
communicated over the network. Instead, the central server is sent only
locally trained model weights provided by the IoT devices. In short, FL begins
with a central server initializing a “global” machine learning model with
random weights, and then the FL training loop begins. The FL training loop,
depicted in Fig. 1, is as follows:  (i) the central server shares the current
global model with the IoT devices;  (ii) the IoT devices receive a copy of the
global model and train it on their private data;  (iii) the IoT devices send
back their locally trained copies of the model to the central server;  (iv)
the central server aggregates all the received locally trained models to
update the global model; and the loop repeats.

Because no raw data are communicated over the network in FL, it is a promising
method for training a model on sensitive or private data. However, since FL
requires that clients send their model updates across a network, the extent to
which FL actually protects client data is unclear.

Client models are trained on local data, meaning they could contain implicit
revealing information about private training samples. Recent work has shown
that gradients from a small batch of samples can be used to produce pixel-wise
reconstructions of images and token-wise reconstructions of text [3]. While
these findings are alarming, their application is limited to settings in which
the gradients being transmitted are from a single, small batch of at most
eight samples [3], or to a limited extent when gradients are averaged over a
larger number of samples [4]. Other work [5] has shown that it is possible to
analytically reconstruct ground-truth labels for training samples from
gradients, but this technique is viable only when gradients are known for
every single training sample. By contrast, in typical FL settings, client
models undergo one or more epochs of training over their local data before
their parameters are sent to the server, rather than just a single gradient
update on a small batch of samples.

Figure 1: The standard steps for FL process: global model distribution, client
model training, client model communication, and aggregation. The steps in red
outline the procedure by which an adversary might conduct an attack.

In this work, we show that even in realistic federated settings, sharing
client model gradients or parameters leaks private information about training
samples. In particular, we show that an adversary could train a deep neural
network to extract a client’s label distribution from its shared model
parameters with near-perfect accuracy. In many federated settings, the ability
of an adversary to predict data distributions of IoT devices can be dangerous.
Consider, for example, a federated language modelling task in which the words
of text messages sent by users of mobile phones serve as labels. In this
setting, it is conceivable that an adversary who is able to reconstruct the
labels of the samples on client devices might be able to deduce identifiable
information about a user. An adversary might, for example, be able to make a
list of clients whose text messages contain information related to banking,
and then target these clients in the future, with more refined attacks.
Moreover, the clients participating in FL might be under the false assumption
that their privacy is assured because training samples remain on their device.
To understand the extent to which user data is protected, we need a deeper
understanding of what kinds of attacks are possible and what kinds of
mitigation strategies are effective.

Our contributions can be summarized as follows: 1) we show empirically that
client model parameters contain information about the labels of samples on
which they were trained; 2) we perform several experiments that suggest
explanations for the presence of this information; 3) we demonstrate that a
deep neural network can be trained to recover information about the
distributions of client labels; and 4) we test the ability of simple defenses,
such as adding Gaussian and Laplacian noise to gradients, to prevent the
leakage of information about client label distributions. Our findings have
significant implications for the use of FL in privacy-sensitive settings.

##  II Related Work

Much of the popularity of FL can be accredited to the intuitive idea that it
improves data privacy in decentralized machine learning since no raw data are
communicated across the network [6, 7, 8]. However, recent works have begun to
challenge this assumption [5]. Zhu et al. coin the notion of “deep leakage
from gradients” to refer their demonstrated optimization-based approach to
reconstruct data from gradient updates [3]. Similarly, Geiping et al.
demonstrate a numerical reconstruction method which reconstructs images from a
neural network gradient update [4]. However, these works have to relax the
problem of data reconstruction substantially to demonstrate the effectiveness
of their techniques. Specifically, their approaches are most effective when
the gradient update from training the model is done on a single image. They
both note that multi-image reconstruction can be done but that effectiveness
of the technique suffers when gradients are averaged across modest batch sizes
$B\geq 8$.

These growing concerns around data leakage from model updates have encouraged
exploration into additional privacy-preserving techniques for FL. Some of the
most popular include homomorphic encryption [9, 10] and differential privacy
[11, 12, 13, 14]. Homomorphic encryption is a technique that encrypts data in
such a way that they are equivariant to many mathematical transformations.
This means that the global model can be updated by an average of encrypted
model weights. These averaged model weights can then be decrypted to train the
models locally without corrupting the local models [15]. A central concern
surrounding homomorphic encryption for FL is that the resource cost of
encrypting data on resource-constrained IoT devices is non-trivial [16, 17].
Differential privacy takes a different, less resource-costly approach. This
technique directly injects an amount of noise, $\epsilon>0$, into client model
updates before they are sent to aggregation node [18, 19]. While cheaper to
perform, there is an obvious trade-off to consider. If $\epsilon$ is too
large, the global model’s accuracy will suffer; if $\epsilon$ is too small,
privacy benefits are reduced.

Our work takes a fundamentally different approach from prior works in the FL
privacy literature. While the works on data reconstruction from model
gradients are critically important in the FL domain, we argue that they rely
on unrealistic simplifying assumptions related to the amount of data model
updates are trained on to support the reconstruction. Instead, we consider a
more realistic scenario where the amount of data on IoT devices is larger than
a handful of data samples (i.e., $\gg 8$). We then show that even with this
more challenging scenario, the local data distributions of IoT devices in an
FL system can still be approximated based on the model weights.

##  III Problem Description & Research Questions

We consider a standard two-tier FL system comprised of decentralized IoT
devices with local, private data and a central aggregation node which manages
the FL process. We seek to learn the extent to which client information is
compromised by sharing client gradients. In particular, could an adversary use
the model parameters shared by clients to deduce information about their
private training data?

Further, we study an adversary who has access to information about the
federated training (e.g., an honest-but-curious server). In the case of a
supervised objective, we assume the adversary has access to the dataset labels
and their order in the output layer of the network. In the case of an
unsupervised objective, as we will show in Section IV-A, models
representations contain relevant information about the input features of the
data, and an adversary could potentially still deduce information about client
data distributions. We also assume the adversary has access to the client
optimization parameters, such as the learning rate $\eta$ and the number of
local training epochs $E$. We focus on the possibility that an adversary could
deduce information about the distribution of a client’s training samples.
Client training data are often not independently and identically distributed
(non-IID), and their data distributions may contain private information.

Specifically, we are interested in answering the following questions:  Q1.
What is the relationship between client model parameters and the data
distribution on which the parameters were learned?  Q2. Can a client model’s
parameters be used directly to predict the distribution on which they were
learned?  We expand on our approach to each of these questions below.

###  III-A Visualizing Client Model Parameters

To understand the relationship between client model parameters and the data
distribution on which they were learned, we train clients on a variety of data
distributions and visualize their model parameters using principal component
analysis (PCA). We refer to clients trained on synthetic data distributions as
“dummy clients,” and to the space of dimension-reduced model parameters as a
“model-latent space.”

In Section IV, we show several experiments which explain the relationship
between the position of a client model in the model-latent space to the
client’s data distribution. By performing a layer-by-layer decomposition of
the client model parameters, and by examining the behavior of client model
parameters when the training objective is unsupervised, we also propose an
explanation for this clustering.

###  III-B Predicting Client Data Distributions

The existence of a relationship between client model parameters and the data
distribution used for learning them raises a second question: can a model’s
position in the model-latent space be used directly to predict a client’s
label distribution? To answer this question, we introduce the notion of a
“meta-dataset.” For client models trained on disjoint training samples
$\mathcal{P}_{k}$ from some original dataset, the meta-dataset is a mapping
from client model parameters $w_{k}$ to label distributions
$\text{Pr}(y=y_{i}\text{ }|\text{ }i\in\mathcal{P}_{k})$. We propose that an
adversary might be capable of training a deep neural network on a synthetic
meta-dataset, enabling them to extract information about clients.

In Algorithm 1, we propose one potential strategy an adversary might use to
obtain information about client data distributions. Later, in Section IV, we
validate empirically the effectiveness of this strategy in the presence of a
simple defense: the addition of Gaussian and Laplacian noise to gradients. The
adversary first intercepts the global model $w_{t}$ sent by the server to the
clients (line 1). They then construct a meta-dataset by creating a large
number of dummy clients whose local training data comes from a synthetic non-
IID distribution (lines 3-6). The client models are initialized with the
intercepted global model parameters and trained on their local training data.
Then, their trained parameters $\hat{w}^{k}$ are reduced to some small number
of dimensions using PCA (lines 8-10). Each dummy client provides one training
sample in the meta-dataset, in which the PCA-reduced model parameters are the
inputs, and the label distributions are the outputs. The adversary then trains
a neural network on the meta-dataset (lines 11-14); the dimension of the PCA-
reduced model parameters $\hat{w}^{k}_{\text{pca}}$ should be sufficiently
small such that they can fit in the input layer. A key assumption is that, in
order to train the “dummy” client models, the adversary can access a set of
proxy training samples that are similar to the training samples on the client
devices. Upon intercepting model parameters sent from the clients to the
server, they could apply the same PCA model to these parameters, and use their
trained label distribution predictor to perform inference (lines 15-17).

1 intercept global model weights, $w_{t}$

2 $\hat{S}^{t}\leftarrow$ set of dummy clients

3 foreach  _dummy client $k\in\hat{S}^{t}$_ do

4 $\hat{\mathcal{P}}^{k}\leftarrow\text{(random subset of proxy data)}$

5 $\hat{w}^{k}\leftarrow\mathit{ClientUpdate}(k,w_{t})$

6 $y^{k}\leftarrow\text{Pr}(y=y_{i}\text{ }|\text{ }i\in\mathcal{P}_{k})$

7

8

9$W_{\text{pca}}\leftarrow\text{(fit PCA model on the $\hat{w}^{k}$)}$

10

11foreach  _dummy client $k\in\hat{S}^{t}$_ do

12 $\hat{w}^{k}_{\text{pca}}\leftarrow\text{(apply $W_{\text{pca}}$ to
$\hat{w}^{k}$)}$

13

14initialize label distribution predictor $f_{\theta}$

15 $\mathcal{B}\leftarrow$ (split the $(\hat{w}^{k}_{\text{pca}},y^{k})$ into
batches of size $\hat{B}$)

16 foreach  _local epoch $e=1,\ldots,\hat{E}$_ do

17 foreach  _batch $b\in\mathcal{B}$_ do

18 $\theta\leftarrow\theta-\eta^{\prime}\nabla\hat{\ell}(\theta;b)$

19

20

21foreach  _client $k\in S_{t}$_ do

22 intercept client model weights $w_{t+1}^{k}$

23 $y^{k}_{\text{(pred)}}\leftarrow f_{\theta}(w_{t+1}^{k})$

24

25

1 Procedure  _ClientUpdate $(k,w)$_

2 $\mathcal{B}\leftarrow\text{(split $\hat{\mathcal{P}}_{k}$ into batches of
size $B$)}$

3 foreach  _local epoch $e=1,\ldots,E$_ do

4 foreach  _batch $b\in\mathcal{B}$_ do

5 $w\leftarrow w-\eta\nabla\ell(w;b)$

6

7

8 return $w$ to server

9

Algorithm 1 Adversarial data distribution prediction

##  IV Experimental Design & Results

###  IV-A Visualizing Client Model Parameters

As a starting point for understanding the relationship between client model
parameters and their label distributions, we replicate and extend an
experiment performed by Wang et al. [20]. We initialize a global model, which
is broadcast to 100 clients with disjoint training sets. For each of the ten
labels in the dataset, ten clients are provided with a synthetic non-IID
training set, according to one of two sampling schemes. In the “80%-20%
sampling” scheme, we provide each client with a data distribution in which 80%
of the samples have one label, and the remaining 20% are uniformly distributed
across the remaining labels. In the “Dirichlet sampling” scheme, client data
distributions are sampled from a Dirichlet distribution. Each client is
trained for one epoch, and all client model parameters are reduced to two
dimensions via PCA for visualization.

We perform experiments using supervised classification tasks on the MNIST
dataset using a multi-layer perceptron (MLP), and on the CIFAR-10 dataset
using a convolutional neural network (CNN); the results of these experiments
are shown in Fig. 2. For both datasets and architectures, we observe the same
phenomenon: client model parameters are clustered according to their dominant
labels, suggesting that they contain implicit information about the samples on
which they were learned. We might hypothesize that this clustering originates
from one or both of two distinct phenomena:

  1. 1.

Images belonging to the same class have similar features, leading clients with
the same dominant label to train models with similar parameters.

  2. 2.

Images belonging to the same class have the same labels, and the nature of the
supervised classification task means client models trained on sample
distributions with the same dominant label will have similar parameters. Zhao
et al. demonstrate a related phenomenon analytically for gradients from a
single training sample [5].

Figure 2: Clustering of PCA-reduced client model parameters on the MNIST and
CIFAR-10 datasets according to the dominant label in their training sets. In
the “80%-20% sampling” scheme, the dominant label accounted for 80% of a
client’s training samples. In the “Dirichlet sampling” scheme, client data
distributions were sampled from a Dirichlet distribution with $\alpha=0.2$;
smaller points denote that a smaller percentage of samples come from the
dominant class, and these points appear to drift from the centers of their
clusters. These experiments used a small learning rate ($\eta=10^{-5}$), which
we found to show the clearest clustering, but the same phenomenon was observed
using all learning rates we tested.

If the clustering is a result of the first phenomenon, we would expect that
client models with semantically similar dominant labels (e.g., car and truck)
would be close together in the model-latent space, and that clustering would
occur even with an unsupervised training objective. Meanwhile, if the
clustering is a result of the second phenomenon, we would expect the positions
of dominant-label clusters to be arbitrary.

To understand which of these phenomena accounts for the clustering, we perform
two more experiments. First, we examine how the model-latent space behaves
when we study each layer of the neural network in isolation. Fig. 3 shows the
PCA-reduced client model parameters for each layer of the CNN used for
CIFAR-10. We find in the earlier convolutional layers, client models trained
on distributions with semantically similar dominant labels (e.g., cat and dog,
car and truck, plane and ship) are close together in the model-latent space.
In the final fully-connected layers, this semantic clustering disappears,
suggesting that the clustering is primarily a result of the labels provided in
the supervised classification task. We hypothesize, therefore, that the
formation of dominant-label clusters is a result of both phenomena we
identified, the first manifesting in the earlier layers, which extract
features from input samples, and the second manifesting in the later layers,
designed for the classification task.

Figure 3: Layer-by-layer clustering of PCA-reduced client model parameters on
the CIFAR-10 dataset. Semantic clustering appears to occur in the earlier
convolutional layers, suggesting that the clustering is the result of the
similarity in features of images belonging to the same class. Note, for
example, that client models trained primarily on images of animals are
concentrated in the bottom-left corner of the model-latent space for the first
convolutional layer. Meanwhile, the clustering appears roughly arbitrary in
the final fully-connected layer, suggesting that the clustering here is a
result of labels.

To validate this hypothesis, we perform a second experiment, in which we train
an autoencoder on the MNIST dataset. If part of the clustering phenomenon can
be accounted for by feature similarity, we would expect some clustering to
occur even with an unsupervised task. Indeed, as shown in Fig. 4, we observe a
strong clustering effect in the model-latent space according to the dominant
label. The problem of label deficiency at the edge suggests the application of
FL in unsupervised and self-supervised settings [21]. The finding that
clustering occurs even when labels are not used in training suggests that in
these applications, and in instances where the adversary does not have access
to the labels used for a supervised objective, client privacy might be leaked.

Figure 4: Clustering of autoencoder parameters trained with an unsupervised
objective on MNIST. The clustering phenomenon occurs even in the absence of
labels, suggesting that on simple datasets, it can be attributed in part to
the similarity of features in training samples with the same label.

###  IV-B Predicting Client Data Distributions

To determine the extent to which the adversarial strategy proposed in
Algorithm 1 might work in a practical federated setting, and to test the
viability of adding noise to gradients as a defense against this strategy, we
run experiments on the MNIST and CIFAR-10 datasets, injecting three different
scales of Gaussian and Laplacian noise into client gradients.

For our experiments, the training set of the meta-dataset contains client
model parameters trained on samples from the training set of the original
dataset. Meanwhile, the test set of the meta-dataset contains client model
parameters trained on samples from the test set of the original dataset. This
train/test split reflects an assumption that the proxy data accessible by the
adversary is different from the training data on IoT devices, but drawn from
the same distribution. To generate synthetic non-IID client label
distributions on which dummy clients are trained, we sample from a symmetric
Dirichlet distribution, where the concentration parameter $\alpha$ governs the
uniformity of the client label distribution. To capture a wide range of
uniformities, we construct the meta-dataset using an equal number of client
label distributions with each of the concentration parameters
$\alpha\in\\{0.1,1,10,100,1000\\}$. Example client label distributions are
shown in Fig. 5. For each experiment, the training split of the meta-dataset
contains 10,000 samples in total (2,000 per value of $\alpha$), and the test
split contains 2,000 samples in total (400 per value of $\alpha$). We use a
simplified model in which the adversary only performs the attack during the
first round of FedAvg. In other words, the adversary intercepts the initial
global model $w_{0}$ and the client models $w_{1}^{k}$. In constructing the
meta-dataset, we reduce the client model parameters to 10 dimensions via PCA.

Figure 5: Example client label distributions sampled from a Dirichlet
distribution, used to create synthetic non-IID data. We constructed a meta-
dataset by training clients which had sample distributions drawn using
concentration parameters $\alpha\in\\{0.1,1,10,100,1000\\}$, which provided a
wide range of uniformities.

For each level of noise, we first measure the accuracy of a global model
trained with FedAvg across 500 rounds, in which training samples are assigned
to clients uniformly at random. In our tests, we use $K=100$ clients, each
receiving 500 training samples. We use client fraction of $C=0.1$ and train
each client for $E=1$ local epoch before aggregation. For MNIST, the client
updates use a batch size of 32 and a learning rate of $10^{-4}$; for CIFAR-10,
the updates use a batch size of 32 and a learning rate of $10^{-3}$. To reduce
the effects of random model initialization and client sampling, we perform
eight trials for each dataset and level of noise.

Then, we train a label distribution predictor on the training set of the meta-
dataset and assess its performance on the test set. The predictor uses an MLP
architecture consisting of an input layer with 10 units, eight hidden layers
with 1000 units each and ReLU activations, and an output layer with 10 units
and a softmax activation to produce a valid probability distribution over the
labels.

The results of the meta-dataset training are shown for the CIFAR-10 dataset in
Fig. 6; experiments on MNIST showed almost identical results. We find that
with $10^{-3}$-scale Gaussian and Laplacian noise, the global model is able to
learn, but the adversary can completely uncover the client data distributions.
Furthermore, even with $10^{-2}$-scale noise, when the model performance is
notably compromised, we found almost no benefit in preventing privacy leakage.
Consequently, our findings suggest that the adversarial strategy proposed in
Algorithm 1 is capable of predicting client label distributions with high
accuracy, and adding noise to client gradients is insufficient in protecting
against the adversarial model.

Figure 6: On the left is the test accuracy of the global model trained on the
CIFAR-10 dataset with FedAvg; on the right is the test loss of the label
distribution predictor trained on the meta-dataset. Adding Gaussian or
Laplacian noise appears to be insufficient in protecting against data
distribution leakage. In both cases, the scale of the noise required
($10^{-1}$) prevents the global model from learning.

##  V Conclusions & Future Directions

We have introduced an algorithm by which, in realistic settings, an adversary
can predict the label distributions of federated clients. Our approach lays
the groundwork for a general class of algorithms in which an adversary creates
a meta-dataset which maps client model parameters to some property of their
local data. Future work might investigate other approaches within this class
of algorithms, such as predicting the presence of a particular type of sample,
rather than a full label distribution. More research is also necessary to
determine the effectiveness of the proposed adversarial technique on more
complex datasets and with self-supervised training objectives. Finally, we
have shown that injecting noise into client gradients is insufficient to
protect privacy. Further work is needed to determine the effectiveness of
other defensive strategies in limiting these forms of attacks.

## Acknowledgment

We want to thank Bhupendra A. Raut, Seongha Park, Robert C. Jackson, Sean
Shahkarami, Yongho Kim, and Ian Foster for their insightful comments, ideas,
and support. This work was supported in part by U.S. Department of Energy,
Office of Science, under contract DE-AC02-06CHI1357, National Science
Foundation’s Mid-Scale Research Infrastructure grant, NSF-OAC-1935984, and
U.S. Department of Energy, Office of Science, under grant DOE-145-SE-
PRJ1009506. This research used resources of the Argonne Leadership Computing
Facility, which is a DOE Office of Science User Facility supported under
Contract DE-AC02-06CH11357.

## References

  * [1] N. Hudson, M. J. Hossain, M. Hosseinzadeh, H. Khamfroush, M. Rahnamay-Naeini, and N. Ghani, “A framework for edge intelligent smart distribution grids via federated learning,” in 2021 International Conference on Computer Communications and Networks (ICCCN), pp. 1–9, IEEE, 2021.
  * [2] N. Rieke, J. Hancox, W. Li, F. Milletarì, H. R. Roth, S. Albarqouni, S. Bakas, M. N. Galtier, B. A. Landman, K. Maier-Hein, S. Ourselin, M. Sheller, R. M. Summers, A. Trask, D. Xu, M. Baust, and M. J. Cardoso, “The future of digital health with federated learning,” npj Digital Medicine, vol. 3, sep 2020.
  * [3] L. Zhu, Z. Liu, and S. Han, “Deep leakage from gradients,” Advances in neural information processing systems, vol. 32, 2019.
  * [4] J. Geiping, H. Bauermeister, H. Dröge, and M. Moeller, “Inverting gradients-how easy is it to break privacy in federated learning?,” Advances in Neural Information Processing Systems, vol. 33, 2020.
  * [5] B. Zhao, K. R. Mopuri, and H. Bilen, “idlg: Improved deep leakage from gradients,” arXiv preprint arXiv:2001.02610, 2020.
  * [6] B. McMahan, E. Moore, D. Ramage, S. Hampson, and B. A. y Arcas, “Communication-efficient learning of deep networks from decentralized data,” in Artificial intelligence and statistics, PMLR, 2017.
  * [7] K. Bonawitz, H. Eichner, W. Grieskamp, D. Huba, A. Ingerman, V. Ivanov, C. Kiddon, J. Konečnỳ, S. Mazzocchi, B. McMahan, et al., “Towards federated learning at scale: System design,” Proceedings of machine learning and systems, vol. 1, pp. 374–388, 2019.
  * [8] T. Yang, G. Andrew, H. Eichner, H. Sun, W. Li, N. Kong, D. Ramage, and F. Beaufays, “Applied federated learning: Improving google keyboard query suggestions,” arXiv preprint arXiv:1812.02903, 2018.
  * [9] C. Fontaine and F. Galand, “A survey of homomorphic encryption for nonspecialists,” EURASIP Journal on Information Security, 2007.
  * [10] F. Wibawa, F. O. Catak, S. Sarp, M. Kuzlu, and U. Cali, “Homomorphic encryption and federated learning based privacy-preserving cnn training: Covid-19 detection use-case,” 2022.
  * [11] C. Dwork, F. McSherry, K. Nissim, and A. Smith, “Calibrating Noise to Sensitivity in Private Data Analysis,” in Theory of Cryptography (S. Halevi and T. Rabin, eds.), (Berlin, Heidelberg), pp. 265–284, Springer Berlin Heidelberg, 2006.
  * [12] M. Adnan, S. Kalra, J. C. Cresswell, G. W. Taylor, and H. R. Tizhoosh, “Federated learning and differential privacy for medical image analysis,” Scientific Reports, vol. 12, p. 1953, Feb. 2022.
  * [13] Z. Ji, Z. C. Lipton, and C. Elkan, “Differential privacy and machine learning: a survey and review,” 2014.
  * [14] G. Andrew, P. Kairouz, S. Oh, A. Oprea, H. B. McMahan, and V. Suriyakumar, “One-shot empirical privacy estimation for federated learning,” 2023.
  * [15] E. Hosseini and A. Khisti, “Secure aggregation in federated learning via multiparty homomorphic encryption,” in 2021 IEEE Globecom Workshops (GC Wkshps), pp. 1–6, 2021.
  * [16] W. Jin, Y. Yao, S. Han, C. Joe-Wong, S. Ravi, S. Avestimehr, and C. He, “Fedml-he: An efficient homomorphic-encryption-based privacy-preserving federated learning system,” 2023.
  * [17] C. Zhang, S. Li, J. Xia, W. Wang, F. Yan, and Y. Liu, “BatchCrypt: Efficient homomorphic encryption for Cross-Silo federated learning,” in 2020 USENIX Annual Technical Conference (USENIX ATC 20), pp. 493–506, USENIX Association, July 2020.
  * [18] K. Wei, J. Li, M. Ding, C. Ma, H. H. Yang, F. Farokhi, S. Jin, T. Q. Quek, and H. V. Poor, “Federated learning with differential privacy: Algorithms and performance analysis,” IEEE Transactions on Information Forensics and Security, vol. 15, pp. 3454–3469, 2020.
  * [19] A. El Ouadrhiri and A. Abdelhadi, “Differential privacy for deep and federated learning: A survey,” IEEE Access, vol. 10, 2022.
  * [20] H. Wang, Z. Kaplan, D. Niu, and B. Li, “Optimizing federated learning on non-IID data with reinforcement learning,” in IEEE INFOCOM 2020 - IEEE Conference on Computer Communications, pp. 1698–1707, 2020.
  * [21] C. He, Z. Yang, E. Mushtaq, S. Lee, M. Soltanolkotabi, and S. Avestimehr, “SSFL: Tackling label deficiency in federated learning via personalized self-supervision,” 2021.