# Putting 3D Spatially Sparse Networks on a Diet

Junha Lee1 Christopher Choy2 Jaesik Park1
POSTECH1 NVIDIA2

###### Abstract

3D neural networks have become prevalent for many 3D vision tasks including
object detection, segmentation, registration, and various perception tasks for
3D inputs. However, due to the sparsity and irregularity of 3D data, custom 3D
operators or network designs have been the primary focus of 3D research, while
the size of networks or efficacy of parameters has been overlooked. In this
work, we perform the first comprehensive study on the weight sparsity of
spatially sparse 3D convolutional networks and propose a compact weight-sparse
and spatially sparse 3D convnet (WS3-ConvNet) for semantic segmentation and
instance segmentation. We employ various network pruning strategies to find
compact networks and show our WS3-ConvNet achieves minimal loss in performance
(2.15% drop) with orders-of-magnitude smaller number of parameters (1/100
compression rate). Finally, we systematically analyze the compression patterns
of WS3-ConvNet and show interesting emerging sparsity patterns common in our
compressed networks to further speed up inference.

(a) Input point cloud

(b) Reference Network,
37.85M Param., 71.57 mIoU(%)

(c) Our WS3-ConvNet,
0.396M Param., 69.42 mIoU(%)

Figure 1:  Visualization of semantic label prediction of reference neural
network (middle) and our WS3-ConvNet that is obtained by pruning 99% of the
weights from the reference neural network (right). While the pruned model has
100 times smaller parameters, the mIoU over the entire ScanNet [6] validation
split is 69.42%, only the 2.15% drops.

##  1 Introduction

Recent advances in 3D neural networks made various 3D vision applications such
as 3D shape classification [31, 32], semantic segmentation [4], object
detection [30, 10], reconstruction [25], 3D registration [5, 3], and many
other tasks accessible. These work focus on broadening the spectrum of the 3D
neural networks, but, over the years, the majority of research on 3D
perception focus more on introducing customized convolution kernels that
preserve the continuity of point cloud [31, 32, 38, 40, 37, 24, 47]. However,
these continuous convolutions sometimes show lower performance or slower
processing time on large scale datasets such as ScanNet [6] or autonomous
driving datasets [2, 34] due to the slow and memory intensive nearest neighbor
search in 3D. Specifically, to accelerate the nearest neighbor search on GPU,
they use pairwise distance computation that is expensive and requires
quadratic $\mathcal{O}(N^{2})$ memory footprint. This limits the resolution of
the point cloud and prevents wider adoption on large-scale 3D perception.

Discretized convolutions, on the other hand, sacrifice small quantization
error ($\sim$1cm) for scalability that can process millions of points with
faster inference using the spatially sparse representation and GPU hash tables
[9, 4, 10]. Thus, they show promising results on many large-scale indoor [1,
6] and outdoor perception tasks [2, 34]. These discrete convolutions are
effective in processing large-scale point clouds, yet these networks are still
bulky due to their 3D convolution kernels, preventing wide adoption on edge
devices, commodity servers, and low-power devices.

In the image domain, we observe significant advances in network pruning
techniques that compress a network by removing redundant parameters [19, 7, 8,
20, 26, 27] and network architecture design methods that reduce the memory and
computation cost of a network by introducing efficient architectures [15, 46,
35]. Compared with more mature research in 2D network compression, the network
pruning and compression of 3D convnets have rarely been studied except for
pruning spatially dense 3D convolution kernels for 2D videos [41].

In this work, to the best of our knowledge, we conduct the first comprehensive
study on weight pruning of spatially sparse 3D convolutional networks and
create Weight-Sparse and Spatially Sparse 3D convolutional networks
(WS3-ConvNet). We extensively employ the network pruning methods [11, 20, 7,
26] on two indoor datasets (ScanNet [6] and S3DIS [1]). For each dataset, we
prune networks for 1) semantic segmentation networks and 2) instance
segmentation networks creating four distinct tasks for network compression.
Specifically, we use six pruning methods with three backbone network
architectures to find compact 3D neural networks. In total, it requires 5800
hours of GPU computation on 30 GPUs for us to find the optimal pruning setting
for 3D network compression. The most effective networks we found show 1/100
compression rate with a negligible drop in performance on 3D perception tasks.
We name these networks WS3-ConvNet.

In addition, we propose the weight-sparse spatially sparse convolution
(WS3-conv) algorithm. Network compression results in smaller networks and
reduction in memory footprint, but faster inference requires dedicated
software. Our WS3-conv algorithm shows up to 34% reduction in inference time.
Lastly, we analyze compression patterns of WS3-ConvNet and propose additional
structural pruning on top of the most performant network to further speed up
inference.

We summarize our contributions as follows:

  * •

We apply network pruning methods on various 3D spatially sparse convolution
networks for the various 3D vision tasks that require 5800 hours of GPU time
and propose WS3-ConvNet.

  * •

We propose an efficient feed-forward algorithm for our weight-sparse spatially
sparse convolution, which makes use of the compressed weights for faster
inference.

  * •

Analysis of the pruned spatially sparse convnets reveal interesting properties
of 3D convnets. 1) Pruning a larger 3D convnet is more effective. 2) The
decoder part of a 3D U-net is less significant than the encoder. 3) the
pruning patterns reveal that a few spatial convolution kernels have no
significance.

  * •

We make use of these observations to apply structural pruning and further
speed up our WS3-ConvNet.

##  2 Related Work

Deep neural networks for 3D data. We limit the scope of our work to spatial 3D
data, although many refer to videos as 3D data, as spatial 3D data requires
special operators to handle sparsity and irregularity. MLP-based methods [31,
32, 47] directly process the continuous 3D coordinates of point clouds from 3D
sensors using a set of MLPs. Recent approaches [24, 38, 40, 37] define
customized convolution kernels that apply a series of local transformations on
each point. Sparse voxel-based approaches [9, 4], on the other hand,
discretize point clouds into discrete voxels and apply spatially sparse
convolution, which is mathematically equivalent to regular dense convolution.

Each group targets different aspects of point cloud perception and there are
pros and cons of each category. For instance, MLP-based methods [31, 32] are
the fastest and the simplest, but often fail to capture local context; others
preserve the continuity of points while requiring expensive nearest neighbor
search with quadratic memory complexity on GPU [32, 38, 47]. Voxel-based
methods introduce a small quantization error, but discretization provides
$\mathcal{O}(1)$ search speed and $\mathcal{O}(N)$ memory footprint using GPU
hash tables. Also, thanks to the mathematical equivalence to 2D convnets, we
could borrow many successful network architectures such as skip-connection
[13], batch normalization [18], and U-shaped network [33] from 2D CNNs to 3D
spatially sparse CNNs.

Network compression. Network pruning or network compression is a technique
that removes a portion of weights from a neural network to shrink the size of
a network and to reduce the amount of computation required for a forward pass.
There are three broad categories of network compression depending on when,
where, and how a network is pruned.

_1) Unstructured pruning_ removes network weights without a predefined
structure. LeCun  _et al_. [19] identifies the redundant parameters by
analyzing sensitivities of parameters using the second-order derivatives of
the objective function. Later work [12, 26] replaced the hessian with the
Taylor expansion to relax the computation burden of computing the second-order
derivatives. Han  _et al_. [11] showed that pruning weights with small
magnitude and fine-tuning the pruned network leads to a compact model without
incurring large loss in accuracy.

_2) Structured pruning_ removes a set of network weights with a specific
structure (e.g., prune channels, filters altogether) that could provide gains
in computation and inference speed without special hardware or software
support. Li  _et al_. [21] applied magnitude-based pruning to network filters.
Other work [22, 16] learn a scaling factor for each channel and prune an
entire channel with a lower magnitude scaling factor. Other approaches [23,
14] find channels that don’t contribute to the feature reconstruction.

_3) Pruning at initialization_ is a recent network pruning method based on the
Lottery Ticket Hypothesis (LTH) [7], which claims that there is a subnetwork
within a network that achieves optimal performance and one can find the
subnetwork by rewinding parameters that survive pruning to their initial
values. Frankle  _et al_. [8] use larger architectures by relaxing the
restriction of reverting the weights to initial values. Zhou  _et al_. [48]
show that the difference between the initial weights and fine-tuned weights
can be another pruning criterion. The proxies for determining lottery tickets
in a data-efficient way have been recently studied extensively as well. SNIP
[20] aims to find performant subnetworks with a few mini-batch iterations.
GraSP [39] and SynFlow [36] suggest that analyzing gradient-flow between
layers enables identifying lottery tickets with a small set of training data
or even without data.

Note that all of the pruning methods we discussed above mainly focus on 2D
perception, i.e., image classification. There is a compression work on dense
3D ConvNets [42], but it is for video understanding to find spatiotemporal
features not for spatial 3D data. Hence, to the best of our knowledge, we use
network compression methods for spatially sparse 3D neural networks for the
first time for 3D perception.

Among the various categories of pruning methods, we mainly use _unstructured
pruning_ since it is the most effective in finding highly compressed networks
that have order-of-magnitude smaller parameters, which is a critical factor
for 3D convnet compression. Also, it can be applied to generic network
architectures due to its simplicity, is simpler to implement, and can provide
insights into structural properties from the final pruning patterns.

##  3 Preliminary

Perception for spatial 3D data such as LiDAR point clouds or RGB-D scans
requires non-standard operators that are not provided in the most off-the-
shelf neural network libraries to handle sparsity and irregularity of 3D data
and to implement fast and efficient local aggregation functions. There are
many representations and local aggregation functions, but in this work we
focus on spatially sparse convnets. So, we cover the basics of a sparse tensor
and spatially sparse convolution in this section.

A sparse tensor is the most common data structure for high-dimensional sparse
data representation. It assumes that the majority of the data is 0 and saves
only the non-zero elements by representing them as their coordinates and
values. For high-dimensional tensors, the coordinate format (COO) is the most
easy-to-use and versatile representation and we follow the convention.
Formally, a $D+1$ dimensional spatial sparse tensor consists of $D$ spatial
dimensions and one dense feature dimension. Such sparse tensor with $N$ non-
zero elements can be expressed as follows:

$\mathcal{T}=(\mathbf{C},\mathbf{F}),$ (1)

where $\mathbf{C}\in\mathbb{Z}^{N\times D}$ is row-wise concatenation of
discretized coordinates, and $\mathbf{F}\in\mathbb{R}^{N\times N_{F}}$ is
corresponding feature vectors of size $N_{F}$. We will use $\mathbf{C}$ to
denote a set of non-zero coordinates and a coordinate matrix, but it will be
clear from the context.

Spatially sparse convolution is a generalized version of the conventional
dense convolution so that it can operate on a sparse tensor. Given an input
sparse tensor at $l$-th layer,
$\mathcal{T}^{l}=(\mathbf{C}^{l},\mathbf{F}^{l})$ where
$\mathbf{C}^{l}\in\mathbb{Z}^{N\times D}$ denotes the $D$-dimensional
coordinates, and $\mathbf{F}^{l}\in\mathbb{R}^{N\times N_{F}^{l}}$ denotes
corresponding features. The spatially spare convolution is defined as follows:

$\mathbf{f}_{\mathbf{c}}^{l+1}=\sum_{\mathbf{i}\in\mathcal{K}(\mathbf{c})}\mathbf{W}_{\mathbf{i}}^{l}\cdot\mathbf{f}_{\mathbf{c+i}}^{l},$
(2)

where $\mathbf{f}^{l}_{\mathbf{c}}\in\mathbf{F}^{l}$ is a feature vector at
coordinate $\mathbf{c}$, $\mathbf{W}\in\mathbb{R}^{K^{D}\times
N_{F}^{l+1}\times N_{F}^{l}}$ is the kernel weights with kernel size $K$,
$\mathbf{c}\in\mathbb{Z}^{D}$ is the current position of kernel center, and
$\mathcal{K}(\mathbf{c})$ denotes a list of offsets in $D$-dimensional
hypercube centered at $\mathbf{c}$. The convolution kernel,
$\mathcal{K}(\cdot)$, is applied only on non-zero elements as the rest of the
input is 0.

Neighbor search, where a convolution kernel finds non-zero neighbors within
the current receptive field around $\mathbf{c}$ ($\mathcal{K}(\mathbf{c})$),
is one of the areas the spatially sparse convolution shows its effectiveness
and scalability. It uses a GPU hash table to speed up the search with
$\mathcal{O}(1)$ complexity, which differentiates itself with other approaches
that use continuous points and $\mathcal{O}(N^{2})$ memory complexity for
neighbor search [31, 32, 38, 47, 24].

