# PatchGraph: In-hand tactile tracking with learned surface normals

Paloma Sodhi1,2, Michael Kaess2, Mustafa Mukadam1, and Stuart Anderson1
1Facebook AI Research, 2Carnegie Mellon University Code and supplementary
material can be found on https://psodhi.github.io/tactile-in-hand

###### Abstract

We address the problem of tracking 3D object poses from touch during in-hand
manipulations. Specifically, we look at tracking small objects using vision-
based tactile sensors that provide high-dimensional tactile image measurements
at the point of contact. While prior work has relied on a-priori information
about the object being localized, we remove this requirement. Our key insight
is that an object is composed of several local surface patches, each
informative enough to achieve reliable object tracking. Moreover, we can
recover the geometry of this local patch online by extracting local surface
normal information embedded in each tactile image. We propose a novel two-
stage approach. First, we learn a mapping from tactile images to surface
normals using an image translation network. Second, we use these surface
normals within a factor graph to both reconstruct a local patch map and use it
to infer 3D object poses. We demonstrate reliable object tracking for over
$100$ contact sequences across unique shapes with four objects in simulation
and two objects in the real-world.

##  I Introduction

We focus on the problem of tracking 3D object poses during in-hand
manipulations using tactile image measurements from vision-based tactile
sensors [1, 2]. Specifically, we look at tracking small objects without prior
geometric models. For instance, a dexterous robot operating in a real-world
household environment will need to manipulate novel household objects for
which CAD models may not be available. We address the question: Can an object
be tracked precisely enough for in-hand manipulation using only local
measurements of its geometry?

Prior work has looked at the object tracking problem primarily in the context
of planar pushing [3, 4, 5, 6]. However, the problem of 3D in-hand
manipulation poses additional challenges. Firstly, the motion is less
constrained such that different object motions can explain the same
measurements. Second, physics priors, such as quasi-static planar pushing
models, are less informative in the 3D case. Hence, prior work on in-hand
object tracking using tactile feedback has relied on a-priori information
about the object being localized, such as 3D global models or a database
generated by a simulator [7, 8, 9, 10].

Our key insight is that the tactile object tracking problem can be efficiently
decomposed in two ways. First, we can decompose an object into many smaller
local surface patches, which can be treated independently. Second, most of the
the information needed to infer the local surface patch geometry is already
embedded in corresponding tactile images.

We find that reliable tracking is achievable with only a local patch—a fused
map created from a sequence of key frame images within a continuous contact
episode. For instance, even though two objects can have very different global
geometries, they can contain very similar local patches that suffice for
tracking.

To both create the local patch and track motion relative to it, we must fuse
multiple tactile image measurements online while inferring the latent object
poses. We formulate this as an inference problem over a factor graph that
offers a flexible and efficient way to fuse such information while
incorporating other priors derived from physical and geometric constraints
[11, 12].

Figure 1: Tracking latent 3D object poses from a tactile image sequence during
in-hand manipulation. We solve this as an inference over a factor graph that
does not rely on having prior object models.

What makes for a good representation for a local patch? An idealized tactile
image captures the surface normals of the gel’s reflective layer, based on the
color and intensity of illumination at each pixel. Hence, it is natural to
learn a mapping from image to gel surface normals. While objects may have
varying global shapes, the learned surface normal mapping can generalize
across these different shapes, because the relationship between pixel
intensities and surface normals depends only on the sensor configuration and
the local contact geometry. Hence, we infer the surface normals at each pixel
in the tactile image, then integrate those normals to create a 3D model of the
visible section of the local patch.

We propose a novel two-stage approach for tracking objects in-hand without any
prior object model information (Fig. 1). First, we learn a mapping from
tactile images to surface normals using an image translation network. Second,
we use these surface normals within a factor graph to both reconstruct a local
patch map and use it to infer 3D object poses. Our key contributions are:

  1. 1.

A factor graph formulation for 3D in-hand tactile tracking that does not rely
on prior object models.

  2. 2.

A factor that works across different global object shapes by relying on local
patches generated from learned surface normals.

  3. 3.

Empirical evaluation on both simulation and real-world trials.

##  II Related Work

Factor graphs for localization Localization and mapping problems are
increasingly formulated as optimization objectives that leverage the inherent
sparsity of the problem to give tractable and more accurate solutions over
filtering approaches [13]. Factor graphs are a popular way for solving such
optimization objectives [11, 12, 3, 14, 15]. They offer a flexible way to fuse
multiple measurements while being computationally efficient to optimize.
Factors in the graph encode local potentials on variables such as observation
models between measurements and states as well as other priors such as physics
and geometry. We formulate our problem using a factor graph based framework in
this paper.

Vision-based touch sensing The advent of vision-based tactile sensors [1, 2,
16, 17] has enabled high-dimensional tactile image measurements that capture
the local deformation at the point of contact. Recent work has also looked at
creating accurate simulation models for such sensors [18, 19, 20]. These
sensors are being explored for various tactile manipulation applications. One
class of approaches use tactile images directly as local feedback to solve for
control actions on tasks such as object insertion [21], box packing [22], and
in-hand manipulations [23, 1]. However, such representations tend to overfit
to the particular task or require significant amount of data to generalize
across tasks. Our work focuses on extracting a state representation like the
global object pose that is easy to use and generalizes across different
downstream control and planning tasks.

Estimation from touch Prior work on estimating states from touch during
manipulation has included filtering methods [24, 25, 26, 27], learning-only
methods [28, 29], methods utilizing prior model information [8, 9, 10], and
graph-based optimization for planar pushing [6, 5, 4, 3]. In particular,
graph-based optimization offers benefits such as being more accurate than
filtering, an ability to incorporate analytic as well as learned models, and
can be solved in real-time making use of efficient, incremental solvers [30,
31, 32] in literature.

Of these different approaches, the work in [10, 9, 8, 7] is most closely
related in terms of the final objective of tracking 3D object poses during in-
hand manipulations.These, however, require prior object model information
either as offline models [7, 10, 9] or a database from a simulator [8]. In
contrast, we do not require a prior model of the object being tracked.
Instead, we build a local patch map on the fly for the current contact episode
and use that within a factor graph framework.

Figure 2: Overview of our factor graph formulation, where variables are object
poses and factors are constraints on poses. We first predict learned surface
normals from tactile images used to reconstruct a 3D point cloud. We then
create a set of factors. Image-to-image factors (blue) encode relative pose
between two point clouds. Image-to-patch factors (red) encode relative pose
between point cloud and a fused local patch. We also include global pose
priors as unary factors on end-effector (yellow) and first object pose (pink).

##  III Problem Formulation

We begin by formalizing the estimation problem as factor graph optimization. A
factor graph is a bipartite graph with two types of nodes: variables
$x\in\mathcal{X}$ and factors $\phi(\cdot):\mathcal{X}\rightarrow\mathbb{R}$.
Variable nodes are the latent states to be estimated, and factor nodes encode
potentials on these variables from observations, physics, or geometry (Fig.
2).

We solve for the maximum a posteriori (MAP) objective $\hat{x}$ by maximizing
product of all factor graph potentials, i.e.,

$\begin{split}\hat{x}&=\underset{{x}}{\operatorname{argmax}}\
\prod_{i}\phi_{i}({x_{i}})\\\ \end{split}$ (1)

To solve this inference objective efficiently, we assume $\phi_{t}(x)$ to be
Gaussian factors corrupted by zero-mean normally distributed noise. Under
Gaussian noise model assumptions, MAP inference is equivalent to a nonlinear
least-squares objective [12], i.e.,