Methods | Neighbor Search |  Large-scale Inference |  Latency (Sec.)
---|---|---|---
Preparation | Inference
KPConv  _deform_ [38] | $\mathcal{O}(N\log N)$ | $\mathcal{O}(KM\log N)$ | Crop & Stitch | 105.15
PAConv [44] | $\mathcal{O}(1)$ | $\mathcal{O}(MN\log K)$ | Crop & Stitch | 28.13
PointTransformer [47] | $\mathcal{O}(1)$ | $\mathcal{O}(MN\log K)$ | Multi-shot | 18.07
gray!20Res16UNet34C | $\mathcal{O}(N)$ | $\mathcal{O}(M)$ | Single-shot | 0.07

Table 1: Run time comparison. We analyze the theoretical time complexity of
neighbor search algorithms and evaluate per-scene wall-time latency of each
network. We denote $N$ as # of dataset points, $M$ as # of query points, and
$K$ as # of neighbors to search. Both $M$ and $N$ are much larger than $K$ in
a large-scale point cloud. We use S3DIS [1] Area 5 test scene with a single
Geforce RTX 3090 GPU for the comparison.

Note that continuous convolution methods [38, 44, 47] require multiple crops
with a sliding window to handle large-scale scenes. Such sliding window makes
them a few orders of magnitude slower than hash-table based neighbor search of
spatially sparse convolution.

Spatial sparsity vs. weight sparsity. The 3D convolutions are bulkier than 2D
networks due to the extra spatial dimension. Note that in 2D, there are
$K^{2}$ weight matrices of size $N_{F}^{l+1}\times N_{F}^{l}$ where $K$ is the
kernel size, but in 3D, it requires $K^{3}$ weight matrices. Thus, to reduce
the number of parameters, we prune some of the weights from the convolution
filters by applying successful pruning approaches [19, 20, 7, 8, 12, 26, 11].
As a result, the pruned network has sparse weights, and we use them for
convolution. To revisit Eq. 2, we prune the weight matrix $W_{i}$ and convert
it to a sparse matrix. However, the convolution is applied only on spatially
sparse coordinates $\mathbf{C}$.

##  4 Spatially Sparse Network Compression

Compressing a neural network requires first training a neural network with
dense weights. We first discuss the target tasks and the training procedure in
Section 4.1. Next, we introduce pruning methods we adopt in Section 4.2 and
discuss the pruning procedure. In Section 4.3, we use the compressed spatially
sparse convnet or weight-sparse spatially sparse convnet (WS3-ConvNet) for
faster inference by introducing our weight-sparse convolution algorithm.

###  4.1 Training Spatially Sparse ConvNets

To achieve high performance and good network compression rate, neural network
has to be fine-tuned on a specific task. We define two most common 3D tasks
for our network compression that are practical and computationally demanding:
semantic segmentation and instance segmentation. For all experiments, we use
U-shaped fully convolutional networks with residual blocks as depicted in
Figure 2 (top). For the last layer, we use logit prediction per voxel for
semantic segmentation. For instance segmentation, we use logits with object
center offset prediction as two outputs. We use a pytorch spatially sparse
convolution implementation for all training [4].

We follow the standard training procedure in [4], but do not apply any
rotation averaging, or more advanced data augmentation techniques [28] to
compare our network fairly with other methods that use the same datasets. We
use cross-entropy loss for semantic segmentation:

$\mathcal{L}_{\mathrm{sem}}=-\frac{1}{N}\sum_{n=1}^{N}\sum_{c=1}^{C}y_{n,c}\log({p_{n,c}}),$
(3)

and for instance segmentation, we use the auxiliary loss for regressing the
object center offset as follows:

$\mathcal{L}_{\mathrm{ins}}=\mathcal{L}_{\mathrm{sem}}+\frac{1}{N}\sum_{n=1}^{N}\left(|\mathbf{c}_{n}-\mathbf{c}_{n}^{*}|-\frac{\mathbf{c}_{n}\cdot\mathbf{c}_{n}^{*}}{\|\mathbf{c}_{n}\|\|\mathbf{c}_{n}^{*}\|}\right).$
(4)

We train total three variants of networks with various depths and widths and
train them with the same training procedure. The detailed configuration of
each network and the training parameters can be found in Section 5.2 and the
supplementary material.

Input : $\mathcal{F}(\cdot,\theta_{\mathrm{init}})$: a pre-trained net. with
param. $\theta_{\mathrm{init}}$.

$i_{\mathrm{train}}$: number of training iterations.

$i_{\mathrm{prune}}$: number of fine-tuning iterations per step.

$N_{\mathrm{prune}}$: number of pruning iteration to repeat.

$p_{\mathrm{target}}$: target pruning ratio.

criteria$(\cdot)$: a score function for pruning.

Output: $\theta^{\prime}$: a pruned network parameter.

1 function IterativePruning

2 $\theta\leftarrow\theta_{\mathrm{init}}$

$m\leftarrow$ ones_like($\theta$) // init. pruning mask

3 $p\leftarrow 1-(1-p_{\mathrm{target}})^{\frac{1}{N_{\mathrm{prune}}}}$

4 for  _$i\leftarrow 1$ to $i_{\mathrm{train}}$_ do

5 if
_$i\mathbin{\mathchoice{\textstyle{\%}}{{\%}}{{\%}}{{\%}}}i_{\mathrm{prune}}==0$_
then

$s\leftarrow$ criteria($\theta\bigodot m$) // score

$\mathrm{idx}\leftarrow$ argsort($s,p$) // sort

$m[\mathrm{idx}]\leftarrow 0$ // update mask

6

7 else

$\theta\leftarrow$ train($\mathcal{F}(\cdot,\theta\bigodot m)$) // fine-tune

8

9 end if

10

11 end for

12 return _$\theta^{\prime}\leftarrow\theta$_

13 end function

14

Algorithm 1 Iterative Pruning Algorithm.

Tag | Approach | Pruning criteria | Scope
---|---|---|---
L1L | L1 | $|\mathbf{w}_{f}|$ | Local
L1G | L1 | $|\mathbf{w}_{f}|$ | Global
FGL | First-order grad. | $|\frac{\partial\mathcal{L}}{\partial\mathbf{w}}\bigodot\mathbf{w}_{f}|$ | Local
FGG | First-order grad. | $|\frac{\partial\mathcal{L}}{\partial\mathbf{w}}\bigodot\mathbf{w}_{f}|$ | Global
L1SL | L1 w/ same sign | $\frac{\mathbf{w}_{0}}{|\mathbf{w}_{0}|}\cdot\mathbf{w}_{f}$ | Local
L1SG | L1 w/ same sign | $\frac{\mathbf{w}_{0}}{|\mathbf{w}_{0}|}\cdot\mathbf{w}_{f}$ | Global

Table 2: Set of pruning strategies applied to the baseline models. We consider
three different pruning criteria with Local and Global scope, which results in
total six combinations. Here, $\mathcal{L}$ denotes the objective function,
$\mathbf{w}_{0}$ is the initial, pretrained weights, and $\mathbf{w}_{f}$
denotes the weights that can be obtained after fine-tuning.

###  4.2 Network Pruning

Once we finish training a spatially sparse convnet with dense weights, we
prune some of the weights while preserving the test accuracy. In this paper,
we utilize several pruning methods that belong to the family of magnitude-
based pruning [20, 7, 48, 11]. Specifically, we use three pruning methods with
different pruning criteria and apply the criteria on the entire weights
(global) or per layer (local). We put all our pruning criteria on Table 2.
Starting from a pre-trained model, we cycle through a series of the network
pruning and fine-tuning iteratively. At each pruning step, we use the pruning
criteria to sort the weights and remove lower $p$ portion of weights and
repeat this $N_{\mathrm{prune}}$ times so the final target pruning rate is
$p_{\mathrm{target}}=1-(1-p)^{N_{\mathrm{prune}}}$, where $p$ and
$N_{\mathrm{prune}}$ are predefined hyperparameters.

We apply each pruning method to each pre-trained network with the same target
compression rate to compare the effect of different pruning criteria. All
hyperparameters other than pruning strategies, e.g., number of fine-tuning
steps between consecutive pruning steps, learning rates, etc., are set to be
the same for every model and a pruning strategy pair. We put our pruning
procedure on Algorithm 1.

###  4.3 Weight-Sparse Spatially Sparse Convolution

Network compression of spatially sparse convnets results in highly compressed
networks with minimal loss in accuracy. This leads to the reduction in memory
footprint and faster transfer, but, faster inference requires dedicated
software to make use of the compressed convolution kernels and the
computational reduction from pruning. To this end, we extend a pytorch
spatially sparse convolution library [4] with a specialized convolution layer
for weight-sparse spatially sparse convolution (WS3-convolution).

The WS3-convolution consists of several steps:

  * •

Kernel mapping: The spatially sparse convolution computes convolution only on
non-zero elements neighbors. This requires extracting appropriate features
from neighbors and placing the result on the center of convolution. This input
neighbor to output pair is known as kernel mapping and we use it to select
input features.

  * •

CSR compression: We compress the dense weight matrices with pruning masks $m$.
We use the Compressed Sparse Row (CSR) format and convert the masked dense
weight matrices to sparse matrices.

  * •

Transposition: Majority of the off-the-shelf BLAS libraries such as Intel
oneMKL and cuSparse commonly use sparse-dense matrix multiplication routine
with the sparse matrix being the left operand. However, the feature maps are
represented in memory with row-major layout [4, 29]. Thus, we transpose the
feature map for sparse dense multiplication. We pre-transpose weight matrices
before the CSR compression.

  * •

Sparse matrix multiplication: Finally, we use off-the-shelf libraries to
compute temporary features and place the output on the correct position
defined by the kernel mapping.

We put the detailed algorithm on Algorithm 2.

Input : $\mathcal{T}_{\mathrm{in}}=\\{\mathbf{C},\mathbf{F}\\}$: an input
sparse tensor with coordinates in COO layout.

$\mathbf{C}\in\mathbb{Z}^{N\times 3}$ is non-zero discrete coordinates.

$\mathbf{F}\in\mathbb{R}^{N\times N_{\mathrm{in}}}$ is a feature map of
$\mathbf{C}$.

$\mathbf{W}\in\mathbb{R}^{K\times N_{\mathrm{out}}\times N_{\mathrm{in}}}$:
the weight tensor of sparse convolution layer with kernel size $K$, input
feature dimension $N_{\mathrm{in}}$, output dimension $N_{\mathrm{out}}$

Output: $\mathcal{T}_{\mathrm{out}}=\\{\mathbf{C},\mathbf{F}^{\prime}\\}$: the
output sparse tensor.

1 function Weight-Sparse Convolution

// construct kernel in-out map

2 $\mathcal{K}(\cdot)\leftarrow$ constructKernelInOutMap($\mathbf{C},K$)

// compress weights to CSR layout

3
$\mathbf{W}_{\mathrm{sp}}^{T}\leftarrow\text{compress}(\text{transpose}(\mathbf{W}))$

4

// initialize output features map

5 $\mathbf{F}_{\mathrm{out}}\leftarrow$
zeros$(shape=\\{N,N_{\mathrm{out}}\\})$

6 for  _$k\leftarrow 1$ to $K^{3}$_ do

// retrieve kernel in-out map

7 $(\mathcal{K}_{\mathrm{in}},\mathcal{K}_{\mathrm{out}})\leftarrow$
retrieveInOutMap$(\mathbf{C},k)$

// transpose matrices

8
$\mathbf{F}_{\mathrm{in}}^{T}\leftarrow\text{transpose}(\mathbf{F}[\mathcal{K}_{\mathrm{in}},\colon])$

// sparse-dense matmul

9 $\mathbf{F}_{\mathrm{out}}[\mathcal{K}_{\mathrm{out}}]$ +=
transpose$(\mathbf{W}_{\mathrm{sp},k}^{T}@\mathbf{F}_{\mathrm{in}}^{T})$

10

11 end for

12
$\mathcal{T}_{\mathrm{out}}\leftarrow\\{\mathbf{C},\mathbf{F}_{\mathrm{out}}\\}$

13

14 return _$\mathcal{T}_{\mathrm{out}}$_

15 end function

16

Algorithm 2 Weight-Sparse Convolution Algorithm.

##  5 Experiment