$\begin{split}&\phi_{i}({x_{i}})\propto\exp\left\\{-\frac{1}{2}||f_{i}(x_{i};z_{i})||_{\Sigma_{i}}^{2}\right\\}\\\
\Rightarrow\ &\hat{x}=\underset{{x}}{\operatorname{argmin}}\
\frac{1}{2}\sum_{i}||f_{i}(x_{i};z_{i})||_{\Sigma_{i}}^{2}\\\ \end{split}$ (2)

For the in-hand object tracking problem, we define states in the graph to be
the 6-DOF object and end-effector poses at every time step $t=1\ldots T$, i.e.
$x_{t}=[o_{t}\ e_{t}]^{T}$, where $o_{t},e_{t}\in SE(3)$. Factors in the graph
include image-to-image factors $f_{im2im}(\cdot)$, image-to-patch factors
$f_{im2pc}(\cdot)$, velocity smoothness priors $f_{vel}(\cdot)$, end-effector
pose priors $f_{eff}(\cdot)$ and vision priors for re-localization at the
beginning of a contact episode $f_{vis}(\cdot)$. At every time step, new
variables and factors are added to the graph. Writing out Eq. LABEL:eq:eq3.2
for our problem,

$\begin{split}\scriptsize&\hat{x}_{1:T}=\underset{x_{1:T}}{\operatorname{argmin}}\sum_{t=1}^{T}\left\\{\right.||f_{im2im}(o_{t-1},o_{t},e_{t-1},e_{t})||^{2}_{\Sigma_{im2im}}+\\\
&||f_{im2pc}(o_{t},e_{t})||^{2}_{\Sigma_{im2pc}}+||f_{vel}(o_{t-2},o_{t-1},o_{t})||^{2}_{\Sigma_{vel}}+\\\
&||f_{eff}(e_{t})||^{2}_{\Sigma_{eff}}+||f_{vis}(o_{t})||^{2}_{\Sigma_{vis}}\left.\right\\}\end{split}$
(3)

Individual cost terms in Eq. LABEL:eq:eq3.3 are described in detail in Section
IV-C. Eq. LABEL:eq:eq3.3 is the optimization objective that we must solve for
every time step. Instead of resolving from scratch every time step, we make
use of efficient, incremental solvers [31] to solve this in real-time.

##  IV Approach

We present a two-stage approach: First, we learn a mapping from tactile images
to surface normals (Section IV-A) which can be integrated to create a 3D
reconstruction (Section IV-B). Second, we use these surface normals to create
a 3D local patch map online within a factor graph for inferring the latent 3D
object poses (Section IV-C).

###  IV-A Learning surface normals

Here we discuss how we learn to predict surface normals using tactile color
images from the Digit sensor [1]. Color images are generated by illuminating
the gel surface with three light sources (Fig. 3(a)), such that each light has
a unique color and direction. When pressing an object into the gel, the gel
surface conforms to the object. Under an idealized model of the sensor, the
color and intensity of light reaching the camera due to diffused reflection
from a point on the surface is directly related to the surface normal of the
gel at that point. Hence, the RGB intensity values of image pixel are expected
to contain significant information about the corresponding local surface
normals of the object.

To infer surface normal images from tactile images we train an image
translation network, pix2pix [33]. The pix2pix model is based on a generator-
discriminator network architecture that enables it to learn mappings from a
low amount of training data. It also enables us to learn a generalized
mapping, i.e. we test on objects unseen during training. To train the model,
we use a dataset of ground truth pairs of color and surface normal images
$z=\\{\mathcal{I}_{c},\mathcal{I}_{n}\\}$. We extract datasets in two ways.
For simulation dataset, we generate surface normals by adding a normal shader
to the Tacto simulator [18] based on Pyrender [34]. For the real dataset, we
follow a similar procedure as [2, 7] of using a ball bearing of known radius
whose ground truth normals can be synthesized. We manually annotate the
circular patches in RGB images, and synthesize ground truth normals in the
foreground annotations. We also augment our dataset with additional simulated
images to improve generalization.

###  IV-B Reconstruction from normals

Given an image with local surface normal information, we generate a 3D
representation of the corresponding surface geometry. This surface normal
image $N(x,y)$ is related to the gradient of a corresponding depth image
$z(x,y)$ as,

$\begin{split}\frac{\partial z}{\partial x}=\frac{n_{x}}{n_{z}},\hskip
20.0pt\frac{\partial z}{\partial y}=\frac{n_{y}}{n_{z}}\end{split}$ (4)

Given depth gradients $\\{\frac{\partial z}{\partial x},\frac{\partial
z}{\partial y}\\}$, we can recover the depth map $z(x,y)$ by integration using
a fast Poisson solver with discrete sine transform (DST) [35] as used in prior
work [2, 7]. For the boundary conditions, we use the mean distance to the
undisturbed gel surface. A consequence of this choice is that contact regions
crossing the edge of the image will not be correctly handled by the solver. We
filter out such images in our trials.

Finally, once we have computed the depth map, we inverse project it to obtain
a 3D point cloud. We use the OpenGL clip projection model [36] to map from
pixel to world coordinates, i.e. $\begin{bmatrix}x_{w}\ y_{w}\
z_{w}\end{bmatrix}=VP\begin{bmatrix}x_{pix}\ y_{pix}\ z_{depth}\end{bmatrix}$,
where $P$ is the projection matrix based on near, far plane values of the
projection frustum, and $V$ is the camera view matrix.

###  IV-C Factor graph optimization

Once we have the surface normal images, we integrate those along with other
priors as factors within a factor graph. The factor graph optimizer then
solves for the joint objective in Eq. LABEL:eq:eq3.3. We look at each of the
cost terms in Eq. LABEL:eq:eq3.3 in detail.

#### Image-to-image factors

In some cases, it suffices to look at consecutive tactile images and infer the
relative transformation between them. Color images from the tactile sensor at
the current and previous time step are converted into point clouds
$\\{\mathcal{P}_{t-1},\mathcal{P}_{t}\\}$ in their respective end-effector or
sensor frame. The two point clouds are registered against each other using a
point-to-plane iterative closest point (ICP) algorithm. The resultant relative
transformation is added as a binary factor between consecutive poses in the
graph. This is expressed as the $f_{im2im}(.)$ term in Eq. LABEL:eq:eq3.3,
i.e.,

$\begin{split}||f_{im2im}(o_{t\text{-}1},e_{t\text{-}1},o_{t},e_{t})||^{2}_{\Sigma_{im2im}}:=||{T}^{graph}_{t\text{-}1,t}\ominus
T^{reg}_{t\text{-}1,t}||^{2}_{\Sigma_{im2im}}\end{split}$ (5)

where
${T}^{graph}_{t\text{-}1,t}=(e_{t\text{-}1}o^{-1}_{t\text{-}1})\ominus(e_{t}o_{t}^{-1})$
are current relative estimates from the graph and
${T}^{reg}_{t\text{-}1,t}=\operatorname*{argmin}_{T}\sum_{i}||p^{(i)}_{t}-Tp^{(i)}_{t\text{-}1}||_{2}^{2}$
is the measured transformation from ICP registration. Both
${T}^{graph}_{t\text{-}1,t}$, ${T}^{reg}_{t\text{-}1,t}$ use gel center as
their origin.
$(p^{(i)}_{t-1},p^{(i)}_{t})\in(\mathcal{P}_{t-1},\mathcal{P}_{t})$ are pairs
of point correspondences in the two point clouds. $\ominus$ denotes difference
between two SE(3) manifold elements.

#### Image-to-patch factors

Image-to-image factors fail whenever the tactile image changes by a non-
trivial amount leading to large registration errors. This happens whenever the
object undergoes a larger transformation or moves in and out of the gel. To
stabilize tracking in these situations we introduce a local patch model by
fusing together multiple point clouds and register new tactile point cloud
data against this patch. The local patch is maintained by fusing together
images at specific key frames within the current contact episode, rather than
fusing all available frames. We choose these key frames at fixed intervals
given uniform motions but one can also select these based on a field-of-view
overlap threshold. The current cloud $\mathcal{P}_{t}$ is registered against
the local patch map cloud $\mathcal{P}_{map}$, and the relative transformation
is added as factors to the graph,

$\begin{split}||f_{im2pc}(o_{t},e_{t})||^{2}_{\Sigma_{im2pc}}:=||{T}^{graph}_{map,t}\ominus{T}^{reg}_{map,t}||^{2}_{\Sigma_{im2pc}}\end{split}$
(6)

where ${T}^{graph}_{map,t}=o_{t}e_{t}^{-1}$ are current estimates from the
graph and
${T}^{reg}_{map,t}=\operatorname*{argmin}_{T}\sum_{i}||p^{(i)}_{t}-Tp^{(i)}_{map}||_{2}^{2}$
is the measured transformation from ICP registration.
$(p^{(i)}_{t},p^{(i)}_{map})\in(\mathcal{P}_{t},\mathcal{P}_{map})$ are pairs
of point correspondences in the current point cloud and local patch map.

#### Constant velocity priors

We add a prior that assumes objects move at constant velocity, which has the
effect of smoothing tracked trajectories. This is a ternary factor between
triplets of object poses,

$\begin{split}||f_{vel}(o_{t\text{-}2},o_{t\text{-}1},o_{t})||^{2}_{\Sigma_{vel}}:=||o^{-1}_{t\text{-}2}o_{t\text{-}1}\ominus
o^{-1}_{t\text{-}1}o_{t}||^{2}_{\Sigma_{vel}}\end{split}$ (7)

#### End-effector priors

We model uncertainty about end-effector locations as unary priors on end-
effector variables,

$\begin{split}||f_{eff}(e_{t})||^{2}_{\Sigma_{eff}}:=||e_{t}\ominus\tilde{e}^{mc}_{t}||^{2}_{\Sigma_{eff}}\\\
\end{split}$ (8)

where, $\tilde{e}^{mc}_{t}={e}^{mc}_{t}\oplus\ \mathcal{N}(0,\Sigma_{eff})$
are poses from the motion capture system with added Gaussian noise. For a
robot end-effector, these pose measurements would instead come from the robot
kinematics.

#### Vision prior

We add a global vision pose prior for only the object pose at the start of an
episode,

$\begin{split}||f_{vis}(o_{t})||^{2}_{\Sigma_{vis}}:=||o_{t}\ominus\tilde{o}^{vis}_{t}||^{2}_{\Sigma_{vis}}\end{split}$
(9)

where, $\tilde{o}^{vis}_{t}$ are poses with added Gaussian noise
$\mathcal{N}(0,\Sigma_{vis})$. For re-localization during multi-contact
episodes, we can add such a factor at the start pose of every new contact
episode.

Figure 3: Real-world experimental setup showing: (a) a Sphere object mounted
on a workbench with a close-up of the Digit sensor, and (b) a contact episode
with the sensor moved around the object.

##  V Results and Evaluation

We evaluate our approach qualitatively and quantitatively on a number of
episodes where an object, unknown a priori, must be tracked from a sequence of
tactile measurements. We compare against a set of baselines on two fronts: a)
on surface normal predictions from images and b) on the final tracking error
of object poses. We use PyTorch [37] for training surface normal models, and
GTSAM C++ library [38] for factor graph optimization. We specifically use the
iSAM2 [31] solver for efficient, real-time optimization.

Figure 4: Training and evaluation of the pix2pix (unet) model that maps input
color images to predicted surface normal images. TABLE I: Learned surface
normal performance (validation loss)

Model type |  |  |  Object Shapes |  |
---|---|---|---|---|---
Sim Sphere |  Cube |  Real sphere |  Toy brick |  Toy human
pix2pix (unet) |  0.4e-3 |  0.5e-3 |  1.0e-3 |  1.1e-3 |  0.9e-3
pix2pix (resnet) |  0.6e-3 |  1.2e-3 |  1.6e-3 |  1.1e-3 |  1.0e-3
MLP 3-layer |  1.4e-3 |  0.5e-3 |  5.8e-3 |  5.1e-3 |  6.7e-3

Figure 5: Qualitative pose tracking performance on simulated trials. (a)
Object being tracked. (b) Local patch map reconstructed online. (c) Estimated
3D pose coordinates (light) against ground truth (dark). Rotations in orange
(x), purple (y), green (z). Translations in grey. Figure 6: Final tracking
error box plots on simulated trials: (top) rotation and (bottom) translation
errors in semi-log scale.

###  V-A Experimental setup

#### Simulator

We collected simulation data using the Tacto [18] simulator where one can load
the Digit sensor, an object and render high-resolution tactile image readings
in real-time. The simulator uses PyBullet [39] as the underlying physics
engine and Pyrender [34] as the back-end rendering engine for generating
images. We generated trials by using a position controller to move the object
on the sensor surface. We collect data for a diverse set of objects: Sphere,
Cube, Toy human and Toy brick.

#### Real-world

For real-world episode, we used a Digit [1] sensor to get tactile
measurements. We mounted the object on a workbench and mounted the Digit on a
movable plate, both of which were tracked by an OptiTrack motion capture
system to get ground truth poses. We collect real-world data for two objects:
a Sphere (ball-bearing) of 1/2” diameter and a Pyramid 1/2” tall with 1.75”
side length.

Figure 7: Qualitative pose tracking performance on real trials. (a) Object
being tracked. (b) Local patch map reconstructed online. (c) Estimated 3D pose
coordinates (light) against ground truth (dark). Rotations in orange (x),
purple (y), green (z). Translations in grey. Figure 8: Final tracking error
box plots on real trials: (top) rotation and (bottom) translation errors.

###  V-B Surface normal reconstructions

We first analyze the accuracy of learned surface normals from color images
collected by the Digit tactile sensor. For simulation episodes, we only train
on two simple objects: 50 images of Sphere and 50 images of Cube. We tested
the model on two different held-out shapes: Toy brick and Toy human. For real
episodes, we only had ground truth normals for the real Sphere. Hence, we
expanded the training dataset to include both the training data from
simulation and the real Sphere data. We tested the model on held-out real
Pyramid object. For baselines, we picked two different pix2pix architectures,
unet and resnet. We also trained a baseline 3-layer MLP 5-32-32-3 with tanh
activation on an L2 loss, mapping a single color pixel (r,g,b,x,y) to a
surface normal (nx,ny,nz), similar to prior work [7].

#### Qualitative reconstructions

Fig. 4 shows qualitative performance of the pix2pix (unet) model on both real
and simulated images. We can see the model generalizes to fairly different
shapes such as the Toy human and Toy brick, even though it is trained on
simple objects. Moreover, fine-tuning this model with very little real data,
i.e., 50 images of the real Sphere, enables it to generalize to unseen
geometries in the real-world such as different local patches of a Pyramid.

#### Quantitative model performance

Table I compares mean-squared pixel loss on the validation dataset for
different model choices. We see that pix2pix, for both unet and resnet
architecture, has a fairly low MSE loss and generalizes to unseen shapes such
as the Toy human and Toy brick. On the other hand, the per-pixel MLP baseline
incurs a high MSE loss. In general, the per-pixel MLP is likely to be
insufficient for non-ideal tactile sensors with effects such as self-shadowing
of the gel. Having convolutional layers, such as with the pix2pix
architectures, can address such confounding effects since these are typically
localized in the image, e.g. shadows are cast by nearby object features since
object depth is small relative to the image size.

###  V-C Factor graph optimization