|  |  |  ScanNet [6] |  S3DIS [1]
---|---|---|---|---
| # Param.(M) | Rel.(%) | mIoU(%) | Rel.(%) | mAcc(%) | Rel.(%) | mIoU(%) | Rel.(%) | mAcc(%) | Rel.(%)
PointNet [31] | 3.5 |  $\downarrow$ 56.6 | 12.2 |  $\downarrow$ 82.2 | 17.9 |  $\downarrow$ 77.0 | 41.1 |  $\downarrow$ 37.8 | 49.0 |  $\downarrow$ 34.6
PointNet++ [32] | 2.0 |  $\downarrow$ 75.2 | 53.5 |  $\downarrow$ 21.9 | - | - | 57.3 |  $\downarrow$ 13.3 | 63.5 |  $\downarrow$ 15.3
TangentConv [37] | 1.0 |  $\downarrow$ 87.4 | 40.9 |  $\downarrow$ 40.3 | 55.1 |  $\downarrow$ 29.1 | 52.8 |  $\downarrow$ 20.1 | 62.2 |  $\downarrow$ 17.0
PointConv [43] | - | - | 61.0 |  $\downarrow$ 11.0 | - | - | - | - | - | -
PointASNL [45] | - | - | 63.5 |  $\downarrow$ 7.31 | - | - | - | - | - | -
KPConv [38] | 14.9 |  $\uparrow$ 85.9 | 69.2 |  $\uparrow$ 1.01 | - | - | 67.1 |  $\uparrow$ 1.53 | 72.8 |  $\downarrow$ 2.87
PAConv [44] | 11.8 |  $\uparrow$ 47.2 | - | - | - | - | 66.6 |  $\uparrow$ 0.77 | 73.0 |  $\downarrow$ 2.60
PointTransformer [47] | 7.8 |  $\downarrow$ 2.67 | - | - | - | - | 70.4 |  $\uparrow$ 6.52 | 76.5 |  $\uparrow$ 2.07
gray!20 Res16UNet14A | 8.0 |  $\pm$0.00 | 68.5 |  $\pm$0.00 | 77.7 |  $\pm$0.00 | 66.1 |  $\pm$0.00 | 75.0 |  $\pm$0.00
$\vdash$ 90% pruning rate | 0.8 |  $\downarrow$ 89.9 | 67.9 |  $\downarrow$ 0.95 | 77.3 |  $\downarrow$ 0.53 | 65.1 |  $\downarrow$ 1.48 | 72.0 |  $\downarrow$ 3.90
$\vdash$ 95% pruning rate | 0.4 |  $\downarrow$ 94.9 | 67.2 |  $\downarrow$ 1.96 | 76.9 |  $\downarrow$ 1.08 | 64.8 |  $\downarrow$ 1.91 | 72.2 |  $\downarrow$ 3.64
$\vdash$ 99% pruning rate | 0.1 |  $\downarrow$ 98.9 | 59.3 |  $\downarrow$ 13.4 | 70.3 |  $\downarrow$ 9.46 | 60.9 |  $\downarrow$ 7.87 | 68.7 |  $\downarrow$ 8.35
Res16UNet18A | 15.5 |  $\uparrow$ 93.4 | 70.3 |  $\uparrow$ 2.58 | 79.7 |  $\uparrow$ 2.54 | 67.6 |  $\uparrow$ 2.30 | 74.6 |  $\downarrow$ 0.51
$\vdash$ 90% pruning rate | 1.6 |  $\downarrow$ 80.6 | 69.9 |  $\uparrow$ 2.03 | 79.0 |  $\uparrow$ 1.72 | 66.5 |  $\uparrow$ 0.57 | 73.0 |  $\downarrow$ 2.60
$\vdash$ 95% pruning rate | 0.8 |  $\downarrow$ 90.2 | 69.3 |  $\uparrow$ 1.09 | 78.6 |  $\uparrow$ 1.20 | 66.3 |  $\uparrow$ 0.30 | 73.2 |  $\downarrow$ 2.32
$\vdash$ 99% pruning rate | 0.2 |  $\downarrow$ 97.9 | 65.7 |  $\downarrow$ 4.10 | 76.1 |  $\downarrow$ 2.03 | 64.2 |  $\downarrow$ 2.84 | 71.4 |  $\downarrow$ 4.68
Res16UNet34C | 37.9 |  $\uparrow$ 372.9 | 71.6 |  $\uparrow$ 4.47 | 80.4 |  $\uparrow$ 3.53 | 68.8 |  $\uparrow$ 4.04 | 75.8 |  $\uparrow$ 1.17
$\vdash$ 90% pruning rate | 3.8 |  $\downarrow$ 52.6 | 71.0 |  $\uparrow$ 3.68 | 79.8 |  $\uparrow$ 2.70 | 66.5 |  $\uparrow$ 0.59 | 73.1 |  $\downarrow$ 2.44
$\vdash$ 95% pruning rate | 1.9 |  $\downarrow$ 76.2 | 71.0 |  $\uparrow$ 3.56 | 79.7 |  $\uparrow$ 2.63 | 66.4 |  $\uparrow$ 0.48 | 73.2 |  $\downarrow$ 2.31
$\vdash$ 99% pruning rate | 0.4 |  $\downarrow$ 95.1 | 69.5 |  $\uparrow$ 1.40 | 78.9 |  $\uparrow$ 1.60 | 65.7 |  $\downarrow$ 0.67 | 72.8 |  $\downarrow$ 2.92

Table 3: Quantitative results of _semantic segmentation_ models on ScanNet [6]
and S3DIS dataset [1]. Gray row indicates the reference network, and we color
green for the positive changes and red for the negative changes w.r.t. the
reference network.

|  |  |  ScanNet [6] |  S3DIS [1]
---|---|---|---|---
| # Param.(M) | Rel.(%) | mAP50(%) | Rel.(%) | mAP25(%) | Rel.(%) | mAP50(%) | Rel.(%) | mAP25(%) | Rel.(%)
gray!20 Res16UNet14A | 8.0 |  $\pm$0.00 | 52.8 |  $\pm$0.00 | 68.1 |  $\pm$0.00 | 54.5 |  $\pm$0.00 | 61.7 |  $\pm$0.00
$\vdash$ 90% pruning rate | 0.8 |  $\downarrow$ 89.9 | 51.2 |  $\downarrow$ 3.03 | 66.4 |  $\downarrow$ 2.50 | 48.7 |  $\downarrow$ 10.6 | 60.7 |  $\downarrow$ 1.62
$\vdash$ 95% pruning rate | 0.4 |  $\downarrow$ 94.9 | 48.9 |  $\downarrow$ 7.39 | 65.4 |  $\downarrow$ 3.96 | 48.7 |  $\downarrow$ 10.6 | 59.4 |  $\downarrow$ 3.73
$\vdash$ 99% pruning rate | 0.1 |  $\downarrow$ 98.9 | 34.6 |  $\downarrow$ 34.5 | 55.9 |  $\downarrow$ 17.9 | 37.7 |  $\downarrow$ 30.8 | 52.2 |  $\downarrow$ 15.4
Res16UNet18A | 15.5 |  $\uparrow$ 93.3 | 55.8 |  $\uparrow$ 5.68 | 71.3 |  $\uparrow$ 4.70 | 55.6 |  $\uparrow$ 2.02 | 62.7 |  $\uparrow$ 1.62
$\vdash$ 90% pruning rate | 1.6 |  $\downarrow$ 80.5 | 55.7 |  $\uparrow$ 5.49 | 71.8 |  $\uparrow$ 5.43 | 48.5 |  $\downarrow$ 11.0 | 59.6 |  $\downarrow$ 3.40
$\vdash$ 95% pruning rate | 0.8 |  $\downarrow$ 90.2 | 55.5 |  $\uparrow$ 5.11 | 71.2 |  $\uparrow$ 4.55 | 49.8 |  $\downarrow$ 8.62 | 58.7 |  $\downarrow$ 4.86
$\vdash$ 99% pruning rate | 0.2 |  $\downarrow$ 97.9 | 47.3 |  $\downarrow$ 10.4 | 65.3 |  $\downarrow$ 4.11 | 45.1 |  $\downarrow$ 17.2 | 56.6 |  $\downarrow$ 8.27
Res16UNet34C | 37.9 |  $\uparrow$ 372.6 | 57.0 |  $\uparrow$ 7.95 | 73.1 |  $\uparrow$ 7.34 | 54.7 |  $\uparrow$ 0.37 | 64.1 |  $\uparrow$ 3.89
$\vdash$ 90% pruning rate | 3.8 |  $\downarrow$ 52.6 | 56.9 |  $\uparrow$ 7.77 | 71.7 |  $\uparrow$ 5.29 | 48.3 |  $\downarrow$ 11.4 | 59.4 |  $\downarrow$ 3.73
$\vdash$ 95% pruning rate | 1.9 |  $\downarrow$ 76.2 | 56.3 |  $\uparrow$ 6.63 | 71.5 |  $\uparrow$ 4.99 | 51.0 |  $\downarrow$ 6.42 | 60.0 |  $\downarrow$ 2.76
$\vdash$ 99% pruning rate | 0.4 |  $\downarrow$ 95.0 | 54.2 |  $\uparrow$ 2.65 | 71.0 |  $\uparrow$ 4.26 | 48.9 |  $\downarrow$ 10.3 | 58.4 |  $\downarrow$ 5.35

Table 4: Quantitative results of _instance segmentation_ models on ScanNet [6]
and S3DIS dataset [1]. Gray row indicates the reference network, and we color
green for the positive changes and red for the negative changes w.r.t. the
reference network.

In this section, we present the results and analysis of our network pruning
experiments that require in total 5800 hours of GPU computation on 30 GPUs.

###  5.1 Experiment setup and metrics

We use semantic segmentation and instance segmentation on two indoor datasets
for the compression of spatially sparse convnets. ScanNet [6] consists of
$1.5$k scenes captured with commercial RGB-D cameras. It provides diverse
indoor scenes but exhibits incomplete reconstruction and occlusion. On the
other hand, the S3DIS [1] dataset provides 271 indoor office rooms captured
with LiDAR 3D scanning devices. The reconstruction is complete and well
aligned, but the scenes are monotonic compared with ScanNet.

We measure the intersection of union (IoU) and per-point classification
accuracy for each semantic class, and report the mean IoUs (mIoU), and mean
accuracy (mAcc) for semantic segmentation. For instance segmentation, we
report average precision per class and mean average precision with varying
accuracies 50% and 25% which we denote as mAP50 and mAP25 respectively.

We use three different network architectures with varying depth and width:
Res16UNet14A, Res16UNet18A, Res16UNet34C [4]. We apply unstructured pruning
[11, 20, 7, 48] and compare the results with representative point-based
approaches [31, 32, 37, 43, 45, 38, 44, 47].

###  5.2 Implementation detail

We train all our networks for $60K$ iterations with an initial learning rate
$0.1$, $2cm$ voxelization, and batch size 8 on both datasets. After the
baseline networks are converged, we apply our iterative pruning with
$i_{\mathrm{train}}=60K$ total training iteration, and $N_{\mathrm{prune}}=10$
pruning and fine-tuning sequences with $i_{\mathrm{prune}}=2K$ fine-tuning
iterations between the consecutive pruning sequences for all experiments. All
the experiments are evaluated on an NVIDIA Titan RTX and Intel Xeon Gold
5220R.

Figure 2: (Top) Res16UNet network architecture. It is a U-shape network with
an encoder (yellow) and decoder (blue), and has skip connections. (Bottom)
Visualization of compressed convolution kernels. Each grid represents the
average density of the first $3\times 3\times 3$ convolution kernel in each
residual block from encoder to decoder from left to right. Note that kernels
in decoder are pruned more aggressively in general, but weights along z-axis
tend to survive.  Figure 3: Ablation study on different pruning criteria on
ScanNet semantic segmentation. For all architectures, L1G (Global L1
magnitude-based pruning) shows the least drop in performance. We use 34C, 18A,
and 14A to denote Res16UNet variants. For pruning method tags, please refer to
Table 2.

###  5.3 Experiment Results and Analysis

We compress various backbone networks with multiple pruning rates and present
the semantic segmentation results on Table 3 and instance segmentation results
on Table 4. We present our analysis of these results in this section.

Conventional 3D ConvNets are bulky. On Table 3, we can see that KPConv [38],
PAConv [44], and PointTransformer [47] achieve high performance on S3DIS while
having as many parameters as Res16UNet14A. However, even when we prune 90% of
parameters from Res16UNet14A, the relative performance drop is only
0.5$\sim$3.9%. On the other hand, PointNet [31], PointNet++ [32], and
TangentConv [37] have almost the same number of parameters as 90% pruned
Res16UNet14A. However, the 90% pruned ResNet outperforms these networks by a
decent margin. Similarly, PAConv [44] having 11.8M parameters shows comparable
results to 90% pruned Res16UNet18A with 1.6M parameters.

Prune overparametrized networks. Another interesting observation is that we
see a small performance drop when we compress heavier networks. For example,
in Table 3, Res16UNet14A with a 95% compression (Network A) and Res16UNet34C
with a 99% compression (Network B) both have a similar number of parameters.
Yet, Network B consistently shows higher accuracy on ScanNet and S3DIS. We
speculate that heavier and deeper networks still retain the same depths
regardless of compression and overparametrization leads to a higher chance of
maintaining more informative connections during pruning.