We now look at the final task performance of tracking 3D object poses using
tactile image measurements. We compare 4 objectives, each of which use
different factors: ConstVel uses only a constant velocity prior, ImageToImage
uses Image-to-image factors, PatchGraph (ours) uses Image-to-patch factors and
GroundtruthPatch which uses a global object model. GroundtruthPatch assumes
the object is known apriori and hence represents the best a method can do.

#### Simulation tracking performance

Fig. 5 shows the qualitative tracking performance of PatchGraph for various
objects in simulated trials. For Cube, Toy brick and Toy human, PatchGraph is
able to reliably track rotations of the object. The local patch, constructed
from online estimates, appears to be consistent with the local object
geometry, thus explaining the good tracking performance. The object that was
most difficult to track was Toy brick, as evidenced by the distortions in the
patch, owing to jerkiness of the in-contact motions. We also note that
tracking rotations for a Sphere has high errors, which is expected since
rotation of a sphere are unobservable from tactile images.

Fig. 6 shows the quantitative rotation and translation tracking performance
against baselines on all objects, with $20$ distinct contact sequences per
object. Overall, PatchGraph has the lowest errors matching closely to that of
GroundtruthPatch that has the global object model. This lends credence to our
claim that a local patch suffices for reliable object tracking. ConstVel
understandably has the highest variance among any baselines. ImageToImage
fails to outperform PatchGraph on any of the datasets. Toy brick appears to be
the most challenging among datasets, where GroundtruthPatch has a clear
performance gap.

#### Real tracking performance

Fig. 7 shows qualitative tracking performance of PatchGraph for various
objects in real trials. Fig. 7(b) shows the local patches which appear
consistent with the local object geometry, and is a key piece to reliable
tracking. Fig. 7(c) shows good rotation tracking for Pyramid and translation
tracking for Sphere.

Fig. 8 shows quantitative rotation and translation tracking performance
against baselines on all objects, with $10$ distinct contact sequences per
object. The main observation is that ImageToImage performs much worse than
ImageToPatch, particularly in translation errors for Sphere. This is primarily
because point clouds generated from individual images are not as geometrically
discriminative as the fused local patch, causing ImageToImage factors alone to
diverge quickly. PatchGraph keeps translation errors under 4mm and rotation
errors under 0.2rad, which looks promising for use in dexterous object
manipulations.

##  VI Conclusion

We presented a factor graph-based approach for tracking 3D object poses from
tactile image sequences during in-hand manipulations. We showed reliable
tracking on $4$ simulated objects and $2$ real objects without relying on any
a priori object information. We achieved this by exploiting two decompositions
of the tracking problem. First, that a complex object can be treated as a
composition of many local patches each of which can be mapped and tracked
largely independently. Second, surface normal information is highly localized
within a tactile image and independent of the global object shape.

A primary limitation that can cause tracking failures is when the local patch
map is not sufficiently discriminative geometrically, e.g. flat or featureless
patches, or when the patch motions are degenerate in the observed image space,
e.g. rotations of a spherical object. As future work, it would be interesting
to explore solutions that take into account geometric degeneracies [40, 41] as
well as approaches that can detect slip and shear [42] to disambiguate motion
degeneracies. Another interesting future direction would be to complement the
tracker with a global first pose re-localization that is able to generalize
across objects, e.g. using visual images to predict a contact location
likelihood.

## References

  * [1] M. Lambeta, P.-W. Chou, S. Tian, B. Yang, B. Maloon, V. R. Most, D. Stroud, R. Santos, A. Byagowi, G. Kammerer, et al., “DIGIT: A novel design for a low-cost compact high-resolution tactile sensor with application to in-hand manipulation,” IEEE Robotics and Automation Letters, vol. 5, no. 3, pp. 3838–3845, 2020.
  * [2] W. Yuan, S. Dong, and E. H. Adelson, “GelSight: High-resolution robot tactile sensors for estimating geometry and force,” Sensors, vol. 17, no. 12, p. 2762, 2017.
  * [3] P. Sodhi, M. Kaess, M. Mukadam, and S. Anderson, “Learning tactile models for factor graph-based estimation,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), 2021.
  * [4] S. Suresh, M. Bauza, K.-T. Yu, J. G. Mangelson, A. Rodriguez, and M. Kaess, “Tactile slam: Real-time inference of shape and pose from planar pushing,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), 2021.
  * [5] A. S. Lambert, M. Mukadam, B. Sundaralingam, N. Ratliff, B. Boots, and D. Fox, “Joint inference of kinematic and force trajectories with visuo-tactile sensing,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), pp. 3165–3171, 2019.
  * [6] K.-T. Yu and A. Rodriguez, “Realtime state estimation with tactile and visual sensing. application to planar manipulation,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), pp. 7778–7785, 2018.
  * [7] S. Wang, Y. She, B. Romero, and E. Adelson, “Gelsight wedge: Measuring high-resolution 3d contact geometry with a compact robot finger,” arXiv preprint arXiv:2106.08851, 2021.
  * [8] M. Bauza, E. Valls, B. Lim, T. Sechopoulos, and A. Rodriguez, “Tactile object pose estimation from the first touch with geometric contact rendering,” Conference on Robot Learning (CoRL), 2020.
  * [9] J. Liang, A. Handa, K. Van Wyk, V. Makoviychuk, O. Kroemer, and D. Fox, “In-hand object pose tracking via contact feedback and gpu-accelerated robotic simulation,” arXiv preprint arXiv:2002.12160, 2020.
  * [10] M. Bauza, O. Canal, and A. Rodriguez, “Tactile mapping and localization from high-resolution tactile imprints,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), pp. 3811–3817, 2019.
  * [11] F. Dellaert, “Factor graphs: Exploiting structure in robotics,” Annual Review of Control, Robotics, and Autonomous Systems, vol. 4, 2020.
  * [12] F. Dellaert and M. Kaess, “Factor graphs for robot perception,” Foundations and Trends in Robotics, vol. 6, no. 1-2, pp. 1–139, 2017.
  * [13] C. Cadena, L. Carlone, H. Carrillo, Y. Latif, D. Scaramuzza, J. Neira, I. Reid, and J. J. Leonard, “Past, present, and future of simultaneous localization and mapping: Toward the robust-perception age,” IEEE Trans. Robotics, vol. 32, no. 6, pp. 1309–1332, 2016.
  * [14] J. Czarnowski, T. Laidlow, R. Clark, and A. J. Davison, “DeepFactors: Real-time probabilistic dense monocular SLAM,” IEEE Robotics and Automation Letters, vol. 5, no. 2, pp. 721–728, 2020.
  * [15] R. Hartley, M. G. Jadidi, L. Gan, J.-K. Huang, J. W. Grizzle, and R. M. Eustice, “Hybrid contact preintegration for visual-inertial-contact state estimation using factor graphs,” in IEEE/RSJ Intl. Conf. on Intelligent Robots and Systems (IROS), pp. 3783–3790, 2018.
  * [16] E. Donlon, S. Dong, M. Liu, J. Li, E. Adelson, and A. Rodriguez, “Gelslim: A high-resolution, compact, robust, and calibrated tactile-sensing finger,” in IEEE/RSJ Intl. Conf. on Intelligent Robots and Systems (IROS), 2018.
  * [17] A. Yamaguchi and C. G. Atkeson, “Combining finger vision and optical tactile sensing: Reducing and handling errors while cutting vegetables,” in IEEE Intl. Conf. on Humanoid Robots (Humanoids), pp. 1045–1051, 2016.
  * [18] S. Wang, M. Lambeta, P.-W. Chou, and R. Calandra, “Tacto: A fast, flexible and open-source simulator for high-resolution vision-based tactile sensors,” arXiv preprint arXiv:2012.08456, 2020.
  * [19] Z. Si and W. Yuan, “Taxim: An example-based simulation model for gelsight tactile sensors,” arXiv preprint arXiv:2109.04027, 2021.
  * [20] A. Agarwal, T. Man, and W. Yuan, “Simulation of vision-based tactile sensors using physics based rendering,” arXiv preprint arXiv:2012.13184, 2020.
  * [21] S. Dong, D. K. Jha, D. Romeres, S. Kim, D. Nikovski, and A. Rodriguez, “Tactile-rl for insertion: Generalization to objects of unknown geometry,” arXiv preprint arXiv:2104.01167, 2021.
  * [22] S. Dong and A. Rodriguez, “Tactile-based insertion for dense box-packing,” in IEEE/RSJ Intl. Conf. on Intelligent Robots and Systems (IROS), 2019.
  * [23] Y. She, S. Wang, S. Dong, N. Sunil, A. Rodriguez, and E. Adelson, “Cable manipulation with a tactile-reactive gripper,” in Robotics: Science and Systems (RSS), 2020.
  * [24] G. Izatt, G. Mirano, E. Adelson, and R. Tedrake, “Tracking objects with point clouds from vision and touch,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), pp. 4000–4007, 2017.
  * [25] B. Saund, S. Chen, and R. Simmons, “Touch based localization of parts for high precision manufacturing,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), pp. 378 – 385, 2017.
  * [26] M. C. Koval, N. S. Pollard, and S. S. Srinivasa, “Pose estimation for planar contact manipulation with manifold particle filters,” Intl. J. of Robotics Research, vol. 34, no. 7, pp. 922–945, 2015.
  * [27] Z. Pezzementi, C. Reyda, and G. D. Hager, “Object mapping, recognition, and localization from tactile geometry,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), pp. 5942–5948, 2011.
  * [28] B. Sundaralingam, A. S. Lambert, A. Handa, B. Boots, T. Hermans, S. Birchfield, N. Ratliff, and D. Fox, “Robust learning of tactile force estimation through robot interaction,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), 2019.
  * [29] R. Li, R. Platt, W. Yuan, A. ten Pas, N. Roscup, M. A. Srinivasan, and E. Adelson, “Localization and manipulation of small parts using gelsight tactile sensing,” in IEEE/RSJ Intl. Conf. on Intelligent Robots and Systems (IROS), pp. 3988–3993, 2014.
  * [30] M. Kaess, A. Ranganathan, and F. Dellaert, “iSAM: Incremental smoothing and mapping,” IEEE Trans. Robotics, vol. 24, no. 6, pp. 1365–1378, 2008\.
  * [31] M. Kaess, H. Johannsson, R. Roberts, V. Ila, J. Leonard, and F. Dellaert, “iSAM2: Incremental smoothing and mapping using the Bayes tree,” Intl. J. of Robotics Research, vol. 31, pp. 216–235, Feb. 2012.
  * [32] P. Sodhi, S. Choudhury, J. G. Mangelson, and M. Kaess, “ICS: Incremental constrained smoothing for state estimation,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), 2020.
  * [33] P. Isola, J.-Y. Zhu, T. Zhou, and A. A. Efros, “Image-to-image translation with conditional adversarial networks,” in Proc. IEEE Int. Conf. Computer Vision and Pattern Recognition, 2017.
  * [34] M. Matl, “Pyrender,” 2019.  https://github.com/mmatl/pyrender.
  * [35] J. Doerner, “Fast poisson reconstruction in python.”  https://gist.github.com/jackdoerner/b9b5e62a4c3893c76e4c.
  * [36] http://www.songho.ca/opengl/gl_projectionmatrix.html.
  * [37] A. Paszke, S. Gross, F. Massa, A. Lerer, J. Bradbury, G. Chanan, T. Killeen, Z. Lin, N. Gimelshein, L. Antiga, A. Desmaison, A. Kopf, E. Yang, Z. DeVito, M. Raison, A. Tejani, S. Chilamkurthy, B. Steiner, L. Fang, J. Bai, and S. Chintala, “PyTorch: An imperative style, high-performance deep learning library,” in Advances in Neural Information Processing Systems 24 (NeurIPS), pp. 8024–8035, 2019.
  * [38] F. Dellaert, “Factor graphs and GTSAM: A hands-on introduction,” tech. rep., Georgia Institute of Technology, 2012.
  * [39] E. Coumans and Y. Bai, “Pybullet, a python module for physics simulation for games, robotics and machine learning,” 2016.  https://pybullet.org/.
  * [40] J. Zhang, M. Kaess, and S. Singh, “On degeneracy of optimization-based state estimation problems,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), pp. 809–816, 2016.
  * [41] E. Westman and M. Kaess, “Degeneracy-aware imaging sonar simultaneous localization and mapping,” Journal of Oceanic Engineering, 2019.
  * [42] W. Yuan, R. Li, M. A. Srinivasan, and E. H. Adelson, “Measurement of shear and slip with a gelsight tactile sensor,” in IEEE Intl. Conf. on Robotics and Automation (ICRA), 2015.


# Technical Report (1)

EPSRC Early Career Fellowship, EP/T00729X/1.

##  1 Summary

Spectral Toolkit of Algorithms for Graphs (STAG) is an open-source C++ and
Python library of efficient spectral algorithms for graphs. Our objective is
to implement advanced graph algorithms developed through algorithmic spectral
graph theory, while making it practical to end users. This series of technical
reports is to document our progress on STAG, including implementation details,
engineering considerations, and the data sets against which our implementation
is tested. The report is structured as follows:

  * •

Section 2 describes the local clustering algorithm, which is the main update
in this STAG release. The discussion is at a high level such that domain
knowledge beyond basic algorithms is not needed.

  * •

Section 3 provides a user guide to the essential features of STAG which allow
a user to apply local clustering.

  * •

Section 4 includes experiments and demonstrations of the functionality of
STAG.

  * •

Finally, Section 5 discusses several technical details; these include our
choice of implemented algorithms, the default setup of parameters, and other
technical choices. We leave these details to the final section, as it’s not
necessary for the reader to understand this when using STAG.

###  1.1 Implemented Algorithms

STAG 1.2 provides an implementation of the following key algorithms.

#### Local Graph Clustering.

Given a large graph and some starting vertex $v$ in the graph, the goal of
local graph clustering is to find some cluster containing $v$. Moreover, the
running time of the algorithm should depend only on the size of the returned
cluster and should be independent of the total size of the graph [11].

STAG provides the first open-source local clustering algorithm which does not
require the entire graph to be loaded into memory. This allows users to apply
local clustering on massive graphs stored on disk or even in a cloud database,
such as Neo4j111https://neo4j.com/. Section 4 demonstrates these applications.

#### Spectral Clustering.

One of the most fundamental algorithms from Spectral Graph Theory is the
spectral clustering algorithm [8, 10, 13]. Spectral clustering is “global” in
the sense that it returns a partition of the entire vertex set of the graph.

#### Generating Graphs from Random Models.

The Stochastic Block Model (SBM) and Erdős-Rényi model are popular random
graph models which are frequently used to evaluate and analyse graph
algorithms. STAG provides several convenient methods to generate graphs from
these models.

##  2 Local Graph Clustering

Graph clustering algorithms are designed to partition an input graph into two
or more clusters. As a basic technique in data science and machine learning,
graph clustering has many applications in numerous areas of computer science
and beyond. Most graph clustering algorithms need to read an entire input
graph for the clustering task, which is computationally expensive if the graph
is massive. If one is interested only in some “local” cluster information,
then local graph clustering provides a more efficient method.

$S$$v$ Figure 1: Given a massive graph containing some small cluster $S$, a
local clustering algorithm takes as input a vertex $v\in S$ and returns an
approximation of $S$ without exploring the whole graph.