Global pruning is better than local pruning. We visualize the performance vs.
network size graph on Figure 3. Global pruning that prunes the entire network
weights perform better than local pruning that prunes the same percentage of
weights from each layer. For instance, ‘34C+L1G’ shows smaller drop in
performance than ‘34C+L1L’. This is common in all backbone networks and we
will explain this in the next paragraph.

Decoder of U-Net is overparametrized. We visualize the pruning patterns of
convolution kernels on Figure 2 (bottom). Note that the encoder retains more
connections while the decoder is pruned aggressively. This indicates that the
decoder has more redundant connections and shows that pruning globally allows
the pruning algorithm to choose the importance freely among all layers leading
to higher performance than local pruning.

Kernels along gravity-axis are more important. Pruning patterns on Figure 2
(bottom) show that kernels along the gravity-axis have higher density than any
other directions. We speculate that the vertical geometric patterns are
especially important in indoor scenes due to structural similarity and height
from the ground is critical in indoor perception.

Figure 4: Visualization of semantic segmentation prediction on ScanNet [17]
using Res16UNet34C with varying sparsity (Dense, 90%, 95%, and 99% pruned
models from left to right).

###  5.4 WS3-ConvNet

Based on the experiments and observations, we pick Res16UNet34C network with
99% pruning rate as our representative result, and name it WS3-ConvNet. As
shown in Table 3, Table 4, Figure 1, and Figure 4, WS3-ConvNet shows very high
accuracy (69.4% mIoU) with only 0.4M parameters, which is even smaller than
PointNet [31] and PointNet++ [32].

In addition, we test our weight-sparse convolution algorithm in Algorithm 2
and measure the speed up on Table 5. Compared with the dense inference time,
the sparse time, which measures the inference time for 99% pruned networks,
shows up to 34% speedup and the speedup is consistent across tasks and
networks.

###  5.5 Structural Pruning of Spatial Kernels

Lastly, we utilize the fact that, on the decoder side, the convolution kernels
along the gravity-axis are more important, and further compress the networks
using structural pruning. We _switch-off_ the kernels that are y off the
z-axis and present the performance drop and speedup in Table 6. Note that
switching off kernels in the last few layers doesn’t harm the accuracy much,
while it enables faster inference.

|  ScanNet [6] |  S3DIS [1]
---|---|---
Arch. | Dense time | Sparse time | Rel.(%) | Dense time | Sparse time | Rel.(%)
| Semantic Segmentation
Res16UNet14A | 2.63 | 1.95 |  $\downarrow$ 25.8 | 1.77 | 1.23 |  $\downarrow$ 30.3
Res16UNet18A | 4.55 | 3.06 |  $\downarrow$ 32.9 | 2.83 | 1.96 |  $\downarrow$ 30.9
Res16UNet34C | 4.57 | 3.25 |  $\downarrow$ 28.8 | 2.84 | 2.14 |  $\downarrow$ 24.7
| Instance Segmentation
Res16UNet14A | 2.83 | 1.95 |  $\downarrow$ 31.1 | 1.78 | 1.33 |  $\downarrow$ 25.1
Res16UNet18A | 4.47 | 3.20 |  $\downarrow$ 28.4 | 3.09 | 2.04 |  $\downarrow$ 33.9
Res16UNet34C | 4.70 | 3.28 |  $\downarrow$ 30.2 | 3.13 | 2.13 |  $\downarrow$ 31.8

Table 5: Comparison on the inference speed of spatially-sparse convolutional
networks and WS3-ConvNet. Times are the average latency per scene and measure
on CPU. Our implementation of WS3-ConvNet exhibits up to 34% reduction of
inference time.

B1 | B2 | B3 | B4 | B5 | B6 | B7 | B8 | mIoU(%) | mAcc(%) | Time(s)
---|---|---|---|---|---|---|---|---|---|---
gray!20 |  |  |  |  |  |  |  |  71.6(+3.02) |  80.4(+1.90) |  4.57(+38.5)
|  |  |  |  |  |  |  |  69.5($\pm$0.00) |  78.9($\pm$0.00) |  3.30($\pm$0.00)
|  |  |  |  |  |  | ✓ |  69.2(-0.42) |  78.8(-0.18) |  2.59(-21.6)
|  |  |  |  |  | ✓ | ✓ |  64.3(-7.39) |  73.5(-6.90) |  2.26(-31.5)
|  |  |  |  | ✓ | ✓ | ✓ |  40.6(-41.6) |  46.0(-41.7) |  2.12(-35.7)
|  |  |  | ✓ | ✓ | ✓ | ✓ |  11.9(-82.9) |  15.7(-80.2) |  2.25(-31.9)
|  |  | ✓ | ✓ | ✓ | ✓ | ✓ |  3.65(-94.8) |  7.10(-91.0) |  2.04(-38.2)
|  | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ |  3.50(-95.0) |  6.92(-91.2) |  1.92(-42.0)
| ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ |  2.61(-96.2) |  6.03(-92.4) |  1.85(-44.1)
✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ |  5.63(-91.9) |  9.73(-87.7) |  1.72(-47.8)

Table 6: Ablation study on the effect of z-axis aligned kernels. Columns B1 to
B8 indicate residual blocks in our model (Res16UNet34C, 99% pruned),
$\checkmark$ denotes that the only weights along the z-axis of $3^{3}$ kernel
are used to compute output features for all conv layers in the block. All
inference time measured on CPU. Gray row indicates the dense (0% pruned)
model.

##  6 Discussion

In this work, we apply various network compression methods on 3D spatially
sparse convolutional networks. We found that 3D convnets can be compressed
effectively to retain only 1% of its original weights while maintaining very
high accuracy. Also, to speed up the inference, we propose weight-sparse
convolution algorithm and show up to 34% speedup. Analysis shows interesting
properties of 3D convnets and we further utilize the findings to apply
additional structural pruning.

We only pruned spatially sparse networks, but would like to incorporate other
networks for compression. In addition, we only used unstructured pruning but
will study the effect of structured pruning and network architecture search
(NAS) to find the optimal compressed network.

## References

  * [1] Iro Armeni, Ozan Sener, Amir R. Zamir, Helen Jiang, Ioannis Brilakis, Martin Fischer, and Silvio Savarese.  3d semantic parsing of large-scale indoor spaces.  In Proceedings of the IEEE International Conference on Computer Vision and Pattern Recognition, 2016.
  * [2] J. Behley, M. Garbade, A. Milioto, J. Quenzel, S. Behnke, J. Gall, and C. Stachniss.  Towards 3D LiDAR-based semantic scene understanding of 3D point cloud sequences: The SemanticKITTI Dataset.  The International Journal on Robotics Research, 40(8-9):959–967, 2021.
  * [3] Christopher Choy, Wei Dong, and Vladlen Koltun.  Deep global registration.  In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, pages 2514–2523, 2020.
  * [4] Christopher Choy, JunYoung Gwak, and Silvio Savarese.  4d spatio-temporal convnets: Minkowski convolutional neural networks.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 3075–3084, 2019.
  * [5] Christopher Choy, Jaesik Park, and Vladlen Koltun.  Fully convolutional geometric features.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 8958–8966, 2019.
  * [6] Angela Dai, Angel X Chang, Manolis Savva, Maciej Halber, Thomas Funkhouser, and Matthias Nießner.  Scannet: Richly-annotated 3d reconstructions of indoor scenes.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 5828–5839, 2017.
  * [7] Jonathan Frankle and Michael Carbin.  The lottery ticket hypothesis: Finding sparse, trainable neural networks.  In International Conference on Learning Representations, 2019.
  * [8] Jonathan Frankle, Gintare Karolina Dziugaite, Daniel M Roy, and Michael Carbin.  Stabilizing the lottery ticket hypothesis.  arXiv preprint arXiv:1903.01611, 2019.
  * [9] Benjamin Graham.  Spatially-sparse convolutional neural networks.  arXiv preprint arXiv:1409.6070, 2014.
  * [10] JunYoung Gwak, Christopher Choy, and Silvio Savarese.  Generative sparse detection networks for 3d single-shot object detection.  In Computer Vision–ECCV 2020: 16th European Conference, Glasgow, UK, August 23–28, 2020, Proceedings, Part IV 16, pages 297–313. Springer, 2020.
  * [11] Song Han, Jeff Pool, John Tran, and William J Dally.  Learning both weights and connections for efficient neural network.  In NIPS, 2015.
  * [12] Babak Hassibi, David G Stork, and Gregory J Wolff.  Optimal brain surgeon and general network pruning.  In IEEE international conference on neural networks, pages 293–299. IEEE, 1993.
  * [13] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 770–778, 2016.
  * [14] Yihui He, Xiangyu Zhang, and Jian Sun.  Channel pruning for accelerating very deep neural networks.  In Proceedings of the IEEE international conference on computer vision, pages 1389–1397, 2017.
  * [15] Andrew G Howard, Menglong Zhu, Bo Chen, Dmitry Kalenichenko, Weijun Wang, Tobias Weyand, Marco Andreetto, and Hartwig Adam.  Mobilenets: Efficient convolutional neural networks for mobile vision applications.  arXiv preprint arXiv:1704.04861, 2017.
  * [16] Zehao Huang and Naiyan Wang.  Data-driven sparse structure selection for deep neural networks.  In Proceedings of the European conference on computer vision (ECCV), pages 304–320, 2018.
  * [17] Forrest N Iandola, Song Han, Matthew W Moskewicz, Khalid Ashraf, William J Dally, and Kurt Keutzer.  Squeezenet: Alexnet-level accuracy with 50x fewer parameters and¡ 0.5 mb model size.  arXiv preprint arXiv:1602.07360, 2016.
  * [18] Sergey Ioffe and Christian Szegedy.  Batch normalization: Accelerating deep network training by reducing internal covariate shift.  In International conference on machine learning, pages 448–456. PMLR, 2015.
  * [19] Yann LeCun, John S Denker, and Sara A Solla.  Optimal brain damage.  In NIPS, pages 598–605, 1990.
  * [20] Namhoon Lee, Thalaiyasingam Ajanthan, and Philip Torr.  Snip: Single-shot network pruning based on connection sensitivity.  In International Conference on Learning Representations, 2018.
  * [21] Hao Li, Asim Kadav, Igor Durdanovic, Hanan Samet, and Hans Peter Graf.  Pruning filters for efficient convnets.  International Conference on Learning Representations, 2017.
  * [22] Zhuang Liu, Jianguo Li, Zhiqiang Shen, Gao Huang, Shoumeng Yan, and Changshui Zhang.  Learning efficient convolutional networks through network slimming.  In Proceedings of the IEEE international conference on computer vision, pages 2736–2744, 2017.
  * [23] Jian-Hao Luo, Jianxin Wu, and Weiyao Lin.  Thinet: A filter level pruning method for deep neural network compression.  In Proceedings of the IEEE international conference on computer vision, pages 5058–5066, 2017.
  * [24] Jiageng Mao, Xiaogang Wang, and Hongsheng Li.  Interpolated convolutional networks for 3d point cloud understanding.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 1578–1587, 2019.
  * [25] Lars Mescheder, Michael Oechsle, Michael Niemeyer, Sebastian Nowozin, and Andreas Geiger.  Occupancy networks: Learning 3d reconstruction in function space.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 4460–4470, 2019.
  * [26] Pavlo Molchanov, Stephen Tyree, Tero Karras, Timo Aila, and Jan Kautz.  Pruning convolutional neural networks for resource efficient inference.  arXiv preprint arXiv:1611.06440, 2016.
  * [27] Michael C Mozer and Paul Smolensky.  Skeletonization: A technique for trimming the fat from a network via relevance assessment.  In Advances in neural information processing systems, pages 107–115, 1989.
  * [28] Alexey Nekrasov, Jonas Schult, Or Litany, Bastian Leibe, and Francis Engelmann.  Mix3d: Out-of-context data augmentation for 3d scenes.  arXiv preprint arXiv:2110.02210, 2021.
  * [29] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Kopf, Edward Yang, Zachary DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala.  Pytorch: An imperative style, high-performance deep learning library.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems 32, pages 8024–8035. Curran Associates, Inc., 2019.
  * [30] Charles R Qi, Or Litany, Kaiming He, and Leonidas J Guibas.  Deep hough voting for 3d object detection in point clouds.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 9277–9286, 2019.
  * [31] Charles R Qi, Hao Su, Kaichun Mo, and Leonidas J Guibas.  Pointnet: Deep learning on point sets for 3d classification and segmentation.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 652–660, 2017.
  * [32] Charles R Qi, Li Yi, Hao Su, and Leonidas J Guibas.  Pointnet++: Deep hierarchical feature learning on point sets in a metric space.  arXiv preprint arXiv:1706.02413, 2017.
  * [33] Olaf Ronneberger, Philipp Fischer, and Thomas Brox.  U-net: Convolutional networks for biomedical image segmentation.  In International Conference on Medical image computing and computer-assisted intervention, pages 234–241. Springer, 2015.
  * [34] Pei Sun, Henrik Kretzschmar, Xerxes Dotiwalla, Aurelien Chouard, Vijaysai Patnaik, Paul Tsui, James Guo, Yin Zhou, Yuning Chai, Benjamin Caine, et al.  Scalability in perception for autonomous driving: Waymo open dataset.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 2446–2454, 2020.
  * [35] Mingxing Tan and Quoc Le.  Efficientnet: Rethinking model scaling for convolutional neural networks.  In International Conference on Machine Learning, pages 6105–6114. PMLR, 2019.
  * [36] Hidenori Tanaka, Daniel Kunin, Daniel LK Yamins, and Surya Ganguli.  Pruning neural networks without any data by iteratively conserving synaptic flow.  arXiv preprint arXiv:2006.05467, 2020.
  * [37] Maxim Tatarchenko, Jaesik Park, Vladlen Koltun, and Qian-Yi Zhou.  Tangent convolutions for dense prediction in 3d.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 3887–3896, 2018.
  * [38] Hugues Thomas, Charles R Qi, Jean-Emmanuel Deschaud, Beatriz Marcotegui, François Goulette, and Leonidas J Guibas.  Kpconv: Flexible and deformable convolution for point clouds.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 6411–6420, 2019.
  * [39] Chaoqi Wang, Guodong Zhang, and Roger Grosse.  Picking winning tickets before training by preserving gradient flow.  arXiv preprint arXiv:2002.07376, 2020.
  * [40] Yue Wang, Yongbin Sun, Ziwei Liu, Sanjay E Sarma, Michael M Bronstein, and Justin M Solomon.  Dynamic graph cnn for learning on point clouds.  Acm Transactions On Graphics (tog), 38(5):1–12, 2019.
  * [41] Zhenzhen Wang, Weixiang Hong, Yap-Peng Tan, and Junsong Yuan.  Pruning 3d filters for accelerating 3d convnets.  IEEE Transactions on Multimedia, 22(8):2126–2137, 2020.
  * [42] Zhenzhen Wang, Weixiang Hong, Yap-Peng Tan, and Junsong Yuan.  Pruning 3d filters for accelerating 3d convnets.  IEEE Transactions on Multimedia, 22(8):2126–2137, 2020.
  * [43] Wenxuan Wu, Zhongang Qi, and Li Fuxin.  Pointconv: Deep convolutional networks on 3d point clouds.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 9621–9630, 2019.
  * [44] Mutian Xu, Runyu Ding, Hengshuang Zhao, and Xiaojuan Qi.  Paconv: Position adaptive convolution with dynamic kernel assembling on point clouds.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 3173–3182, 2021.
  * [45] Xu Yan, Chaoda Zheng, Zhen Li, Sheng Wang, and Shuguang Cui.  Pointasnl: Robust point clouds processing using nonlocal neural networks with adaptive sampling.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 5589–5598, 2020.
  * [46] Xiangyu Zhang, Xinyu Zhou, Mengxiao Lin, and Jian Sun.  Shufflenet: An extremely efficient convolutional neural network for mobile devices.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 6848–6856, 2018.
  * [47] Hengshuang Zhao, Li Jiang, Jiaya Jia, Philip HS Torr, and Vladlen Koltun.  Point transformer.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 16259–16268, 2021.
  * [48] Hattie Zhou, Janice Lan, Rosanne Liu, and Jason Yosinski.  Deconstructing lottery tickets: Zeros, signs, and the supermask.  arXiv preprint arXiv:1905.01067, 2019.