Typically, the objective of local graph clustering is to find some highly-
connected vertex set (cluster) in an input graph. Let’s assume that $S$ is a
highly-connected vertex set of an underlying undirected graph $G$, i.e., $S$
forms a cluster. Then, a local clustering algorithm is given some vertex $v\in
S$ as input, and returns some set $S^{\prime}$ such that $S^{\prime}$ is a
reasonable approximation of the target set $S$. Moreover, the running time of
the algorithm is proportional to the size of $S$ and independent of the size
of $G$. In applications, local clustering can be viewed as a search for
related objects: given a query vertex, a local clustering algorithm returns a
set of closely related vertices. Figure 1 illustrates local graph clustering.

Andersen, Chung, and Lang [1] introduced a key local clustering algorithm
which we refer to as the ACL algorithm. At a high level, the algorithm finds a
local cluster by analysing the behaviour of random walks on the graph,
beginning at the starting vertex. The ACL algorithm has proved extremely
useful and has inspired extensive further research and applications [2, 7, 12,
14]. The local_cluster method of STAG provides an implementation of the ACL
local clustering algorithm.

Existing open-source local clustering methods require that the entire graph is
loaded into RAM in order to apply the local algorithm. In this sense, they are
not truly “local” since they cannot be applied to graphs larger than the
available memory and the total running time depends on the size of the graph.
STAG provides the first open-source local clustering algorithm which can be
applied to massive graphs without loading them into RAM. Moreover, the
provided interface is simple and the algorithm can be applied to graphs stored
in memory, on disk, or in a Neo4j database.

##  3 User’s Guide to Local Graph Clustering with STAG

This section provides a guide to the essential features of STAG which allow a
user to apply local clustering. Section 3.1 describes how to install the STAG
C++ and STAG Python libraries. Then, Section 3.2 introduces the graph file
formats supported by STAG and demonstrates the methods for reading and writing
graphs to disk. Finally, Section 3.3 explains the graph classes provided by
STAG and Section 3.4 documents the local_cluster method for local clustering.

Although most of the examples in this section use C++, the functionality of
STAG C++ is also available in STAG Python. Appendix B includes example code
demonstrating how to perform local clustering with STAG Python. The full
documentation of STAG C++ and STAG Python is available on the STAG library
website.

###  3.1 Installation of STAG

This section describes how to install STAG for use with C++ and Python.

#### Installing STAG for C++.

STAG is built on the Eigen and Spectra C++ libraries, and these must be
installed before STAG. For information on installing Eigen and Spectra, please
refer to their documentation. For convenience, Appendix A provides a bash
script which, at the time of writing, will install Eigen and Spectra on a
standard Linux system. Then, the latest version of STAG should be downloaded
from

https://github.com/staglibrary/stag/releases.

After downloading and extracting the source code, STAG can be compiled and
installed with cmake.

⬇

1 mkdir build_dir

2 cd build_dir

3 cmake ..

4 sudo make install

Once STAG has been installed, it is available for use in C++ projects built
with the cmake build tools. The following cmake code will link a C++ project
with STAG.

⬇

1 find_package(stag REQUIRED)

2 include_directories(${STAG_INCLUDE_DIRS})

3 target_link_libraries(YOUR_PROJECT stag)

An example STAG project demonstrating the full cmake configuration is
available at

https://github.com/staglibrary/example-stag-project.

#### Installing STAG for Python.

STAG Python can be installed from the Python Package Index with the pip tool.

⬇

1 python -m pip install stag

Then, the modules of STAG can be directly imported into any Python script.

###  3.2 File Formats

STAG supports two simple file formats for storing graphs on disk: EdgeList and
AdjacencyList. Many graph datasets are provided in EdgeList format [6], and
will work directly with STAG.

#### EdgeList File Format.

In an EdgeList file, each line corresponds to one edge in the graph. A line
consists of two integer node IDs and an optional edge weight, all separated
with spaces. Here is an example of a simple EdgeList file.

⬇

1 # This is a comment

2 0 1 0.5

3 1 2 1

4 2 0 3

In this example, line 2 defines an edge between nodes $0$ and $1$ with weight
$0.5$.

#### AdjacencyList File Format.

In an AdjacencyList file, each line corresponds to one node in the graph. A
line consists of the node ID, followed by a list of adjacent nodes. The node
IDs at the beginning of each line must be sorted in increasing order. Here is
an example of a simple AdjacencyList file.

⬇

1 # This is a comment

2 0: 1 2

3 1: 0 3 2

4 2: 0 1

5 3: 1

In this example, node $1$ has edges to nodes $0$, $2$, and $3$.

#### Working with Files.

STAG provides several methods for reading, writing, and converting between
EdgeList and AdjacencyList files, as demonstrated in the following example.

⬇

1 #include <stag/graphio.h>

2 ...

3 // Read an AdjacencyList graph

4 std::string filename = "mygraph.adjacencylist";

5 stag::Graph myGraph = stag::load_adjacencylist(filename);

6

7 // Save an EdgeList graph

8 std::string new_filename = "mygraph.edgelist";

9 stag::save_edgelist(myGraph, new_filename);

10

11 // Convert an AdjacencyList to EdgeList directly

12 stag::adjacencylist_to_edgelist(filename, new_filename);

13 ...

###  3.3 Graph Classes

STAG provides several graph classes which can be applied for a wide variety of
applications. Figure 2 summarises the available graph classes.

$S$$v$DiskRAMDatabase$\Longleftrightarrow$$\Longleftrightarrow$$\Longleftrightarrow$
Figure 2: The graph classes provided by STAG. Every class inherits from the
abstract LocalGraph class, which provides the local clustering method. The key
difference between the different classes is the location of the graph data.
The AdjacencyListLocalGraph class reads node adjacency data from disk, the
Graph class stores the entire graph in RAM, and the Neo4jGraph class queries
node adjacency data from a Neo4j database.

#### The LocalGraph Class.

STAG provides an abstract stag::LocalGraph class which defines the data
structure necessary to apply local clustering. The only required method on the
data structure is neighbors(v) which returns a list of the neighbors of node
$v$. Every graph class provided by STAG inherits from LocalGraph.

#### The Graph Class.

The stag::Graph class is the basic graph object within the STAG library. The
class stores the adjacency matrix of the graph in memory as a sparse matrix.

#### The AdjacencyListLocalGraph Class.

The stag::AdjacencyListLocalGraph class provides an implementation of the
stag::LocalGraph interface for a graph stored on disk as an AdjacencyList. The
graph is loaded into memory in a local way only. This allows for local
algorithms to be executed on very large graphs stored on disk without loading
the whole graph into memory. The following example demonstrates how to create
an AdjacencyListLocalGraph with STAG C++.

⬇

1 #include <stag/graph.h>

2 ...

3 // Create an AdjacencyListLocalGraph

4 std::string filename = "mygraph.adjacencylist";

5 stag::AdjacencyListLocalGraph myGraph(filename);

6

7 // Get the neighbours of node 0

8 std::vector<long long> neighbors = myGraph.neighbors_unweighted(0);

9 ...

#### The Neo4jGraph Class.

STAG Python additionally provides the Neo4jGraph class, which provides an
implementation of the LocalGraph interface for a graph stored in a Neo4j
database. The following example shows how to create the Neo4jGraph using the
database connection information.

⬇

1 import stag.neo4j

2

3 # Connect to the database

4 uri = "<Database URI>"

5 username = "neo4j"

6 password = "<password>"

7 my_graph = stag.neo4j.Neo4jGraph(uri, username, password)

8

9 # Print the neighbors of node 0

10 print(my_graph.neighbors_unweighted(0))

The Neo4jGraph class provides additional methods for querying the properties
of the nodes in the database. The details are available in the full STAG
documentation.

###  3.4 Local Clustering

STAG provides the following local_cluster method.

⬇

1 std::vector<long long> stag::local_cluster(stag::LocalGraph* graph,

2 long long seed_vertex,

3 double target_volume)

Given a graph and a starting vertex, the local_cluster method finds a cluster
close to the starting vertex. The running time of the algorithm is
proportional to the size of the returned cluster and independent of the size
of the entire graph. The parameters of the method are described as follows:

  * •

graph \- a LocalGraph object. This could be a Graph, an
AdjacencyListLocalGraph, or a Neo4jGraph.

  * •

seed_vertex \- the starting vertex in the graph.

  * •

target_volume \- an estimate of the volume of the target cluster. This
parameter does _not_ impose a hard constraint on the algorithm and so an
approximate volume is sufficient.

When working with very large graphs, it is recommended to use the
AdjacencyListLocalGraph object for local clustering in order to avoid the
overhead of reading the entire graph into memory. Section 4.1 demonstrates the
advantage of using the AdjacencyListLocalGraph for local clustering. The
following code demonstrates a complete program which uses STAG C++ to find a
local cluster in a graph stored in an AdjacencyList file on disk.

⬇

1 #include <iostream>

2 #include <stag/graph.h>

3 #include <stag/cluster.h>

4

5 int main() {

6 // Create the graph backed by a large file on disk

7 std::string filename = "mygraph.adjacencylist";

8 stag::AdjacencyListLocalGraph mygraph(filename);

9

10 // Perform local clustering

11 int start_vertex = 1;

12 double target_volume = 100;

13 auto cluster = stag::local_cluster(&mygraph, start_vertex, target_volume);

14

15 // Print the returned cluster

16 for (auto v : cluster) std::cout << v << ", ";

17 std::cout << std::endl;

18 }

##  4 Showcase studies

STAG makes local clustering straightforward for a variety of applications, and
this section presents some examples of local clustering with STAG. The code
used to produce all experimental results is available at

https://github.com/staglibrary/local-clustering-case-study.

All experiments are performed on an HP ZBook laptop with an 11th Gen Intel(R)
Core(TM) i7-11800H @ 2.30GHz processor and 32 GB RAM.

###  4.1 Example 1

The advantage of local clustering over other clustering algorithms is that the
running time of local clustering is proportional to the size of the returned
cluster and independent of the total size of the graph. If we first load the
entire graph into memory before applying local clustering, then we lose the
advantage of the sub-linear running time. For this reason, STAG provides the
AdjacencyListLocalGraph class which provides local access to a graph stored on
disk without reading the entire graph. In this example, we compare the running
time of local clustering on an AdjacencyListLocalGraph object, which accesses
the graph locally on disk, and a Graph object, which loads the entire graph
into memory.

We generate graphs of various sizes from the stochastic block model as
follows. Given parameters $k$, $p$, and $q$, we create a graph with $k$
clusters $C_{1},\ldots C_{k}$, each containing 1,000 vertices. For every pair
of vertices $(u,v)\in V\times V$, we add the edge $(u,v)$ with probability $p$
if $u$ and $v$ are in the same cluster and with probability $q$ otherwise. We
always set $p=0.01$ and $q=0.001/k$. This ensures that the conductance of the
constructed clusters is always close to $0.1$.

We perform local clustering on the constructed graphs for a random starting
node and target volume 20,000, and compare the following two methods:

  * •

In memory: the entire graph is loaded into memory as a Graph object before
applying local clustering.

  * •

On disk: the graph is read locally from a file on disk with an
AdjacencyListLocalGraph object.

Figure 3 shows the running time of the local clustering algorithm for each
method across a range of graph sizes. These results demonstrate that for large
graphs, the overhead of reading the entire graph into memory dominates the
running time of the algorithm and reading the graph directly from disk is
significantly more efficient.

Figure 3: The comparison of the running time of local clustering on graphs in
memory and on disk. For massive graphs, reading the graph locally from disk is
significantly faster.

###  4.2 Example 2

In the second example, we demonstrate the applicability of local clustering
for finding sets of related nodes in a real-world graph. We use the wiki-
topcats dataset [14] which is a graph of Wikipedia hyperlinks constructed in
2011. The graph includes the pages in the top 100 Wikipedia categories and
includes 1,791,489 vertices and 28,511,807 edges. The dataset is available on
the SNAP datasets page [6] as an EdgeList file. We convert the EdgeList to an
AdjacencyList with the edgelist_to_adjacencylist method, and use the
AdjacencyListLocalGraph object for local clustering.

With a few lines of code, STAG allows us to create a “related pages” search
using local clustering. By providing a search page and setting the target
volume to be 100, the local clustering returns a set of pages which are
closely connected to the search page. Figure 4 shows some example of local
clustering results from the Wikipedia graph.

$S$$v$DiskRAMDatabase$\Longleftrightarrow$$\Longleftrightarrow$$\Longleftrightarrow$
Search page: Emacs Emacs Robert J. Chassell Zmacs Climacs Aquamacs Dunnet
(game) Bernard Greenberg Macro recorder GNU Manifesto TNT (instant messenger)
Agda (theorem prover)

$S$$v$DiskRAMDatabase$\Longleftrightarrow$$\Longleftrightarrow$$\Longleftrightarrow$
Search page: Stag Hound Stag Hound Over-canvassed sailing Sea Serpent
(clipper) Medium Clipper Donald McKay Extreme clipper Memnon (clipper) Flying
Cloud (clipper) Eleanor Creesy Weigh anchor Stowage

$S$$v$DiskRAMDatabase$\Longleftrightarrow$$\Longleftrightarrow$$\Longleftrightarrow$
Search page: Stagflation Stagflation Agflation Biflation Differential
accumulation Embedded liberalism Supply shock Shimshon Bichler Douglas Harper
Online savings account Jonathan Nitzan

Figure 4: Example of local clusters found in the Wikipedia dataset.

###  4.3 Example 3

In this example, we demonstrate local clustering on a Neo4j database in the
cloud. We first follow the Neo4j documentation to create a cloud database
using the AuraDB service [9]. We use the “Movies” example dataset provided by
Neo4j. Then, by creating a Neo4jGraph object with STAG Python, we are able to
search for related movies using local clustering. Listing 1 shows the complete
Python script used to perform this search, and Figure 5 shows some of the
search results.

⬇

1import stag.neo4j

2import stag.cluster

3

4# Initialise the graph object with the Neo4j database credentials

5uri = "<Database URI>"

6username = "neo4j"

7password = "<password>"

8g = stag.neo4j.Neo4jGraph(uri, user, password)

9

10# Perform local clustering

11seed_id = g.query_id("title", "The Matrix")

12cluster = stag.cluster.local_cluster(g, seed_id, 5)

13

14# Print the names of the returned movies.

15for node_id in cluster:

16 title = g.query_property(node_id, ’title’)

17 print(title)

Listing 1: Python code for local clustering with a Neo4j database.

$S$$v$DiskRAMDatabase$\Longleftrightarrow$$\Longleftrightarrow$$\Longleftrightarrow$
Search: The Matrix The Matrix The Matrix Reloaded The Matrix Revolutions V for
Vendetta Ninja Assassin

$S$$v$DiskRAMDatabase$\Longleftrightarrow$$\Longleftrightarrow$$\Longleftrightarrow$
Search: You’ve Got Mail You’ve Got Mail RescueDawn When Harry Met Sally As
Good as It Gets Sleepless in Seattle Joe Versus the Volcano

$S$$v$DiskRAMDatabase$\Longleftrightarrow$$\Longleftrightarrow$$\Longleftrightarrow$
Search: A Few Good Men A Few Good Men Hoffa As Good as It Gets Stand By Me