##  A Appendix

In this appendix, we provide additional details and results of our method.

###  A.1 Architecture details

For all experiments, we used three baseline networks for network compression
and analysis. We provide the architectural details of the these networks. As
illustrated in Table a.1, the networks share the same backbone but have
varying depths and widths. The backbone has four steps of downsampling and
upsampling using convolution and transposed convolution layers with stride 2.
Hence, at the very bottom of the bottleneck, the resolution of an input voxel
reduces 16 times than its initial value. In total, there are 8 residual blocks
between upsampling and downsampling. Between the same upsampled and
downsampled features with the same resolution, it has skip connections. For
semantic segmentation, the features are projected to logits using $1\times
1\times 1$ convolution layer with $N_{C}$ output dimension, where $N_{C}=20$
for ScanNet [6] and 13 for S3DIS [1]. Note that instance segmentation models
have an additional head for predicting the offset coordinates for instance
centroids.

###  A.2 Additional analysis of pruning criteria.

In Figure a.1, we visualize the performance vs. model size graphs on S3DIS
semantic segmentation (left), ScanNet instance segmentation (middle), and
S3DIS instance segmentation (right). For semantic segmentation model, we use
mIoU for our main performance metric; for instance segmentation, we use mAP50
for our evaluation metric.

We found that the global pruning strategies show better pruning efficiency
than the local pruning methods across all datasets and tasks. Note that in the
S3DIS instance segmentation experiment, some pruning strategies show the best
performance on a 95% pruned model, not on 90% pruned models. We speculate that
the 95% pruned model generalizes better than the larger 90% pruned models
because the instance segmentation task on S3DIS is monotonic and lack
diversity compared with the ScanNet dataset.

Figure a.1: Ablation study on different pruning criteria on S3DIS [1] semantic
segmentation (Left), ScanNet [6] instance segmentation (Middle), and S3DIS
instance segmentation (Right).

layer name | Res16UNet14A | Res16UNet18A | Res16UNet34C
---|---|---|---
conv0 |  $3\times 3\times 3$, 32
conv1 |  $2\times 2\times 2$, 32, stride 2
block1 | $\begin{bmatrix}3\times 3\times 3,\text{32}\\\ 3\times 3\times 3,32\end{bmatrix}\times 1$ | $\begin{bmatrix}3\times 3\times 3,\text{32}\\\ 3\times 3\times 3,32\end{bmatrix}\times 2$ | $\begin{bmatrix}3\times 3\times 3,\text{32}\\\ 3\times 3\times 3,32\end{bmatrix}\times 2$
conv2 |  $2\times 2\times 2$, 32, stride 2
block2 | $\begin{bmatrix}3\times 3\times 3,\text{64}\\\ 3\times 3\times 3,\text{64}\end{bmatrix}\times 1$ | $\begin{bmatrix}3\times 3\times 3,\text{64}\\\ 3\times 3\times 3,\text{64}\end{bmatrix}\times 2$ | $\begin{bmatrix}3\times 3\times 3,\text{64}\\\ 3\times 3\times 3,\text{64}\end{bmatrix}\times 3$
conv3 |  $2\times 2\times 2\times$, 64, stride 2
block3 | $\begin{bmatrix}3\times 3\times 3,\text{128}\\\ 3\times 3\times 3,\text{128}\end{bmatrix}\times 1$ | $\begin{bmatrix}3\times 3\times 3,\text{128}\\\ 3\times 3\times 3,\text{128}\end{bmatrix}\times 2$ | $\begin{bmatrix}3\times 3\times 3,\text{128}\\\ 3\times 3\times 3,\text{128}\end{bmatrix}\times 4$
conv4 |  $2\times 2\times 2$, 128, stride 2
block4 | $\begin{bmatrix}3\times 3\times 3,\text{256}\\\ 3\times 3\times 3,\text{256}\end{bmatrix}\times 1$ | $\begin{bmatrix}3\times 3\times 3,\text{256}\\\ 3\times 3\times 3,\text{256}\end{bmatrix}\times 2$ | $\begin{bmatrix}3\times 3\times 3,\text{256}\\\ 3\times 3\times 3,\text{256}\end{bmatrix}\times 6$
conv4_tr |  $2\times 2\times 2$, 128, stride 2 |  $2\times 2\times 2$, 256, stride 2
block5 | $\begin{bmatrix}3\times 3\times 3,\text{128}\\\ 3\times 3\times 3,\text{128}\end{bmatrix}\times 1$ | $\begin{bmatrix}3\times 3\times 3,\text{128}\\\ 3\times 3\times 3,\text{128}\end{bmatrix}\times 2$ | $\begin{bmatrix}3\times 3\times 3,\text{256}\\\ 3\times 3\times 3,\text{256}\end{bmatrix}\times 2$
conv5_tr |  $2\times 2\times 2$, 128, stride 2
block6 | $\begin{bmatrix}3\times 3\times 3,\text{128}\\\ 3\times 3\times 3,\text{128}\end{bmatrix}\times 1$ | $\begin{bmatrix}3\times 3\times 3,\text{128}\\\ 3\times 3\times 3,\text{128}\end{bmatrix}\times 2$ | $\begin{bmatrix}3\times 3\times 3,\text{128}\\\ 3\times 3\times 3,\text{128}\end{bmatrix}\times 2$
conv6_tr |  $2\times 2\times 2$, 96, stride 2
block7 | $\begin{bmatrix}3\times 3\times 3,\text{96}\\\ 3\times 3\times 3,\text{96}\end{bmatrix}\times 1$ | $\begin{bmatrix}3\times 3\times 3,\text{96}\\\ 3\times 3\times 3,\text{96}\end{bmatrix}\times 2$ | $\begin{bmatrix}3\times 3\times 3,\text{96}\\\ 3\times 3\times 3,\text{96}\end{bmatrix}\times 2$
conv7_tr |  $2\times 2\times 2$, 96, stride 2
block8 | $\begin{bmatrix}3\times 3\times 3,\text{96}\\\ 3\times 3\times 3,\text{96}\end{bmatrix}\times 1$ | $\begin{bmatrix}3\times 3\times 3,\text{96}\\\ 3\times 3\times 3,\text{96}\end{bmatrix}\times 2$ | $\begin{bmatrix}3\times 3\times 3,\text{96}\\\ 3\times 3\times 3,\text{96}\end{bmatrix}\times 2$
final | $1\times 1\times 1,N_{C}$
offset∗ | $\begin{bmatrix}1\times 1\times 1,96\\\ 1\times 1\times 1,3\end{bmatrix}$
Number of Params | $8.02\times 10^{6}$ | $15.5\times 10^{6}$ | $37.9\times 10^{6}$

Table a.1: Architectures of Res16UNet variants. We denote a convolution layer
with its kernel size, output channel size, and convolution stride size. If a
stride side is not specified, it is assumed that there is no stride. All
convolution layers except for the last layer have a Batch Normalization and a
ReLU layer after them. The layers with the tag ”conv_tr” indicates the
transposed convolution layers. We use a square bracket to denote a residual
block, with the number of blocks stacked. For instance segmentation, the
additional head (offset∗) for offset estimation exists, as depicted in Figure
2 (Top) in our main paper.

###  A.3 Additional kernel visualization.

In Figure a.2, we visualize the compressed convolution kernels with different
pruning strategies, or different compression rates. All models are
Res16UNet34C trained and pruned on ScanNet semantic segmentation.

As shown in the figure, all pruning strategies tend to prune the weights of
decoder more aggressively. While the models pruned with global pruning methods
have relatively dense encoder, the local pruning methods make the encoder
sparser. Furthermore, we identify that the pattern of convolution kernels,
alignment along the gravity axis, reveals as the model is compressed more
significantly.

Criteria |  | Prune
---
rate
Kernel Visualization |  | mIoU
---
(%)
L1L | 99% |  |
---
54.7
FGL | 99% |  |
---
34.4
FGG | 99% |  |
---
68.5
L1SL | 99% |  |
---
58.1
L1SG | 99% |  |
---
68.9
L1G | 99% |  |
---
69.5
L1G | 95% |  |
---
71.0
L1G | 90% |  |
---
71.1

Figure a.2: Visualizations of compressed convolution kernels. We use
Res16UNet34C for all kernel visualizations. We present convolution kernels
compressed with different pruning criteria on the first 6 rows, and kernels
compressed with different compression rate on the last three rows.

###  A.4 Additional qualitative results.

We provide additional visualizations of semantic segmentation results on the
S3DIS dataset (Figure a.3), instance segmentation on the ScanNet dataset
(Figure a.4), and instance segmentation on the S3DIS dataset (Figure a.5). We
use compressed Res16UNet34C models for all experiments.

|  |  |
---|---|---|---
|  |  |
|  |  |
| (a) Reference Network
---
(Res16UNet34C)
(b) 90% pruned | (c) 95% pruned | (d) 99% pruned

Figure a.3: Semantic segmentation qualitative results on S3DIS [1].

|  |  |
---|---|---|---
|  |  |
|  |  |
| (a) Reference Network
---
(Res16UNet34C)
(b) 90% pruned | (c) 95% pruned | (d) 99% pruned

Figure a.4: Instance segmentation qualitative results on ScanNet [6].