Figure 5: Example of local clusters found in the Neo4j movies dataset.

##  5 Technical Considerations

In this section we discuss the technical choices made in the design and
implementation of STAG. We discuss our choice of implemented algorithm, the
setting of the default parameters, and the implementation of the
AdjacencyListLocalGraph class for reading a graph locally from a file on disk.

###  5.1 Implemented Algorithm

Many algorithms have been proposed for local clustering, including those based
on PageRank [1, 2], the evolving set process [3, 4], and network flows [5]. We
chose to implement the algorithm based on PageRank presented by Andersen,
Chung, and Lang [1], and we refer to this as the ACL algorithm. We chose this
algorithm because it is relatively simple, easy to understand, and effective
in practice. Furthermore, the theoretical guarantees for the ACL algorithm are
optimal up to constant factors.222The original analysis by Andersen et al. [1]
has an extra factor of $\log(n)$ in the approximation guarantee. This factor
is not necessary and has been removed in later analysis using the same
technique [7, 12]. The ACL algorithm requires two parameters:

  * •

the $\alpha$ parameter controls the “teleport probability” of the personalised
PageRank; and

  * •

the $\epsilon$ parameter controls the approximation error of the approximate
PageRank calculation.

STAG provides the local_cluster_acl method which allows the user to specify
the parameters $\alpha$ and $\epsilon$ directly.

⬇

1 std::vector<long long> stag::local_cluster_acl(stag::LocalGraph* graph,

2 long long seed_vertex,

3 double alpha,

4 double epsilon)

For convenience, STAG also provides the local_cluster method which requires
only an estimate of the volume of the target cluster. Given a volume $\gamma$,
the local_cluster method uses the parameters $\alpha=1/2000$ and
$\epsilon=1/(20\gamma)$ for the ACL algorithm.

###  5.2 Reading Graphs Locally From Disk

A key feature of STAG is the AdjacencyListLocalGraph class which reads the
neighbourhood information of a graph in a local way from an AdjacencyList file
on disk. Since the data in an AdjacencyList file is sorted according to the
node ID, we can query the neighbors of any node in $O(\log(n))$ time by binary
search of the AdjacencyList file. As demonstrated in Section 4.1, this
additional logarithmic factor in the running time is much preferable to the
cost of reading the entire graph into memory when applying local algorithms to
massive graphs.

## References

  * [1] Reid Andersen, Fan Chung, and Kevin Lang.  Local graph partitioning using Pagerank vectors.  In 47th Annual IEEE Symposium on Foundations of Computer Science (FOCS’06), pages 475–486, 2006.
  * [2] Reid Andersen, Fan Chung, and Kevin Lang.  Local partitioning for directed graphs using Pagerank.  In 5th International Workshop on Algorithms and Models for the Web-Graph (WAW’07), pages 166–178, 2007.
  * [3] Reid Andersen, Shayan Oveis Gharan, Yuval Peres, and Luca Trevisan.  Almost optimal local graph clustering using evolving sets.  Journal of the ACM, 63(2):1–31, 2016.
  * [4] Reid Andersen and Yuval Peres.  Finding sparse cuts locally using evolving sets.  In 41st Annual ACM Symposium on Theory of Computing (STOC’09), pages 235–244, 2009.
  * [5] Kimon Fountoulakis, Di Wang, and Shenghao Yang.  $p$-Norm flow diffusion for local graph clustering.  In 37th International Conference on Machine Learning (ICML’20), page 3222–3232, 2020.
  * [6] Jure Leskovec and Andrej Krevl.  SNAP Datasets: Stanford large network dataset collection.  http://snap.stanford.edu/data, June 2014.
  * [7] Peter Macgregor and He Sun.  Local algorithms for finding densely connected clusters.  In 38th International Conference on Machine Learning (ICML ’21), pages 7268–7278, 2021.
  * [8] Peter Macgregor and He Sun.  A tighter analysis of spectral clustering, and beyond.  In 39th International Conference on Machine Learning (ICML ’22), pages 14717–14742, 2022.
  * [9] Neo4j.  Neo4j graph database.  http://neo4j.org/, 2023.
  * [10] Andrew Y Ng, Michael I Jordan, and Yair Weiss.  On spectral clustering: Analysis and an algorithm.  In 15th Advances in Neural Information Processing Systems (NeurIPS’01), pages 849–856, 2001.
  * [11] Daniel A. Spielman and Shang-Hua Teng.  A local clustering algorithm for massive graphs and its application to nearly linear time graph partitioning.  SIAM Journal on Computing, 42(1):1–26, 2013.
  * [12] Yuuki Takai, Atsushi Miyauchi, Masahiro Ikeda, and Yuichi Yoshida.  Hypergraph clustering based on Pagerank.  In 26th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (KDD’20), pages 1970–1978, 2020.
  * [13] Ulrike von Luxburg.  A tutorial on spectral clustering.  Statistics and Computing, 17(4):395–416, 2007.
  * [14] Hao Yin, Austin R Benson, Jure Leskovec, and David F Gleich.  Local higher-order graph clustering.  In 23rd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (KDD’17), pages 555–564, 2017.

##  Appendix A Installing STAG Dependencies

For convenience, we provide the following bash script for installing the STAG
C++ dependencies. At the time of writing, this will download and install the
Eigen and Spectra libraries.

⬇

1 # Create a directory to work in

2 mkdir libraries

3 cd libraries

4

5 # Install Eigen

6 wget https://gitlab.com/libeigen/eigen/-/archive/3.4.0/eigen-3.4.0.tar.gz

7 tar xzvf eigen-3.4.0.tar.gz

8 cd eigen-3.4.0

9 mkdir build_dir

10 cd build_dir

11 cmake ..

12 sudo make install

13 cd ../..

14

15 # Install Spectra

16 wget https://github.com/yixuan/spectra/archive/v1.0.1.tar.gz

17 tar xzvf v1.0.1.tar.gz

18 cd spectra-1.0.1

19 mkdir build_dir

20 cd build_dir

21 cmake ..

22 sudo make install

23 cd ../..

##  Appendix B User Guide Examples using STAG Python

This section includes example code omitted from Section 3 demonstrating how to
use STAG Python for local clustering.

###  B.1 Working with Files

The following example demonstrates how to read and write AdjacencyList and
EdgeList files with STAG Python.

⬇

1 import stag.graphio

2

3 # Read an AdjacencyList graph

4 filename = "mygraph.adjacencylist"

5 myGraph = stag.graphio.load_adjacencylist(filename)

6

7 # Save an EdgeList graph

8 new_filename = "mygraph.edgelist"

9 stag.graphio.save_edgelist(myGraph, new_filename)

10

11 # Convert an AdjacencyList to EdgeList directly

12 stag.graphio.adjacencylist_to_edgelist(filename, new_filename)

###  B.2 Graph Classes

The following example shows how to create an AdjacencyListLocalGraph object
with STAG Python.

⬇

1 import stag.graph

2

3 # Create the graph object

4 filename = "mygraph.adjacencylist"

5 my_graph = stag.graph.AdjacencyListLocalGraph(filename)

6

7 # Show the neighbors of node 0

8 print(my_graph.neighbors_unweighted(0))

###  B.3 Local Clustering

The following example gives a complete program for finding a local cluster in
a graph stored in an AdjacencyList file with STAG Python.

⬇

1 import stag.graph

2 import stag.cluster

3

4 # Create the graph object

5 filename = "mygraph.adjacencylist"

6 mygraph = stag.graph.AdjacencyListLocalGraph(filename)

7

8 # Find a local cluster

9 start_vertex = 1

10 target_volume = 100

11 cluster = stag.cluster.local_cluster(mygraph, start_vertex, target_volume)

12

13 # Display the result

14 print(cluster)