|  |  |
---|---|---|---
|  |  |
|  |  |
|  |  |
| (a) Reference Network
---
(Res16UNet34C)
(b) 90% pruned | (c) 95% pruned | (d) 99% pruned

Figure a.5: Instance segmentation qualitative results on S3DIS [1].


# Multiclass Diffuse Interface Models for
Semi-Supervised Learning on Graphs

Cristina Garcia-Cardona1, Arjuna Flenner2 and Allon G. Percus1
1Institute of Mathematical Sciences, Claremont Graduate University, Claremont,
CA 91711, USA
2Naval Air Warfare Center, Physics and Computational Sciences, China Lake, CA
93555, USA
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>

###### Abstract

We present a graph-based variational algorithm for multiclass classification
of high-dimensional data, motivated by total variation techniques. The energy
functional is based on a diffuse interface model with a periodic potential. We
augment the model by introducing an alternative measure of smoothness that
preserves symmetry among the class labels. Through this modification of the
standard Laplacian, we construct an efficient multiclass method that allows
for sharp transitions between classes. The experimental results demonstrate
that our approach is competitive with the state of the art among other graph-
based algorithms.

##  1 INTRODUCTION

Many tasks in pattern recognition and machine learning rely on the ability to
quantify local similarities in data, and to infer meaningful global structure
from such local characteristics [Coifman et al., 2005]. In the classification
framework, the desired global structure is a descriptive partition of the data
into categories or classes. Many studies have been devoted to the binary
classification problems. The multiple-class case, where the data is
partitioned into more than two clusters, is more challenging. One approach is
to treat the problem as a series of binary classification problems [Allwein et
al., 2000]. In this paper, we develop an alternative method, involving a
multiple-class extension of the diffuse interface model introduced in
[Bertozzi and Flenner, 2012].

The diffuse interface model by Bertozzi and Flenner combines methods for
diffusion on graphs with efficient partial differential equation techniques to
solve binary segmentation problems. As with other methods inspired by physical
phenomena [Bertozzi et al., 2007, Jung et al., 2007, Li and Kim, 2011], it
requires the minimization of an energy expression, specifically the Ginzburg-
Landau (GL) energy functional. The formulation generalizes the GL functional
to the case of functions defined on graphs, and its minimization is related to
the minimization of weighted graph cuts [Bertozzi and Flenner, 2012]. In this
sense, it parallels other techniques based on inference on graphs via
diffusion operators or function estimation [Coifman et al., 2005, Chung, 1997,
Zhou and Schölkopf, 2004, Szlam et al., 2008, Wang et al., 2008, Bühler and
Hein, 2009, Szlam and Bresson, 2010, Hein and Setzer, 2011].

Multiclass segmentation methods that cast the problem as a series of binary
classification problems use a number of different strategies: (i) deal
directly with some binary coding or indicator for the labels [Dietterich and
Bakiri, 1995, Wang et al., 2008], (ii) build a hierarchy or combination of
classifiers based on the one-vs-all approach or on class rankings [Hastie and
Tibshirani, 1998, Har-Peled et al., 2003] or (iii) apply a recursive
partitioning scheme consisting of successively subdividing clusters, until the
desired number of classes is reached [Szlam and Bresson, 2010, Hein and
Setzer, 2011]. While there are advantages to these approaches, such as
possible robustness to mislabeled data, there can be a considerable number of
classifiers to compute, and performance is affected by the number of classes
to partition.

In contrast, we propose an extension of the diffuse interface model that
obtains a simultaneous segmentation into multiple classes. The multiclass
extension is built by modifying the GL energy functional to remove the
prejudicial effect that the order of the labelings, given by integer values,
has in the smoothing term of the original binary diffuse interface model. A
new term that promotes homogenization in a multiclass setup is introduced. The
expression penalizes data points that are located close in the graph but are
not assigned to the same class. This penalty is applied independently of how
different the integer values are, representing the class labels. In this way,
the characteristics of the multiclass classification task are incorporated
directly into the energy functional, with a measure of smoothness independent
of label order, allowing us to obtain high-quality results. Alternative
multiclass methods minimize a Kullback-Leibler divergence function [Subramanya
and Bilmes, 2011] or expressions involving the discrete Laplace operator on
graphs [Zhou et al., 2004, Wang et al., 2008].

This paper is organized as follows. Section 2 reviews the diffuse interface
model for binary classification, and describes its application to semi-
supervised learning. Section 3 discusses our proposed multiclass extension and
the corresponding computational algorithm. Section 4 presents results obtained
with our method. Finally, section 5 draws conclusions and delineates future
work.

##  2 DATA SEGMENTATION WITH THE GINZBURG-LANDAU MODEL

The diffuse interface model [Bertozzi and Flenner, 2012] is based on a
continuous approach, using the Ginzburg-Landau (GL) energy functional to
measure the quality of data segmentation. A good segmentation is characterized
by a state with small energy. Let $u(\boldsymbol{x})$ be a scalar field
defined over a space of arbitrary dimensionality, and representing the state
of the system. The GL energy is written as the functional

$E_{GL}(u)=\frac{\epsilon}{2}\int\\!|\nabla
u|^{2}\;d\boldsymbol{x}+\frac{1}{\epsilon}\int\\!F(u)\;d\boldsymbol{x},$ (1)

with $\nabla$ denoting the spatial gradient operator, $\epsilon>0$ a real
constant value, and $F$ a double well potential with minima at $\pm 1$:

$F(u)=\frac{1}{4}\left(u^{2}-1\right)^{2}.$ (2)

Segmentation requires minimizing the GL functional. The norm of the gradient
is a smoothing term that penalizes variations in the field $u$. The potential
term, on the other hand, compels $u$ to adopt the discrete labels of $+1$ or
$-1$, clustering the state of the system around two classes. Jointly
minimizing these two terms pushes the system domain towards homogeneous
regions with values close to the minima of the double well potential, making
the model appropriate for binary segmentation.

The smoothing term and potential term are in conflict at the interface between
the two regions, with the first term favoring a gradual transition, and the
second term penalizing deviations from the discrete labels. A compromise
between these conflicting goals is established via the constant $\epsilon$. A
small value of $\epsilon$ denotes a small length transition and a sharper
interface, while a large $\epsilon$ weights the gradient norm more, leading to
a slower transition. The result is a diffuse interface between regions, with
sharpness regulated by $\epsilon$.

It can be shown that in the limit $\epsilon\to 0$ this function approximates
the total variation (TV) formulation in the sense of functional ($\Gamma$)
convergence [Kohn and Sternberg, 1989], producing piecewise constant solutions
but with greater computational efficiency than conventional TV minimization
methods. Thus, the diffuse interface model provides a framework to compute
piecewise constant functions with diffuse transitions, approaching the ideal
of the TV formulation, but with the advantage that the smooth energy
functional is more tractable numerically and can be minimized by simple
numerical methods such as gradient descent.

The GL energy has been used to approximate the TV norm for image segmentation
[Bertozzi and Flenner, 2012] and image inpainting [Bertozzi et al., 2007,
Dobrosotskaya and Bertozzi, 2008]. Furthermore, a calculus on graphs
equivalent to TV has been introduced in [Gilboa and Osher, 2008, Szlam and
Bresson, 2010].

### Application of Diffuse Interface Models to Graphs

An undirected, weighted neighborhood graph is used to represent the local
relationships in the data set. This is a common technique to segment classes
that are not linearly separable. In the $N$-neighborhood graph model, each
vertex $z_{i}\in Z$ of the graph corresponds to a data point with feature
vector $\boldsymbol{x}_{i}$, while the weight $w_{ij}$ is a measure of
similarity between $z_{i}$ and $z_{j}$. Moreover, it satisfies the symmetry
property $w_{ij}=w_{ji}$. The neighborhood is defined as the set of $N$
closest points in the feature space. Accordingly, edges exist between each
vertex and the vertices of its $N$-nearest neighbors. Following the approach
of [Bertozzi and Flenner, 2012], we calculate weights using the local scaling
of Zelnik-Manor and Perona [Zelnik-Manor and Perona, 2005],

$w_{ij}=\exp\left(-\frac{||\boldsymbol{x}_{i}-\boldsymbol{x}_{j}||^{2}}{\tau(\boldsymbol{x}_{i})\;\tau(\boldsymbol{x}_{j})}\right).$
(3)

Here, $\tau(\boldsymbol{x}_{i})=||\boldsymbol{x}_{i}-\boldsymbol{x}^{M}_{i}||$
defines a local value for each $\boldsymbol{x}_{i}$, where
$\boldsymbol{x}^{M}_{i}$ is the position of the $M$th closest data point to
$\boldsymbol{x}_{i}$, and $M$ is a global parameter.

It is convenient to express calculations on graphs via the graph Laplacian
matrix, denoted by $\boldsymbol{L}$. The procedure we use to build the graph
Laplacian is as follows.

  1. 1.

Compute the similarity matrix $\boldsymbol{W}$ with components $w_{ij}$
defined in (3). As the neighborhood relationship is not symmetric, the
resulting matrix $\boldsymbol{W}$ is also not symmetric. Make it a symmetric
matrix by connecting vertices $z_{i}$ and $z_{j}$ if $z_{i}$ is among the
$N$-nearest neighbors of $z_{j}$ or if $z_{j}$ is among the $N$-nearest
neighbors of $z_{i}$ [von Luxburg, 2006].

  2. 2.

Define $\boldsymbol{D}$ as a diagonal matrix whose $i$th diagonal element
represents the degree of the vertex $z_{i}$, evaluated as

$d_{i}=\sum_{j}w_{ij}.$ (4)

  3. 3.

Calculate the graph Laplacian: $\boldsymbol{L}=\boldsymbol{D}-\boldsymbol{W}$.

Generally, the graph Laplacian is normalized to guarantee spectral convergence
in the limit of large sample size [von Luxburg, 2006]. The symmetric
normalized graph Laplacian $\boldsymbol{L_{s}}$ is defined as

$\boldsymbol{L_{s}}=\boldsymbol{D}^{-1/2}\;\boldsymbol{L}\;\boldsymbol{D}^{-1/2}=\boldsymbol{I}-\boldsymbol{D}^{-1/2}\;\boldsymbol{W}\;\boldsymbol{D}^{-1/2}.$
(5)

Data segmentation can now be carried out through a graph-based formulation of
the GL energy. To implement this task, a fidelity term is added to the
functional as initially suggested in [Dobrosotskaya and Bertozzi, 2010]. This
enables the specification of a priori information in the system, for example
the known labels of certain points in the data set. This kind of setup is
called semi-supervised learning (SSL). The discrete GL energy for SSL on
graphs can be written as [Bertozzi and Flenner, 2012]:

$\displaystyle E_{GL_{\mathrm{SSL}}}(\boldsymbol{u})$ $\displaystyle=$
$\displaystyle\frac{\epsilon}{2}\langle\boldsymbol{u},\boldsymbol{L_{s}}\boldsymbol{u}\rangle+\frac{1}{\epsilon}\sum_{z_{i}\in
Z}F(u(z_{i}))$ (6) $\displaystyle+\sum_{z_{i}\in
Z}\frac{\lambda(z_{i})}{2}\;\left(u(z_{i})-u_{0}(z_{i})\right)^{2}$

In the discrete formulation, $\boldsymbol{u}$ is a vector whose component
$u(z_{i})$ represents the state of the vertex $z_{i}$, $\epsilon>0$ is a real
constant characterizing the smoothness of the transition between classes, and
$\lambda(z_{i})$ is a fidelity weight taking value $\lambda>0$ if the label
$u_{0}(z_{i})$ (i.e. class) of the data point associated with vertex $z_{i}$
is known beforehand, or $\lambda(z_{i})=0$ if it is not known (semi-
supervised).

Equation (6) may be understood as an example of the more general form of an
energy functional for data classification,

$\displaystyle
E(\boldsymbol{u})=||\boldsymbol{u}||_{a}+\frac{\lambda}{2}||\boldsymbol{u}-\boldsymbol{f}||_{b}^{p},$
(7)

where the norm $||u||_{a}$ is a regularization term and $||u-f||_{b}$ is a
fidelity term. The choice of the regularization norm $||\cdot||_{a}$ has non-
trivial consequences in the final classification accuracy. Attractive
qualities of the norm $||\cdot||_{a}$ include allowing classes to be close in
a metric space, and obtain segmentations for nonlinearly separable data. Both
of these goals are addressed using the GL energy functional for SSL.

Minimizing the functional simulates a diffusion process on the graph. The
information of the few labels known is propagated through the discrete
structure by means of the smoothing term, while the potential term clusters
the vertices around the states $\pm 1$ and the fidelity term enforces the
known labels. The energy minimization process itself attempts to reduce the
interface regions. Note that in the absence of the fidelity term, the process
could lead to a trivial steady-state solution of the diffusion equation, with
all data points assigned the same label.

The final state $u(z_{i})$ of each vertex is obtained by thresholding, and the
resulting homogeneous regions with labels of $+1$ and $-1$ constitute the two-
class data segmentation.

##  3 MULTICLASS EXTENSION

The double-well potential in the diffuse interface model for SSL flows the
state of the system towards two definite labels. Multiple-class segmentation
requires a more general potential function $F(u)$ that allows clusters around
more than two labels. For this purpose, we use the periodic-well potential
suggested by Li and Kim [Li and Kim, 2011],

$F(u)=\frac{1}{2}\,\\{u\\}^{2}\,(\\{u\\}-1)^{2},$ (8)

where $\\{u\\}$ denotes the fractional part of $u$,

$\\{u\\}=u-\lfloor u\rfloor,$ (9)

and $\lfloor u\rfloor$ is the largest integer not greater than $u$.

This periodic potential well promotes a multiclass solution, but the graph
Laplacian term in Equation (6) also requires modification for effective
calculations due to the fixed ordering of class labels in the multiple class
setting. The graph Laplacian term penalizes large changes in the spatial
distribution of the system state more than smaller gradual changes. In a
multiclass framework, this implies that the penalty for two spatially
contiguous classes with different labels may vary according to the (arbitrary)
ordering of the labels.

This phenomenon is shown in Figure 1. Suppose that the goal is to segment the
image into three classes: class 0 composed by the black region, class 1
composed by the gray region and class 2 composed by the white region. It is
clear that the horizontal interfaces comprise a jump of size 1 (analogous to a
two class segmentation) while the vertical interface implies a jump of size 2.
Accordingly, the smoothing term will assign a higher cost to the vertical
interface, even though from the point of view of the classification, there is
no specific reason for this. In this example, the problem cannot be solved
with a different label assignment. There will always be an interface with
higher costs than others independent of the integer values used.

Thus, the multiclass approach breaks the symmetry among classes, influencing
the diffuse interface evolution in an undesirable manner. Eliminating this
inconvenience requires restoring the symmetry, so that the difference between
two classes is always the same, regardless of their labels. This objective is
achieved by introducing a new class difference measure.

Figure 1: Three class segmentation. Black: class 0. Gray: class 1. White:
class 2.

###  3.1 Generalized Difference Function

The final class labels are determined by thresholding each vertex $u(z_{i})$,
with the label $y_{i}$ set to the nearest integer:

$y_{i}=\left\lfloor u(z_{i})+\frac{1}{2}\right\rfloor.$ (10)

The boundaries between classes then occur at half-integer values corresponding
to the unstable equilibrium states of the potential well. Define the function
$\hat{r}(x)$ to represent the distance to the nearest half-integer:

$\hat{r}(x)=\left|\frac{1}{2}-\\{x\\}\right|.$ (11)

A schematic of $\hat{r}(x)$ is depicted in Figure 2. The $\hat{r}(x)$ function
is used to define a generalized difference function between classes that
restores symmetry in the energy functional. Define the generalized difference
function $\rho$ as:

$\rho(u(z_{i}),u(z_{j}))=\left\\{\begin{array}[]{lll}\hat{r}(u(z_{i}))+\hat{r}(u(z_{j}))&&y_{i}\neq
y_{j}\\\ &&\\\
\left|\hat{r}(u(z_{i}))-\hat{r}(u(z_{j}))\right|&&y_{i}=y_{j}\end{array}\right.$
(12)

Thus, if the vertices are in different classes, the difference $\hat{r}(x)$
between each state’s value and the nearest half-integer is added, whereas if
they are in the same class, these differences are subtracted. The function
$\rho(x,y)$ corresponds to the tree distance (see Fig. 2). Strictly speaking,
$\rho$ is not a metric since it does not satisfy $\rho(x,y)=0\Rightarrow x=y$.
Nevertheless, the cost of interfaces between classes becomes the same
regardless of class labeling when this generalized distance function is
implemented.

Half-integerInteger$\hat{r}(x)$ Figure 2: Schematic interpretation of
generalized difference: $\hat{r}(x)$ measures distance to nearest half-
integer, and $\rho$ then corresponds to distance on tree.

The GL energy functional for SSL, using the new generalized difference
function $\rho$, is expressed as

$\displaystyle E_{MGL_{\mathrm{SSL}}}(\boldsymbol{u})$ $\displaystyle=$
$\displaystyle\frac{\epsilon}{2}\sum_{z_{i}\in Z}\sum_{z_{j}\in
Z}\frac{w_{ij}}{\sqrt{d_{i}d_{j}}}\,\left[\rho(u(z_{i}),u(z_{j}))\,\right]^{2}$
(13) $\displaystyle+\frac{1}{2\epsilon}\sum_{z_{i}\in
Z}\\{u(z_{i})\\}^{2}\,(\\{u(z_{i})\\}-1)^{2}$ $\displaystyle+\sum_{z_{i}\in
Z}\frac{\lambda(z_{i})}{2}\;\left(u(z_{i})-u_{0}(z_{i})\right)^{2}.$

Note that $\rho$ could also be used in the fidelity term, but for simplicity
this modification is not included. In practice, this has little effect on the
results.

###  3.2 Computational Algorithm

The GL energy functional given by (13) may be minimized iteratively, using
gradient descent:

$u_{i}^{m+1}=u_{i}^{m}-dt\,\left[\frac{\delta E_{MGL_{\mathrm{SSL}}}}{\delta
u_{i}}\right],$ (14)

where $u_{i}$ is a shorthand for $u(z_{i})$, $dt$ represents the time step and
the gradient direction is given by:

$\frac{\delta E_{MGL_{\mathrm{SSL}}}}{\delta u_{i}}=\epsilon
G(u_{i}^{m})+\frac{1}{\epsilon}F^{\prime}(u_{i}^{m})+\lambda_{i}\left(u_{i}^{m}-{u_{i}}_{0}\right)$
(15)

$G(u_{i}^{m})=\sum_{j}\frac{w_{ij}}{\sqrt{d_{i}d_{j}}}\left[\hat{r}(u_{i}^{m})\pm\hat{r}(u_{j}^{m})\right]\hat{r}^{\prime}(u_{i}^{m})$
(16)

$F^{\prime}(u_{i}^{m})=2\;\\{u_{i}^{m}\\}^{3}-3\;\\{u_{i}^{m}\\}^{2}+\\{u_{i}^{m}\\}$
(17)

Algorithm 1 Calculate $\boldsymbol{u}$

0: $\epsilon>0,dt>0,m_{\mathrm{max}}>0,K\mathrm{~{}given}$

0: $\mathrm{out}=\boldsymbol{u}^{m_{\mathrm{max}}}$

$\boldsymbol{u}^{0}\leftarrow rand((0,K))-\frac{1}{2},\ m\leftarrow 0$

for $m<m_{\mathrm{max}}$ do

$i\leftarrow 0$

for $i<n$ do

$u_{i}^{m+1}\leftarrow
u_{i}^{m}-dt\left(\epsilon\>G(u_{i}^{m})+\frac{1}{\epsilon}\>F^{\prime}(u_{i}^{m})+\lambda_{i}\left(u_{i}^{m}-{u_{i}}_{0}\right)\right)$

if $\mathrm{Label}(u_{i}^{m+1})\neq\mathrm{Label}(u_{i}^{m})$ then

$(v_{i})_{k}\leftarrow k+\\{u_{i}^{m+1}\\}$

$u_{i}^{m+1}\leftarrow(v_{i})_{k}\mathrm{~{}where~{}}k=\arg\min_{\;0\leq
k<K}\;\sum_{j}\frac{w_{ij}}{\sqrt{d_{i}d_{j}}}\,\left[\rho((v_{i})_{k},u_{j})\,\right]^{2}$

end if

$i\leftarrow i+1$

end for

$m\leftarrow m+1$

end for

The gradient of the generalized difference function $\rho$ is not defined at
half integer values. Hence, we modify the method using a greedy strategy:
after detecting that a vertex changes class, the new class that minimizes the
smoothing term is selected, and the fractional part of the state computed by
the gradient descent update is preserved. Consequently, the new state of
vertex $i$ is the result of gradient descent, but if this causes a change in
class, then a new state is determined.

Specifically, let $k$ represent an integer in the range of the problem, i.e.
$k\in[0,K-1]$, where $K$ is the number of classes in the problem. Given the
fractional part $\\{u\\}$ resulting from the gradient descent update, define
$(v_{i})_{k}=k+\\{u_{i}\\}$. Find the integer $k$ that minimizes
$\sum_{j}\frac{w_{ij}}{\sqrt{d_{i}d_{j}}}\,\left[\rho((v_{i})_{k},u_{j})\,\right]^{2}$,
the smoothing term in the energy functional, and use $(v_{i})_{k}$ as the new
vertex state. A summary of the procedure is shown in Algorithm 1 with
$m_{\mathrm{max}}$ denoting the maximum number of iterations.

##  4 RESULTS

The performance of the multiclass diffuse interface model is evaluated using a
number of data sets from the literature, with differing characteristics. Data
and image segmentation problems are considered on synthetic and real data
sets.

###  4.1 Synthetic Data

A synthetic three-class segmentation problem is constructed following an
analogous procedure used in [Bühler and Hein, 2009] for “two moon” binary
classification, using three half circles (“three moons”). The half circles are
generated in $\mathbb{R}^{2}$. The two top circles have radius $1$ and are
centered at $(0,0)$ and $(3,0)$. The bottom half circle has radius $1.5$ and
is centered at $(1.5,0.4)$. We sample 1500 data points (500 from each of these
half circles) and embed them in $\mathbb{R}^{100}$. The embedding is completed
by adding Gaussian noise with $\sigma^{2}=0.02$ to each of the 100 components
for each data point. The dimensionality of the data set, together with the
noise, make this a nontrivial problem.

Figure 3: Three-class segmentation. Left: spectral clustering. Right:
multiclass GL (adaptive $\epsilon$).

The difficulty of the problem is illustrated in Figure 3, where we use both
spectral clustering decomposition and the multiclass GL method. The same graph
structure is used for both methods. The symmetric graph Laplacian is computed
based on edge weights given by (3), using $N=10$ nearest neighbors and local
scaling based on the $M=10$ closest point. The spectral clustering results are
obtained by applying a $k$-means algorithm to the first $3$ eigenvectors of
the symmetric graph Laplacian. The average error obtained, over 100 executions
of spectral clustering, is 20% ($\pm 0.6\%$). The figure displays the best
result obtained, corresponding to an error of $18.67\%$.

The multiclass GL method was implemented with the following parameters:
interface scale $\epsilon=1$, step size $dt=0.01$ and number of iterations
$m_{\mathrm{max}}=800$. The fidelity term is determined by labeling 25 points
randomly selected from each class (5% of all points), and setting the fidelity
weight to $\lambda=30$ for those points. Several runs of the procedure are
performed to isolate effects from the random initialization and the arbitrary
selection of fidelity points. The average error obtained, over 100 runs with
four different fidelity sets, is 5.2% ($\pm 1.01\%$). In general terms, the
system evolves from an initially inhomogeneous state, rapidly developing small
islands around fidelity points that become seeds for homogeneous regions and
progressing to a configuration of classes forming nearly uniform clusters.

The multiclass results were further improved by incrementally decreasing
$\epsilon$ to allow sharper transitions between states as in [Bertozzi and
Flenner, 2012]. With this approach, the average error obtained over 100 runs
is reduced to 2.6% ($\pm 0.3\%$). The best result obtained in these runs is
displayed in Figure 3 and corresponds to an average error of 2.13%. In these
runs, $\epsilon$ is reduced from $\epsilon_{0}=2$ to $\epsilon_{f}=0.1$ in
decrements of 10%, with $40$ iterations performed per step. The average
computing time per run in this adaptive technique is 1.53s in an Intel Quad-
Core @ 2.4 GHz, without any parallel processing.

For comparison, we note the results from the literature for the simpler two
moon problem ($\mathbb{R}^{100}$, $\sigma^{2}=0.02$ noise). The best errors
reported include: 6% for p-Laplacian [Bühler and Hein, 2009], 4.6% for ratio-
minimization relaxed Cheeger cut [Szlam and Bresson, 2010], and 2.3% for
binary GL [Bertozzi and Flenner, 2012]. While these are not SSL methods the
last of these does involve other prior information in the form of a mass
balance constraint. It can be seen that both of our procedures, fixed and
adaptive $\epsilon$, produce high-quality results even for the more complex
three-class segmentation problem. Calculation times are also competitive with
those reported for the binary case (0.5s - 50s).

###  4.2 Image Segmentation

As another test setup, we use a grayscale image of size $191\times 196$, taken
from [Jung et al., 2007, Li and Kim, 2011] and composed of 5 classes: black,
dark gray, medium gray, light gray and white. This image contains structure,
such as an internal hole and junctions where multiple classes meet. The image
information is represented through feature vectors defined as
$(x_{i},y_{i},\mathrm{pix}_{i})$, with $x_{i}$ and $y_{i}$ corresponding to
$(x,y)$ coordinates of the pixel and $\mathrm{pix}_{i}$ equal to the intensity
of the pixel. All of these are normalized so as to obtain values in the range
$[0,1]$.

The graph is constructed using $N=30$ nearest neighbors and local scaling
based on the $M=30$ closest point. We use parameters $\epsilon=1$, $dt=0.01$
and $m_{\mathrm{max}}=800$. We then choose 1500 random points (4% of the
total) for the fidelity term, with $\lambda=30$. Figure 4 displays the
original image with the randomly selected fidelity points (top left), and the
five-class segmentation. Each class image shows in white the pixels identified
as belonging to the class, and in black the pixels of the other classes. In
this case, all the classes are segmented perfectly with an average run time of
59.7s. The method of Li and Kim [Li and Kim, 2011] also segments this image
perfectly, with a reported run time of 0.625s. However, their approach uses
additional information, including a pre-assignment of specific grayscale
levels to classes, and the overall densities of each class. Our approach does
not require these.

Figure 4: Image Segmentation Results. Top left: Original five-class image,
with randomly chosen fidelity points displayed. Other panels: the five
segmented classes, shown in white.

###  4.3 MNIST Data

The MNIST data set available at http://yann.lecun.com/exdb/mnist/ is composed
of 70,000 images of size $28\times 28$, corresponding to a broad sample of
handwritten digits 0 through 9. We use the multiclass diffuse interface model
to segment the data set automatically into 10 classes, one per handwritten
digit. Before constructing the graph, we preprocess the data by normalizing
and projecting into 50 principal components, following the approach in [Szlam
and Bresson, 2010]. No further steps, such as smoothing convolutions, are
required. The graph is computed with $N=10$ nearest neighbors and local
scaling based on the $M=10$ closest points.

An adaptive $\epsilon$ variant of the algorithm is implemented, with
parameters $\epsilon_{0}=2$, $\epsilon_{f}=0.01$, $\epsilon$ decrement 10%,
$dt=0.01$, and 40 iterations per step. For the fidelity term, 7,000 images
(10% of total) are chosen, with weight $\lambda=30$. The average error
obtained, over 20 runs with four different fidelity sets, is 7% ($\pm
0.072\%$). The confusion matrix for the best result obtained, corresponding to
a 6.86% error, is given in Table 1: each row represents the segmentation
obtained, while the columns represent the true digit labels. For reference,
the average computing time per run in this adaptive technique is 132s. Note
that, in the segmentations, the largest mistakes made are in trying to
distinguish digits 4 from 9 and 7 from 9.

For comparison, errors reported using unsupervised clustering algorithms in
the literature are: 12.9% for p-Laplacian [Bühler and Hein, 2009], 11.8% for
ratio-minimization relaxed Cheeger cut [Szlam and Bresson, 2010], and 12.36%
for the multicut version of the normalized 1-cut [Hein and Setzer, 2011]. A
more sophisticated graph-based diffusion method applied in a semi-supervised
setup (transductive classification), with function-adapted eigenfunctions, a
graph constructed with 13 neighbors, and self-tuning with the 9th neighbor
reported in [Szlam et al., 2008] obtains an error of 7.4%. Results with
similar errors are reported in [Liu et al., 2010]. Thus, the performance of
the multiclass GL on this data set improves upon other published results,
while requiring less preprocessing and a simpler regularization of the
functions on the graph.

Table 1: Confusion Matrix for the MNIST Data Segmentation. Obtained / True | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9
---|---|---|---|---|---|---|---|---|---|---
0 | 6712 | 3 | 39 | 10 | 6 | 36 | 57 | 10 | 61 | 28
1 | 1 | 7738 | 7 | 15 | 9 | 1 | 9 | 23 | 36 | 12
2 | 24 | 50 | 6632 | 95 | 65 | 17 | 16 | 63 | 65 | 30
3 | 13 | 16 | 84 | 6585 | 8 | 218 | 5 | 42 | 153 | 84
4 | 5 | 6 | 27 | 8 | 6279 | 32 | 13 | 59 | 43 | 305
5 | 21 | 6 | 13 | 128 | 27 | 5736 | 57 | 3 | 262 | 34
6 | 91 | 26 | 50 | 11 | 35 | 91 | 6693 | 0 | 45 | 1
7 | 6 | 6 | 31 | 97 | 26 | 15 | 0 | 6689 | 24 | 331
8 | 27 | 15 | 86 | 156 | 21 | 110 | 25 | 16 | 6065 | 66
9 | 3 | 11 | 21 | 36 | 348 | 57 | 1 | 388 | 71 | 6067

##  5 CONCLUSIONS

We have proposed a new multiclass segmentation procedure, based on the diffuse
interface model. The method obtains segmentations of several classes
simultaneously without using one-vs-all or alternative sequences of binary
segmentations required by other multiclass methods. The local scaling method
of Zelnik-Manor and Perona, used to construct the graph, constitutes a useful
representation of the characteristics of the data set and is adequate to deal
with high-dimensional data.

Our modified diffusion method, represented by the non-linear smoothing term
introduced in the Ginzburg-Landau functional, exploits the structure of the
multiclass model and is not affected by the ordering of class labels. It
efficiently propagates class information that is known beforehand, as
evidenced by the small proportion of fidelity points (4% \- 10% of dataset)
needed to perform accurate segmentations. Moreover, the method is robust to
initial conditions. As long as the initialization represents all classes
uniformly, different initial random configurations produce very similar
results. The main limitation of the method appears to be that fidelity points
must be representative of class distribution. As long as this holds, such as
in the examples discussed, the long-time behavior of the solution relies less
on choosing the “right” initial conditions than do other learning techniques
on graphs.

State-of-the-art results with small classification errors were obtained for
all classification tasks. Furthermore, the results do not depend on the
particular class label assignments. Future work includes investigating the
diffuse interface parameter $\epsilon$. We conjecture that the proposed
functional converges (in the $\Gamma$-convergence sense) to a total
variational type functional on graphs as $\epsilon$ approaches zero, but the
exact nature of the limiting functional is unknown.

## ACKNOWLEDGEMENTS

This research has been supported by the Air Force Office of Scientific
Research MURI grant FA9550-10-1-0569 and by ONR grant N0001411AF00002.

## REFERENCES

  * [Allwein et al., 2000] Allwein, E. L., Schapire, R. E., and Singer, Y. (2000).  Reducing multiclass to binary: A unifying approach for margin classifiers.  Journal of Machine Learning Research, 1:113–141.
  * [Bertozzi et al., 2007] Bertozzi, A., Esedoḡlu, S., and Gillette, A. (2007).  Inpainting of binary images using the Cahn-Hilliard equation.  IEEE Transactions on Image Processing, 16(1):285–291.
  * [Bertozzi and Flenner, 2012] Bertozzi, A. L. and Flenner, A. (2012).  Diffuse interface models on graphs for classification of high dimensional data.  Multiscale Modeling and Simulation, 10(3):1090–1118.
  * [Bühler and Hein, 2009] Bühler, T. and Hein, M. (2009).  Spectral clustering based on the graph $p$-Laplacian.  In Bottou, L. and Littman, M., editors, Proceedings of the 26th International Conference on Machine Learning, pages 81–88. Omnipress, Montreal, Canada.
  * [Chung, 1997] Chung, F. R. K. (1997).  Spectral graph theory.  In Regional Conference Series in Mathematics, volume 92. Conference Board of the Mathematical Sciences (CBMS), Washington, DC.
  * [Coifman et al., 2005] Coifman, R. R., Lafon, S., Lee, A. B., Maggioni, M., Nadler, B., Warner, F., and Zucker, S. W. (2005).  Geometric diffusions as a tool for harmonic analysis and structure definition of data: Diffusion maps.  Proceedings of the National Academy of Sciences, 102(21):7426–7431.
  * [Dietterich and Bakiri, 1995] Dietterich, T. G. and Bakiri, G. (1995).  Solving multiclass learning problems via error-correcting output codes.  Journal of Artificial Intelligence Research, 2(1):263–286.
  * [Dobrosotskaya and Bertozzi, 2008] Dobrosotskaya, J. A. and Bertozzi, A. L. (2008).  A wavelet-laplace variational technique for image deconvolution and inpainting.  IEEE Trans. Image Process., 17(5):657–663.
  * [Dobrosotskaya and Bertozzi, 2010] Dobrosotskaya, J. A. and Bertozzi, A. L. (2010).  Wavelet analogue of the Ginzburg-Landau energy and its gamma-convergence.  Interfaces and Free Boundaries, 12(2):497–525.
  * [Gilboa and Osher, 2008] Gilboa, G. and Osher, S. (2008).  Nonlocal operators with applications to image processing.  Multiscale Modeling and Simulation, 7(3):1005–1028.
  * [Har-Peled et al., 2003] Har-Peled, S., Roth, D., and Zimak, D. (2003).  Constraint classification for multiclass classification and ranking.  In S. Becker, S. T. and Obermayer, K., editors, Advances in Neural Information Processing Systems 15, pages 785–792. MIT Press, Cambridge, MA.
  * [Hastie and Tibshirani, 1998] Hastie, T. and Tibshirani, R. (1998).  Classification by pairwise coupling.  In Advances in Neural Information Processing Systems 10. MIT Press, Cambridge, MA.
  * [Hein and Setzer, 2011] Hein, M. and Setzer, S. (2011).  Beyond spectral clustering - tight relaxations of balanced graph cuts.  In Shawe-Taylor, J., Zemel, R., Bartlett, P., Pereira, F., and Weinberger, K., editors, Advances in Neural Information Processing Systems 24, pages 2366–2374.
  * [Jung et al., 2007] Jung, Y. M., Kang, S. H., and Shen, J. (2007).  Multiphase image segmentation via Modica-Mortola phase transition.  SIAM J. Appl. Math, 67(5):1213–1232.
  * [Kohn and Sternberg, 1989] Kohn, R. V. and Sternberg, P. (1989).  Local minimizers and singular perturbations.  Proc. Roy. Soc. Edinburgh Sect. A, 111(1-2):69–84.
  * [Li and Kim, 2011] Li, Y. and Kim, J. (2011).  Multiphase image segmentation using a phase-field model.  Computers and Mathematics with Applications, 62:737–745.
  * [Liu et al., 2010] Liu, W., He, J., and Chang, S.-F. (2010).  Large graph construction for scalable semi-supervised learning.  Proceedings of the 27th International Conference on Machine Learning.
  * [Subramanya and Bilmes, 2011] Subramanya, A. and Bilmes, J. (2011).  Semi-supervised learning with measure propagation.  Journal of Machine Learning Research, 12:3311–3370.
  * [Szlam and Bresson, 2010] Szlam, A. and Bresson, X. (2010).  Total variation and cheeger cuts.  In Fürnkranz, J. and Joachims, T., editors, Proceedings of the 27th International Conference on Machine Learning, pages 1039–1046. Omnipress, Haifa, Israel.
  * [Szlam et al., 2008] Szlam, A. D., Maggioni, M., and Coifman, R. R. (2008).  Regularization on graphs with function-adapted diffusion processes.  Journal of Machine Learning Research, 9:1711–1739.
  * [von Luxburg, 2006] von Luxburg, U. (2006).  A tutorial on spectral clustering.  Technical Report TR-149, Max Planck Institute for Biological Cybernetics.
  * [Wang et al., 2008] Wang, J., Jebara, T., and Chang, S.-F. (2008).  Graph transduction via alternating minimization.  Proceedings of the 25th International Conference on Machine Learning.
  * [Zelnik-Manor and Perona, 2005] Zelnik-Manor, L. and Perona, P. (2005).  Self-tuning spectral clustering.  In Saul, L. K., Weiss, Y., and Bottou, L., editors, Advances in Neural Information Processing Systems 17. MIT Press, Cambridge, MA.
  * [Zhou et al., 2004] Zhou, D., Bousquet, O., Lal, T. N., Weston, J., and Schölkopf, B. (2004).  Learning with local and global consistency.  In Thrun, S., Saul, L. K., and Schölkopf, B., editors, Advances in Neural Information Processing Systems 16, pages 321–328. MIT Press, Cambridge, MA.
  * [Zhou and Schölkopf, 2004] Zhou, D. and Schölkopf, B. (2004).  A regularization framework for learning from graph data.  In Workshop on Statistical Relational Learning. International Conference on Machine Learning, Banff, Canada.