# LiDAR-assisted Large-scale Privacy Protection in Street-view
Cycloramas

Clint Sebastian, Bas Boom, Egor Bondarev, Peter H.N. de With;
Eindhoven University of Technology, Eindhoven, The Netherlands;
Cyclomedia B.V, Zaltbommel, The Netherlands

###### Abstract

Recently, privacy has a growing importance in several domains, especially in
street-view images. The conventional way to achieve this is to automatically
detect and blur sensitive information from these images. However, the
processing cost of blurring increases with the ever growing resolution of
images. We propose a system that is cost-effective even after increasing the
resolution by a factor of 2.5. The new system utilizes depth data obtained
from LiDAR to significantly reduce the search space for detection, thereby
reducing the processing cost. Besides this, we test several detectors after
reducing the detection space and provide an alternative solution based on
state-of-the-art deep learning detectors to the existing HoG-SVM-Deep system
that is faster and has a higher performance.

## Introduction

Several street-view services such as Google Street View, Bing Maps Streetside,
Mapillary have systematically collected and hosted millions of images.
However, the privacy of an individual is vital and hence concealing this
information is important for privacy protection. Typically, this is achieved
by detecting the privacy-sensitive content and then blurring it. However, this
becomes increasingly challenging with growing image resolution. In a
commercial setting, it is preferred to exploit ultra-high resolution images
which are updated frequently. A higher resolution offers a clearer view of the
important content and allows to extract specific details, while the frequent
updates provide better tracking of possible changes. However, these benefits
come with higher production effort of images that incurs further processing
cost for privacy protection. In this research for street-view imaging, since
the commercial imaging vendor Cyclomedia transitions from a 100-Megapixel
(Mpix) to a 250-Mpix camera system with LiDAR, we present a privacy protection
solution, which offers better performance and similar processing cost to the
existing system.

Currently, Cyclomedia offers 100-Mpix 360$\degree$ panoramic images
(cycloramas), which are taken at a driving interval of every 5 meters. The new
system will have 250-Mpix cycloramas along with LiDAR data acquired at the
same capturing frequency. However, the current detection algorithm has a
quadratic computation cost that makes it infeasible for the processing of the
250-Mpix images. The proposed system exploits the LiDAR to preprocess the
250-Mpix cycloramas such that the search space is largely reduced for the
detection algorithm. To account for the quadratic cost of the current
detector, we replace it with a Convolutional Neural Network (CNN) detector
that has a linear cost for computation. We extensively test several detectors
to reach the best trade-off between speed and performance. Finally, the
proposed method achieves a higher detection performance and a computation cost
similar to the 100-Mpix system.

## Related Work

Various approaches are available for privacy protection in street-view images.
The typical way to achieve this is to detect and blur objects in images. One
of the first approach for privacy protection explores a two-stage detector
[1]. The first stage comprises of two detectors running in parallel, one with
high recall and another with high precision. This is followed by a post-
processing step that converts the output regions in the bounding box to a
feature vector that is passed to a neural network. If the output from the
neural network exceeds a certain threshold, then it is blurred. Apart from
detection and blurring methods, there are only a few approaches that
simultaneously detect, remove and inpaint pedestrians [2, 3]. However,
inpainting based on the context leads to inaccurate content and is not
preferred for commercial applications.

With the advances in deep learning, object detection has improved
significantly. Object detection provides a reliable way to localize and detect
objects of interest. Several deep learning-based object detection algorithms
are available with trade-offs between speed and performance. They are built on
top of CNN architectures such as ResNet and its variants [4, 5]. Each of the
CNN architectures has a meta-architecture on top that completes the object
detection framework. Popular meta-architectures include Faster-RCNN, Region-
based Fully Convolutional Network (RFCN) and Single Shot Multi-box detector
(SSD) [6, 7, 8]. Both Faster-RCNN and RFCN are two-stage detectors. Faster-
RCNN has a region proposal network (RPN) that utilizes predefined anchor boxes
which are later refined, whereas RFCN has a position-sensitive RoI pooling
that encodes class-specific positional information into specific feature maps.
SSD is a single-stage detector that directly predicts a bounding box and its
score at each spatial position using a small convolution layer. Compared to
both Faster-RCNN and RFCN, SSD is faster due to its single-stage approach.
However, the two-stage networks offer better performance. In this research, we
compare the performance of these meta-architectures on top of different CNN
architectures, in order to find a detector that is suited to replace the
existing face and license plate detectors.

Figure 1: Overview of the new detection system. The input RGB-D image is
converted to a point cloud and the planes are extracted using RANSAC. The
found horizontal planes are combined and are reprojected to an image plane.
The red region indicates the removed areas of the ego-vehicle, the yellow
region shows the extracted ground planes and the green region depicts the
buffer added to the top boundary of the extracted planes (bottom-right image).
The extracted areas (green and yellow regions) are processed by the object
detector.

Nowadays, LiDAR-based system are becoming more prevalent in the autonomous
driving industry. RGB-D-based object detection has been explored in [9, 10].
However, to the best of our knowledge, this is the first work that utilizes
depth obtained from LiDAR to reduce the search space of a detection algorithm.
Plane fitting in point clouds using Random Sample Consensus (RANSAC) is well
studied and several improvements have been proposed [11]. In this research, we
utilize depth data obtained from LiDAR to create a point cloud which is later
used to estimate planes. The estimated planes are further applied to reduce
the search space in the cycloramas.

## Overview of Architecture

The final objective of the ideal 250-Mpix system is to have a processing cost
(on the cloud) similar to that of the current 100-Mpix system. The current
system is deployed on CPU, whereas the proposed system utilizes the GPU
architecture. To provide a fair comparison, we compare the processing time
with a specific constraint. The constraint is that cost of processing on the
GPU is 6 times more expensive than the CPU. Hence, the new system should
process a single image in one-sixth of the current system’s processing time.
The emphasis is on achieving high recall since privacy protection is crucial.
False positives that may lower the precision can be manually removed.

### Camera system and dataset

Currently, the cycloramas (360$\degree$ panoramas) are captured at the
100-Mpix resolution. They are processed by the current detection system on
CPUs in the cloud. The new camera system captures 250-Mpix cycloramas along
with the LiDAR data. The LiDAR data is processed to produce a depth map
corresponding to every pixel in the RGB image. The dataset consist of
approximately 3000 cycloramas with 50,000 objects (faces and license plates).
The images are taken from various parts of the Netherlands, United States and
Germany. The dataset is divided in the ratio of 70:30 for training and test
sets.

### Current detection system

The current system at Cyclomedia is composed of an ensemble of eight cascaded
HoG-SVM detectors. Five of the detectors are trained to detect faces, whereas
three of them are used for the detection of license plates. Each of them is a
HoG feature-based SVM object detector that is stacked three times forming a
cascade, followed by a CNN classifier that is used to remove false positives.
The detectors try to achieve a high recall at the cost of higher false
positives. Each detector is also trained with hard negative mining to reduce
false positives. The detectors follow a multi-scale approach for training and
inference. Each image is processed over 50 scales to maximize the recall. For
simplicity, we refer to this complex pipeline as HoG-SVM-Deep. The
disadvantage of this system is that with growing resolution, the number of
scales that needs to be scanned for the objects of interest also increases.
Empirically, we have found that the cost of processing increases approximately
by 6.5 times in the transition from 100-Mpix to the 250-Mpix system.

### Proposed detection system

Figure 2: Precision recall curves for different detectors for face and license
plate detectors.

The proposed system is composed of two parts, a pre-processing step to reduce
the search space for the detection algorithm, and a deep learning detector
step which offers a good trade-off between speed and performance. An overview
of the proposed system is shown in Figure 1. A cyclorama may contain objects
of interest at different parts of an image. However, some regions such as sky
or the area covered by the ego-vehicle do not inherently contain any objects.
In the current system, these regions are removed prior to the detection and
only 66% of a cyclorama is processed by the detector. However, the regions
where the object of interest is present are dependent on the ground plane. We
expect objects such as faces and license plates to be present at most two
meters above the ground level. We utilize the depth maps (calibrated and
corrected) from the LIDAR scanners to estimate the ground planes. We first
project the depth map to a point cloud, then fit the planes using RANSAC.

Figure 3: Outputs from plane extraction algorithm. The red and the dark
regions of the image are not processed. The bands of yellow and green area is
processed by the detector.

Projecting a 250-Mpix image to a point cloud and fitting planes are
computationally expensive operations. Therefore, the RGB-D image is
downsampled prior to projecting to a point cloud. The image is downsampled by
a factor of 10 in both $x-$ and $y-$ directions. The planar segmentation
results in vertical and horizontal planes. We filter out the vertical planes
and retain the horizontal planes, as they may contain people or vehicles. The
horizontal surfaces are reprojected into the image plane, to obtain a mask
that possibly contains the objects of interest. We select a rectangular region
around the boundaries of the mask, with the region size at about one-third of
the cyclorama height.

The region obtained from the mask has high resolution that cannot be directly
processed by a CNN. Therefore, it is converted into patches of 1200 $\times$
600 pixels which is later supplied to the detection algorithm. Due to the
dynamic nature of the ground planar segmentation, there is a significant
reduction in the area for object detection. We have found that on average only
31% of the cyclorama needs to be processed by the detection algorithm. This is
less than one half of the manually estimated processing area in the current
detection system.

Table 1: Table 1: Time vs. Performance of 100-Mpix and 250-Mpix camera systems on face and license plate detectors. Detector | Time (sec) | Max Recall
---|---|---
current | new | Face | L.P.
SSD Inception v2 | 16.50 | 16.98 | 44.8 | 67.8
FRCNN Inception v2 | 21.15 | 21.54 | 79.6 | 89.0
RFCN ResNet-101 | 44.64 | 45.94 | 86.9 | 94.1
FRCNN ResNet-101 | 86.54 | 89.06 | 82.6 | 92.4
FRCNN In. ResNet v2 | 378.68 | 389.62 | 88.7 | 95.3
HoG-SVM-Deep | 215 | 879.12 | 57.8 | 75.2

## Experiments and Results

We have experimented with several deep learning detectors and compared them
with the currently deployed existing HoG-SVM-Deep system. We have tested deep
learning detectors with three meta-architectures, namely SSD, RFCN, and
Faster-RCNN (abbreviated as FRCNN in the table) with Inception v2, ResNet-101
and Inception ResNet v2, as the backbone architecture. We train on the above-
mentioned dataset with a model pretrained on the COCO dataset. We train each
network for 200,000 iterations with a batch size of 4 for all networks except
SSD Inception v2. The SSD Inception v2 network uses a batch size of 24 with a
learning rate of $1^{-05}$. The learning rates for the rest of the networks
are set to $1^{-06}$. The computational costs of the different detectors are
presented in Table 1. Note that the reported times are GPU execution times for
the deep learning-based detectors, whereas the HoG-SVM-Deep detector uses a
CPU. The HoG-SVM-Deep system is executed on the CPU as it is cheaper to
process on the cloud services, whereas the deep learning-based detectors are
more efficiently processed on the GPU. HoG-SVM-Deep is executed on a Xeon
E5-2673 v3 CPU, whereas the rest of the detectors are executed on a Tesla K80
GPU.

The plane extraction is performed using RANSAC. Three points are selected and
the corresponding parameters are computed. Depending on the given threshold,
the outliers are removed and a plane is fitted. This procedure is repeated
until a plane with a maximum number of points is obtained. Once a plane is
estimated, the points corresponding to that plane are removed and a new plane
is computed. This process is repeated until the top-10 planes are found.
Usually, the largest plane is the ground plane. However, we combine the
largest 10 planes in the output image after reprojection of the point cloud,
to maximize the area of ground surface. We have tuned the RANSAC-based plane
fitting to our dataset. In our experiments, we have found that setting the
distance threshold (distance at which a point is considered as an inlier for
the model) to 0.5 meters results in the best estimation of planes. After the
planes are estimated and reprojected to the image plane, we add a small buffer
of 350 pixels (green region in Figure 3) to the top boundary of the extracted
plane (yellow region in Figure 3) on the reprojected image. We have observed
that such a small buffer is needed to account for people standing on the
ground farther away from the image. It was empirically found that adding 350
pixels can account for all the objects of interest present in our dataset.

After the planar segmentation and subsequent patching, each patch of 1200
$\times$ 600 pixels is fed to the detector. The detector performance on our
dataset is presented in Figure 2 and Table 1. All the networks except SSD
Inception v2, offer better performance than HoG-SVM-Deep for both face and
license plate detection. Only Faster-RCNN with Inception ResNet v2 has higher
computation cost than HoG-SVM-Deep, but offers significantly higher
performance. Comparing the execution times of the current and new system of
all the detectors, it can be seen that they yield very similar computation
time. This is due to the search-space reduction achieved using depth data on
the new camera system. The processing cost of executing the current system on
the CPU is 215 seconds. Following the constraint mentioned in the previous
section, we need a system that utilizes one-sixth of the current system ( 36
seconds). From Table 1, only Faster-RCNN with Inception v2 satisfies that
requirement. However, it offers much higher performance than the current
system. On the other hand, RFCN with ResNet-101 have higher performance at the
expense of 28% longer processing time. Both options form feasible solutions
for the existing system.

### Limitations

We have utilized the RANSAC plane extraction algorithm, as it is effective and
fast. However, this may be challenging in locations where the surfaces are not
flat. The underlying assumption that objects such as faces can be found only
on the ground plane have exceptions in scenarios where people are present on
top of a building. We found that such cases occur once in every 50,000 images.
In such scenarios, the objects are manually blurred.

## Conclusion

We have presented a framework that is an alternative to the existing system
used at Cyclomedia for privacy protection in street-view images. The current
system utilizes HoG-SVM-Deep to process a 100-Mpix image. We propose a new
detection system for the 250-Mpix images that utilizes the depth data obtained
from LiDAR to reduce the search space. The reduction results in processing
only 31% of the 250-Mpix cyclorama. We obtain detectors that are suited for
our problem, namely Faster-RCNN with Inception ResNet v2 and RFCN with
ResNet-101. Both of these detectors offer a higher recall than the existing
system for detecting both face and license plates, providing 37-50% and 18-25%
improvement in recall on face and license plate detection. On the other hand,
the proposed system using the new detectors results in processing costs that
are comparable to the existing 100-Mpix camera system.

## References

  * [1] Frome, Andrea, German Cheung, Ahmad Abdulkader, Marco Zennaro, Bo Wu, Alessandro Bissacco, Hartwig Adam, Hartmut Neven, and Luc Vincent. ”Large-scale privacy protection in Google Street View.” In ICCV, pp. 2373-2380. 2009.
  * [2] Nodari, Angelo, Marco Vanetti, and Ignazio Gallo. ”Digital privacy: Replacing pedestrians from google street view images.” In Pattern Recognition (ICPR), 2012 21st International Conference on, pp. 2889-2893. IEEE, 2012.
  * [3] Flores, Arturo, and Serge Belongie. ”Removing pedestrians from google street view images.” In Computer Vision and Pattern Recognition Workshops (CVPRW), 2010 IEEE Computer Society Conference on, pp. 53-58. IEEE, 2010.
  * [4] He, Kaiming, 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, pp. 770-778. 2016.
  * [5] Szegedy, Christian, Sergey Ioffe, Vincent Vanhoucke, and Alexander A. Alemi. ”Inception-v4, inception-resnet and the impact of residual connections on learning.” In AAAI, vol. 4, p. 12. 2017.
  * [6] Ren, Shaoqing, Kaiming He, Ross Girshick, and Jian Sun. ”Faster R-CNN: Towards real-time object detection with region proposal networks.” In Advances in neural information processing systems, pp. 91-99. 2015.
  * [7] Dai, Jifeng, Yi Li, Kaiming He, and Jian Sun. ”R-fcn: Object detection via region-based fully convolutional networks.” In Advances in neural information processing systems, pp. 379-387. 2016.
  * [8] Liu, Wei, Dragomir Anguelov, Dumitru Erhan, Christian Szegedy, Scott Reed, Cheng-Yang Fu, and Alexander C. Berg. ”Ssd: Single shot multibox detector.” In European conference on computer vision, pp. 21-37. Springer, Cham, 2016.
  * [9] Rahman, Mohammad Muntasir, Yanhao Tan, Jian Xue, and Ke Lu. ”RGB-D object recognition with multimodal deep convolutional neural networks.” In Multimedia and Expo (ICME), 2017 IEEE International Conference on, pp. 991-996. IEEE, 2017.
  * [10] Schwarz, Max, Anton Milan, Arul Selvam Periyasamy, and Sven Behnke. ”RGB-D object detection and semantic segmentation for autonomous manipulation in clutter.” The International Journal of Robotics Research 37, no. 4-5 (2018): 437-451.
  * [11] Fischler, Martin A., and Robert C. Bolles. ”Random sample consensus: a paradigm for model fitting with applications to image analysis and automated cartography.” Communications of the ACM 24.6 (1981): 381-395


# Slow decay of infection in the inhomogeneous SIR model

Hidetsugu Sakaguchi and Yuta Nakao  Department of Applied Science for
Electronics and Materials, Interdisciplinary Graduate School of Engineering
Sciences, Kyushu University, Kasuga, Fukuoka 816-8580, Japan

###### Abstract

The SIR model with spatially inhomogeneous infection rate is studied with
numerical simulations in one, two, and three dimensions, considering the case
that the infection spreads inhomogeneously in densely populated regions or hot
spots. We find that the total population of infection decays very slowly in
the inhomogeneous systems in some cases, in contrast to the exponential decay
of the infected population $I(t)$ in the SIR model of the ordinary
differential equation. The slow decay of the infected population suggests that
the infection is locally maintained for long and it is difficult for the
disease to disappear completely.

##  I Introduction

Various nonlinear phenomena such as limit-cycle, chaos, and pattern formation
have been studied with numerical simulations. Population dynamics is an
important research field for nonlinear phenomena. The temporal oscillation in
ecosystems owing to the prey-predator interaction is reproduced by the Lotka-
Volterra equation LV . Chaotic dynamics was found in the logistic map and it
is applied to explain annual number fluctuations of bean weevils by May lg .
It is important to understand the spread of diseases from population dynamics.
Kermack and McKendrick proposed a mathematical model for epidemics in 1927 KM
. The model describes the population dynamics of $S$ (susceptible), $I$
(infected), and $R$ (recovered). Many authors studied the Kermack-McKendrick
model and generalized models, and applied it to the analysis of the spread of
diseases Murray ; Diekmann ; Capasso . Nobel performed numerical simulation of
the geographic and temporal development of plagues using the mathematical
model with diffusion terms Nobel . The model was also used to study
vaccination strategy Shulgin ; Laguzet . The recent COVID-19 pandemic have
induced intensive study of mathematical models of epidemic Science ; Costa .

In statistical physics, simplified stochastic models of infection on lattices
called contact process have been intensively studied. It was shown that the
critical behavior of the contact process is in the same universality class of
directed percolation Grassberger ; Kinzel ; Dickman ; Henkel ; Tome . In the
contact process, there is a phase transition from an infection phase to no
infection phase on finite-dimensional lattices. Pastor-Satorras and Vesignani
studied the spread of infection in scale-free networks Ves . They showed that
the infection does not disappear even for sufficiently small infection rate
when the exponent of the power law of link number is below a critical value.
Several authors studied the contact process with quenched disorder Noest ;
Cafiero . In the disordered contact process below the critical value, there is
a state called Griffiths phase, where the total infection decays to zero but
the decay occurs very slowly. Griffiths proposed a mechanism of the slow decay
originally in random spin systems Griffiths . There are compact infection
clusters of $N$ sites whose lifetime is O$(e^{aN})$ in the random system,
however, the probability of the clusters of $N$ sites is estimated as
O$(e^{-bN})$. As a result, the total infection decays as $t^{-b/a}$ Griffiths
; Noest . The Griffiths phases were studied in various complex networks. Munoz
; Moretti ; Cota ; Cota2  Even slower logarithmic decay was also reported in
some complex networks. Lee ; Cota3

In this paper, we study the decay process of the infection in the
deterministic Kermack-McKendrick model with spatially inhomogeneous infection
rate in one, two, and three dimensions, considering that infection often
spreads rapidly in densely populated zones or hot spots. The hot spots are
inhomogeneously distributed. In this paper, we will show slow decay of
infection such as a power-law decay with numerical simulations. In Sec. II, we
make a brief review of the Kermack-McKendrick model as an ordinary
differential equation. The SIR model has a unique property that there are
infinitely large number of stationary solutions in the SIR model depending on
the initial conditions. In Sec. III, we study the Kermack-McKendrick model
with one hot spot. We show the power-law decay of exponent $1/2$ in the one-
dimension model with one hot spot. In Sec. IV, we perform some numerical
simulations of the Kermack-McKendrick model with quenched randomness in one,
two, and three dimensions. In Sec.V, the results are summarized. Slow dynamics
appears in various complex systems such as glassy soft matter glass ; glass2 .
In the glassy soft matter, the dynamic heterogeneities of clusters and
trapping in a random energy landscape are closely related to the slows
dynamics Adams . The mechanism of slow dynamics in our model is different from
that of the glassy state, in that the freezing of particle motions does not
occur in our system. The slow decay in our model might be related with the
Griffiths phase, in that quenched randomness is important in both systems. It
is important that large clusters with large lifetimes appear with rare
probabilities in the Griffiths phase of the contact process. On the other
hand, the power-law decay of exponent $1/2$ in one dimension is due to the
diffusion process as shown in Sec. III in our model, and a coarsening process
is observed as shown in Sec. IV, which are not so important in the slow
dynamics in the Griffiths phase. It is somewhat similar to the power-law decay
found in the phase separation or ordering process near the first-order phase
transition, where both the diffusion and coarsening processes are important
Gunton ; Onuki . We think that the diffusion, coarsening, and quenched
randomness are important for the slow decay in our model, however, the
theoretical understanding is not sufficient yet and the details are left to
future study.

##  II SIR model

We make a very brief review of the dynamics of the Kermack-McKendrick model as
an ordinary differential equation KM ; Murray . The Kermack-McKendrick model
is expressed as

$\displaystyle\frac{dS}{dt}$ $\displaystyle=$ $\displaystyle-\beta SI,$ (1)
$\displaystyle\frac{dI}{dt}$ $\displaystyle=$ $\displaystyle\beta SI-\gamma
I,$ (2)

where $S$ and $I$ denote susceptible and infected populations. The parameters
$\beta$ and $\gamma$ denote infection and recovering rates. The recovered
population $R$ is calculated from

$\frac{dR}{dt}=\gamma I.$

Since the three variables $S$, $I$, and $R$ are used, the Kermack-McKendrick
model is called the SIR model. If the population $E$ of exposed persons before
the appearance of symptoms is included, the SIR model is generalized to the
SEIR model:

$\displaystyle\frac{dS}{dt}$ $\displaystyle=$ $\displaystyle-\beta SI,$
$\displaystyle\frac{dE}{dt}$ $\displaystyle=$ $\displaystyle\beta SI-\delta
E,$ $\displaystyle\frac{dI}{dt}$ $\displaystyle=$ $\displaystyle\delta
E-\gamma I$ (3)

where $\delta$ denotes the incidence rate. In the SIR model, recovered persons
are assumed not to be infected again owing to the immunity. However, there are
diseases such as Malaria, for which there is a possibility that recovered
persons are infected again. For such diseases, the SIS model where
$dS/dt=-\beta SI+\gamma I$ is used instead of Eq. (1). In the SIR and SEIR
models, $I$ decays to zero finally, but $I$ does not decay to zero in the SIS
model. The contact process is a stochastic version of SIS model. In this
paper, we consider mainly the SIR model, however, slow dynamics is observed
even in the SEIR model as shown in Sec. IV.

Figure 1: (a) Trajectories in $(S,I)$ space starting from $S(0)=4$, $2.5$, and
$1.3$ at $\beta=\gamma=1$. (b) Time evolutions of $I(t)$ for $S(0)=4$ and 2.5
at $\beta=\gamma=1$ in the semi-logarithmic scale. (c) $S_{\infty}$ as a
function of $S(0)$ at $\gamma=\beta=1$.

There is a stationary solution $S=S_{0}$ and $I=0$ to Eqs. (1) and (2). The
stationary state is unstable for $S_{0}\beta>\gamma$. Figure 1(a) shows
trajectories in $(S,I)$ space starting from the initial condition $S(0)=4$,
$2.5$, and $1.3$ at $\beta=\gamma=1$. The initial condition for $I$ is fixed
to be 0.00001. Figure 1(b) shows time evolutions of $I(t)$ for $S(0)=4$ and
2.5 at $\beta=\gamma=1$ in the semi-logarithmic scale. The infection $I(t)$
spreads initially and then decays to zero exponentially since the susceptible
population $S(t)$ becomes below $\gamma/\beta$. This is a state that the herd
immunity is attained. The final state is $(S,I)=(S_{\infty},0)$ where
$S_{\infty}<\gamma/\beta$ depends on the initial value $S(0)$. The uninfected
population $S_{\infty}$ takes a smaller value for a larger initial value
$S(0)$ of $S(t)$. $S_{\infty}$ can be calculated for the conserved quantity
$Q$ of this equation:

$Q=S+I-(\gamma/\beta)\log(S)$ (4)

If $I(0)$ is sufficiently small, $S_{\infty}$ is a solution of

$S_{\infty}-(\gamma/\beta)\log(S_{\infty})=S(0)-(\gamma/\beta)\log(S(0)).$

Figure 1(c) shows $S_{\infty}$ as a function of $S(0)$ at $\gamma=\beta=1$.
The ratio $S_{\infty}/S_{0}$ takes any value between 0 and 1, depending on the
initial value $S_{0}$.

##  III Slow decay of infection in SIR models with one hot spot

Hereafter, we consider spatially extended systems. The spread of infection
occurs in densely populated regions. In this paper, we consider mainly SIR
models with inhomogeneous infection rate on one-, two-, ang three-dimensional
lattices. In one dimension, the model equation is written:

$\displaystyle\frac{dS_{i}}{dt}$ $\displaystyle=$
$\displaystyle-\beta_{i}S_{i}I_{i}+D_{S}(S_{i+1}-2S_{i}+S_{i-1}),$ (5)
$\displaystyle\frac{dI_{i}}{dt}$ $\displaystyle=$
$\displaystyle\beta_{i}S_{i}I_{i}-\gamma I_{i}+D_{I}(I_{i+1}-2I_{i}+I_{i-1}),$
(6)

where $D_{S}$ and $D_{I}$ are diffusion constants. Periodic boundary
conditions are imposed and the system size is $N$. The initial conditions were
set to be $S_{i}=1$ and $I_{i}=0.00001$ for numerical simulations. Figure 2(a)
shows time evolutions of the total population of infection:
$SI=\sum_{i=1}^{N}I_{i}$ for $N=1000$ (solid line) and $N=200$ (dotted line)
at $D_{S}=1,D_{I}=1,\gamma=1$ in the double-logarithmic scale. The infection
rate is $\beta_{i}=\beta_{o}=0.9$ for $i\neq N/2$ and $\beta_{i}=\beta_{p}=3$
at $i=N/2$. The total infected population $SI$ decays exponentially for
$N=200$, however, it decays slowly with a power law $0.4/t^{1/2}$ for
$N=1000$. Even in the case of $N=200$, $SI$ decays in a power law of exponent
$1/2$ until $t\sim 2000$. There is a crossover time from the power law decay
to the exponential decay, and the crossover time increases with $N$. Figure
2(b) shows time evolutions of $SI=\sum_{i=1}^{N}I_{i}$ for $N=1000$ (solid
line) and $N=200$ (dotted line) in the double-logarithmic scale, when the
infection rate is $\beta_{i}=3$ for $N/2-7\leq i\leq N/2+7$ and
$\beta_{o}=0.9$ for the other region. Other parameters are the same as in Fig.
2(a). The total infected population $SI$ decays with the power law of exponent
$1/2$ for $N=1000$. The damping oscillation is observed near $t=100$ in
contrast to Fig. 2(a).

Figure 2: (a) Time evolutions of $SI=\sum_{i=1}^{N}I_{i}$ for $N=1000$ (solid
line) and $N=200$ (dotted line) at $D_{S}=1,D_{I}=1$, and $\gamma=1$ in the
double-logarithmic scale. The infection rate is $\beta_{i}=\beta_{o}=0.9$ for
$i\neq N/2$ and $\beta_{i}=3$ at $i=N/2$. The straight dashed line denotes a
power law of $1/t^{1/2}$. (b) Time evolutions of $SI=\sum_{i=1}^{N}I_{i}$ for
$N=1000$ (solid line) and $N=200$ (dotted line) in the double-logarithmic
scale. The infection rate is $\beta_{i}=3$ for $N/2-7\leq i\leq N/2+7$ and
$\beta_{0}=0.9$ for the other region. Other parameters are the same as in
Fig.2(a). Figure 3: (a) Time evolutions of $SI=\sum_{i=1}^{N}I_{i}$ at
$\beta_{o}=0.9$ (solid line), (b) 1.2 (dotted line) and 1.4 (dotted line) for
$N=1000$, $D_{S}=1,D_{I}=0.5,\gamma=1$ in the double-logarithmic scale. (b)
Critical value $\beta_{c}$ for the power-law decay as a function of $D_{I}$
for $N=1000$.

The power-law decay does not appear for large $\beta_{o}$, i.e.,the infection
rate of $i\neq N/2$ under a fixed value of $\beta_{i}=\beta_{p}$ at $i=N/2$.
Figure 3(a) shows the time evolutions of $SI=\sum_{i=1}^{N}I_{i}$ at
$\beta_{o}=0.9$ (solid line), (b) 1.2 (dashed line) and 1.4 (dotted line) for
$N=1000$. $D_{S}=1,D_{I}=0.5$, and $\gamma=1$ in the double-logarithmic scale.
The power law of exponent $1/2$ is observed at $\beta_{o}=0.9$ and 1.2 but
$SI$ decays rapidly at $\beta_{o}=1.4$. Figure 3(b) shows the critical value
of $\beta_{o}$ as a function of $D_{I}$ for $\gamma=1$ and $D_{S}=1$. In Fig.
2, we have shown numerical results for the case $\beta_{o}S_{i}(0)<\gamma$,
however, the critical value is larger than $\gamma/S_{i}(0)=1$, that is, the
power-law decay appears even when the infection spreads in the surrounding
region. The critical value of $\beta_{o}$ decreases with $D_{I}$. The slow
decay occurs in a wider parameter range for smaller $D_{I}$.

Although the total infection decays to zero finally under the periodic
boundary conditions, a stationary localized structure of infection is realized
if the fixed boundary conditions $S_{0}=S_{N+1}=S_{B}$ are imposed. Similar
localized structure is maintained for a long time even in systems of periodic
boundary conditions. Figures 4(a) and (b) show stationary profiles of (a)
$I_{i}$ and (b) $S_{i}$ for $N=200$ and $S_{B}=1$. The other parameters are
$D_{S}=1$, $D_{I}=1$, $\gamma=1$, and the infection rate is $\beta_{i}=0.9$
for $i\neq N/2$ and $\beta_{i}=3$ at $i=N/2$, which are the same as in Fig.
2(a). The infection occurs locally near $i=N/2$, $S_{i}$ diffuses into the
infection region, and is infected near $i=N/2$. Outside of the infection
region, the profile of $S_{i}$ has a nearly constant slope, that is,
$S_{i}=S_{B}-\alpha i$ for $i<N/2$ and $S_{i}=S_{B}-\alpha(N-i)$ for $i>N/2$.
The profile of $I_{i}$ is approximated as $I_{i}=I_{0}e^{-\lambda|i-N/2|}$.
Figure 4(c) shows the relationship between $\beta_{N/2}$ and $S_{N/2}$. The
localized structure disappears for $\beta_{N/2}<1.5$ and
$S_{i}=S_{N/2}=S_{B}=1$ is satisfied.

Figure 4: Stationary profiles of (a) $I_{i}$ and (b) $S_{i}$ for $N=200$ and
$S_{B}=1$. The other parameters are $D_{S}=1,D_{I}=1,\gamma=1$, and the
infection rate is $\beta_{i}=0.9$ for $i\neq N/2$ and $\beta_{i}=3$ for
$i=N/2$. (c) Relationship between $\beta_{N/2}$ and $S_{N/2}$.

Next, we make an analysis of the localized state and the power-law decay of
exponent 1/2. If the continuum approximation is taken, Eqs. (5) and (6) are
rewritten as

$\displaystyle\frac{\partial S}{\partial t}$ $\displaystyle=$
$\displaystyle-\beta SI+D_{S}\nabla^{2}S,$ (7) $\displaystyle\frac{\partial
I}{\partial t}$ $\displaystyle=$ $\displaystyle\beta SI-\gamma
I+D_{I}\nabla^{2}I,$ (8)

where $\nabla^{2}=\partial^{2}/\partial x^{2}$ in one dimension, however, the
model will be later extended to two and three dimensions. For the stationary
state in one dimension, $S(x)$ and $I(x)$ satisfy

$\displaystyle D_{S}\frac{\partial^{2}S}{\partial x^{2}}$ $\displaystyle=$
$\displaystyle\beta(x)S(x)I(x),$ (9) $\displaystyle
D_{I}\frac{\partial^{2}I}{\partial x^{2}}$ $\displaystyle=$
$\displaystyle\\{\gamma-\beta(x)S(x)\\}I(x).$ (10)

If $I(x)$ and $\beta(x)$ are assumed to be $I(x)=I(N/2)e^{-\lambda|x-N/2|}$
and $\beta(x)=\beta_{0}+\beta_{1}\delta(x-N/2)$, $\lambda$ satisfies

$\lambda=\sqrt{(\gamma-\beta_{0}S(N/2))/D_{I}}.$ (11)

Equation (10) yields

$D_{I}\left\\{\left(\frac{\partial I}{\partial
x}\right)_{(N/2)_{+}}-\left(\frac{\partial I}{\partial
x}\right)_{(N/2)_{-}}\right\\}=-\beta_{1}S(N/2)I(N/2)$ (12)

at $x=N/2$, which leads to

$2D_{I}\lambda=\beta_{1}S(N/2).$ (13)

Equations (11) and (13) yield

$S(N/2)=\frac{\sqrt{4D_{I}^{2}\beta_{0}^{2}+4D_{I}\gamma\beta_{1}^{2}}-2D_{I}\beta_{0}}{\beta_{1}^{2}}.$
(14)

The dashed line in Fig. 4(c) shows the relationship between
$\beta_{N/2}=\beta_{0}+\beta_{1}$ and $S(N/2)$ by Eq. (14). Fairly good
agreement with direct numerical results is seen. If the infection region is
sufficiently small, Eq. (9) gives

$D_{S}\left\\{\left(\frac{\partial S}{\partial
x}\right)_{(N/2)_{+}}-\left(\frac{\partial S}{\partial
x}\right)_{(N/2)_{-}}\right\\}=\beta_{1}S(N/2)I(N/2).$ (15)

From Eq. (15), $I(N/2)$ is approximated at

$I(N/2)=\frac{2D_{S}(S_{B}-S(N/2))/(N/2)}{\beta_{1}S(N/2)}$ (16)

If the boundary conditions are not fixed to $S_{B}$, the power-law decay is
observed. Similar power-law decay is observed even for $D_{I}=0$. Figure 5(a)
shows the time evolution of $I_{N/2}$ for $N=1000$ at $D_{S}=1,\;D_{I}=0$, and
$\gamma=1$ in the double-logarithmic scale. The infection rate is
$\beta_{i}=0.9$ for $i\neq N/2$ and $\beta_{i}=3$ at $i=N/2$. Figure 5(b)
shows two snapshots of $S_{i}$ at $t=100$ and 5000. The population $S_{i}$ at
$i=N/2$ is $1/3$ and the width of the depression increases as $t^{1/2}$.

Figure 5: (a) Time evolution of $I(N/2)$ for $N=1000$ at
$D_{S}=1,\;D_{I}=0,\;\gamma=1$ in the double-logarithmic scale. The infection
rate is $\beta_{i}=0.9$ for $i\neq N/2$ and $\beta_{i}=3$ for $i=N/2$. (b) Two
snapshots of $S_{i}$ at $t=100$ and 5000.

In the continuum approximation, the solution satisfies

$\frac{\partial S}{\partial t}=D_{S}\frac{\partial^{2}S}{\partial x^{2}}$ (17)

for $x\neq N/2$, since $I(x,t)$ is almost zero for $x\neq N/2$, and
$2D_{S}\partial S/\partial x=\beta_{N/2}S(N/2)I(N/2)$ at $x=(N/2)_{+}$.
Because $S(N/2)$ is fixed to be $\gamma/\beta_{N/2}$ in case of $D_{I}=0$,
$I(N/2)$ is determined to be $2D_{S}\partial S/\partial x/\gamma$. That is,
the population of infection is determined by the diffusion process of
susceptible population for large $t$. The slope $\partial S/\partial x$ at
$x=(N/2)_{+}$ can be calculated from the solution of $S(x,t)$ to the diffusion
equation. The diffusion equation can be solved by the Fourier transform

$S(x,t)=S(N/2)+\frac{1-S(N/2)}{N/2}(x-N/2)+\sum_{n=1}^{\infty}A_{n}(t)\sin\\{2\pi
n(x-N/2)/N\\},$ (18)

for $x>N/2$. The Fourier amplitude $A_{n}(t)$ satisfies
$A_{n}(t)=A_{n}(0)e^{-D_{S}4\pi^{2}n^{2}t/N^{2}}$, where $A_{n}(0)$ is
determined from the initial condition $S(x,0)=(S(N/2)-1)\delta(x-N/2)+1$ as

$A_{n}(0)=(1-S(N/2))\frac{2}{\pi n}.$

$I(N/2)$ is evaluated as

$I(N/2)=2D_{S}\frac{\beta_{N/2}}{\gamma}\frac{\partial S}{\partial
x}=2D_{S}\frac{\beta_{N/2}}{\gamma}\sum_{n=1}^{\infty}(1-S(N/2))\frac{2}{N}e^{-4D_{S}\pi^{2}n^{2}t/N^{2}}.$

Since the summation is negligible for $n$ satisfying
$4D_{S}\pi^{2}n^{2}t/N^{2}>>1$ or $n>>N/(4D_{S}\pi^{2}t)^{1/2}$, $I(N/2)$ is
approximated as

$I(N/2)=\frac{4(1-\gamma/\beta_{N/2})\sqrt{D_{s}}}{\pi\gamma\sqrt{t}},$

if the contribution $(1-S(N/2))/(N/2)$ in the second term of Eq. (18) is
neglected for large $N/2$. This is a reason of the power law of exponent
$1/2$. The dashed line in Fig. 5(a) denotes this relation, which is good
approximation to the direct numerical simulation.

Two- and three dimensional models are expressed with Eqs. (7) and (8) if
$\nabla^{2}$ is rewritten as $\partial^{2}/\partial
x^{2}+\partial^{2}/\partial y^{2}$ in two dimensions and
$\partial^{2}/\partial x^{2}+\partial^{2}/\partial y^{2}+\partial^{2}/\partial
z^{2}$ in three dimensions. If $S$, $I$ and $\beta$ depends only the radius
$r$, Eqs. (7) and (8) are expressed as

$\displaystyle\frac{\partial S}{\partial t}$ $\displaystyle=$
$\displaystyle-\beta SI+D_{S}\left\\{\frac{\partial^{2}S}{\partial
r^{2}}+\frac{d-1}{r}\frac{\partial S}{\partial r}\right\\},$ (19)
$\displaystyle\frac{\partial I}{\partial t}$ $\displaystyle=$
$\displaystyle\beta SI-\gamma I+D_{I}\left\\{\frac{\partial^{2}I}{\partial
r^{2}}+\frac{d-1}{r}\frac{\partial I}{\partial r}\right\\},$ (20)

where $d$ denotes the dimension 2 or 3. We have performed numerical simulation
of Eqs. (19) and (20) as a one-dimensional discrete system similar to Eqs. (5)
and (6) with additional terms of $D_{S}(d-1)/r\partial S/\partial r$ and
$D_{I}(d-1)/r\partial I/\partial r$. The system size is $N=1000$, and $r=0$ is
set in the middle of $i=N/2$ and $N/2+1$. That is, $r=i-(N+1)/2$. Parameters
are $D_{S}$, $\gamma=1$, and $\beta_{i}=0.9$ for $i\neq N/2,N/2+1$ and
$\beta_{i}=3$ for $i=N/2$ and $N/2+1$. Figure 6(a) shows the time evolution of
$I(N/2)$ at $D_{I}=1$ (solid line) and $D_{I}=0$ (dashed line) for $d=2$.
Figure 6(b) shows $1/I(N/2)$ as a function of $\log(t)$ at $D_{I}=1$ (solid
line) and $D_{I}=0$ (dashed line). The solid line can be approximated at
$5.8\log(t)+1.2$ and the dashed at $0.71\log(t)+0.5$. That is, $I(N/2)$ decays
as $1/\log(t)$, The decay of infection in two dimensions is slower than in one
dimension. Figure 6(c) shows a snapshot of $S(x,t)$ at $t=5000$ for $D_{I}=0$.
The value of $S(r)$ is almost equal to be $\gamma/\beta(N/2)=1/3$, however,
the slope of $S(r)$ seems to be divergent at $r=0$ in contrast to the one-
dimensional case. The analysis using the expansion by the Bessel and Neumann
functions instead of the Fourier series expansion might explain the
logarithmic decay. Figure 6(d) shows the time evolution of $I(N/2)$ at
$D_{I}=1$ (solid line) and $D_{I}=0$ (dashed line) for $d=3$. The other
parameter values are the same as in Fig. 6(a). The decay of $I(N/2)$ is even
slower in three dimensions. The functional form of the decay is not known yet.
The infection hardly disappears, probably because the susceptible population
in the surrounding region diffuses into the center of the three-dimensional
hot spot and the site number of the surrounding region is large in the three-
dimensional system compared to the one- and two-dimensional systems.

Figure 6: (a) Time evolution of $I(N/2)$ at $D_{I}=1$ (solid line) and
$D_{I}=0$ (dashed line) for $d=2$. Parameters are $D_{S}=1$, $\gamma=1$, and
$\beta_{i}=0.9$ for $i\neq N/2,N/2+1$ and $\beta_{i}=3$ for $i=N/2$ and
$N/2+1$. (b) $1/I(N/2)$ as a function of $\log(t)$ at $D_{I}=1$ (solid line)
and $D_{I}=0$ (dashed line). (c) Snapshot profile of $S(x,t)$ at $t=5000$ for
$D_{I}=0$. (d) Time evolution of $I(N/2)$ at $D_{I}=1$ (solid line) and
$D_{I}=0$ (dashed line) for $d=3$.

##  IV Slow decay of infection in random SIR models

In this section, we study random SIR models in one, two and three dimensions.
The one-dimensional model is again expressed as

$\displaystyle\frac{dS_{i}}{dt}$ $\displaystyle=$
$\displaystyle-\beta_{i}S_{i}I_{i}+D_{S}(S_{i+1}-2S_{i}+S_{i-1}),$
$\displaystyle\frac{dI_{i}}{dt}$ $\displaystyle=$
$\displaystyle\beta_{i}S_{i}I_{i}-\gamma I_{i}+D_{I}(I_{i+1}-2I_{i}+I_{i-1}).$
(21)

where $\beta_{i}$ is a uniform random number between 0 and $\beta_{m}$. Figure
7(a) shows time evolutions of $SI=\sum_{i}I_{i}$ for Eq. (21) at $D_{I}=1$ and
$\beta_{m}=1.5$ (solid line) and $D_{I}=0$ and $\beta_{m}=1.1$ (dashed line).
The other parameters are $D_{S}=1$ and $\gamma=1$. The system size is
$N=5000$. The initial condition is $S_{i}=1$ and $I_{i}=0.00001$. Figure 7(a)
shows that $SI$ decays roughly with a power law for large $t$. The exponent is
around 1 for $D_{I}=1$ and $\beta_{m}=1.5$, and 1.35 for $D_{I}=0$ and
$\beta_{m}=1.1$. The exponent depends on the parameters, Figure 7(b) shows
three snapshot profiles of $I_{i}$ at $t=200,300$ and 2000 for $D_{I}=1$ and
$\beta_{m}=1.5$. The other parameters are the same as the ones in Fig. 7(a).
Figure 7(c) shows four snapshot profiles of $S_{i}$ at $t=200,300,2000$, and
10000 in the same frame. The infection occurs at many hot spots at $t=200$,
and the spatial profile is intermittent as shown in Fig. 7(b). The number of
hot spots decrease with time. The cusp points in Fig, 7(c) correspond the hot
spots. $S_{i}$ decreases with time by the diffusion to the hot spots and
infection at the hot spots. Figure 7(d) shows the local average of $\beta_{i}$
(red dotted line), that is,
$\bar{\beta}_{i}=(1/2)\sum_{j=i-10}^{i+10}\beta_{i}$, $100I_{i}$ (blue dashed
line) at $t=1000$ and the maximum value of $I_{i}(t)$ (green solid line) for
$0<t<10000$ in the range of $2000<i<3000$. In most cases, hot spots appear
near the points where the locally averaged infection rate $\bar{\beta}_{i}$ is
large. Infection is stamped out at some hot spots, however, strong hot spots
survive for long. The lifetime of a hot spot with an interval $L$ is estimated
as O$(L^{2})$, because the width of the diffusion field increases as $t^{1/2}$
as shown in Fig. 5(b) and the power-law decay changes to an exponential decay
when the width reaches the size of interval. If the intervals between strong
hot spots become longer after the burnout of some hot spots, the lifetime of
the survived hot spot becomes even longer. The power-law decay of $t^{-1/2}$
by the diffusion effect and the increase of lifetime by the coarsening process
might be the origin of the slow decay in one dimension.

Figure 7: (a) Time evolutions of $SI=\sum_{i}I_{i}$ for Eq. (21) at $D_{I}=1$
and $\beta_{m}=1.5$ (solid line) and $D_{I}=0$ and $\beta_{m}=1.1$ (dashed
line) for $D_{S}=1$,and $\gamma=1$. (b) Three snapshot profiles of $I_{i}$ at
$t=200,300$ and 2000 for $D_{I}=1$ and $\beta_{m}1.5$. (c) Four snapshot
profiles of $S_{i}$ at $t=200,300,2000$, and 10000 from above. (d) Local
average of $\beta_{i}$ (red dotted line), $100I_{i}$ (blue dashed line) at
$t=1000$ and the maximum value of $I_{i}(t)$ for $0<t<10000$ (green solid
line).

We study the SEIR model to check the generality of the slow dynamics: The
model equation is

$\displaystyle\frac{dS_{i}}{dt}$ $\displaystyle=$
$\displaystyle-\beta_{i}S_{i}I_{i}+D_{S}(S_{i+1}-2S_{i}+S_{i-1}),$
$\displaystyle\frac{dE_{i}}{dt}$ $\displaystyle=$
$\displaystyle\beta_{i}S_{i}I_{i}-\delta_{i}E_{i}+D_{E}(E_{i+1}-2E_{i}+E_{i-1}),$
$\displaystyle\frac{dI_{i}}{dt}$ $\displaystyle=$
$\displaystyle\delta_{i}E_{i}-\gamma I_{i}+D_{I}(I_{i+1}-2I_{i}+I_{i-1}),$
(22)

where $S_{i}$, $E_{i}$ and $I_{i}$ denote respectively the susceptible,
exposed, and infected populations, and $\delta_{i}$ denotes the incidence
rate. The parameters are set to be $D_{S}=D_{E}=1$, $D_{I}=0.5$, $\gamma=1$,
and $\delta=2$. Figure 8(a) shows the time evolutions of $SI=\sum I_{i}$ in a
system with one hot spot for $N=1000$, where $\beta_{i}=\beta_{o}=0.9$ (solid
line), 1.2 (dashed line), and 1.4 (dotted line) for $i\neq N/2$ and
$\beta_{i}=3$ at $i=N/2$. The power law decay of exponent $1/2$ is observed
for $\beta_{o}=0.9$ and 1.2. An exponential decay is observed at
$\beta_{0}=1.4$. Figure 8(b) shows the time evolution of $SI$ in a random SEIR
model of $N=5000$ where $\beta_{i}$ takes a uniform random number between 0
and 1.4. $SI$ exhibits a power-law decay of exponent around 0.85. These
results are similar to those in the SIR model.

Figure 8: (a) Time evolutions of $SI=\sum I_{i}$ in a system with one hot spot
for $\beta_{o}=0.9$ (solid line), 1.2 (dashed line), and 1.4 (dotted line) for
$N=1000$, $D_{S}=D_{E}=1$, $D_{I}=0.5$, $\gamma=1$, and $\delta=2$. (b) Time
evolution of $SI$ in a random SEIR model of $N=5000$ where $\beta_{i}$ takes a
uniform random number between 0 and 1.4. The dashed line denotes a power law
of exponent 0.85.

The two-dimensional random SIR model is expressed as

$\displaystyle\frac{dS_{i,j}}{dt}$ $\displaystyle=$
$\displaystyle-\beta_{i,j}S_{i,j}I_{i,j}+D_{S}(S_{i+1,j}+S_{i-1,j}+S_{i,j}+S_{i,j-1}-4S_{i,j}),$
$\displaystyle\frac{dI_{i,j}}{dt}$ $\displaystyle=$
$\displaystyle\beta_{i,j}S_{i,j}I_{i,j}-\gamma
I_{i,j}+D_{I}(I_{i+1,j}+I_{i-1,j}+I_{i,j+1}+I_{i,j-1}-4I_{i,j}),$ (23)

Figure 9(a) shows time evolutions of $SI=\sum_{i,j}I_{i,j}$ for Eq. (23) at
$D_{S}=D_{I}=1,\gamma=2$ in random systems, when $\beta_{i,j}$ takes a uniform
random value between 0 and $\beta_{m}=3.5$ (solid line) or between 0 and
$\beta_{m}=3.2$ (dashed line). The dotted line is $SI\propto 1/t^{0.85}$. The
system size is $600\times 600$. The initial condition is $S_{i,j}=1$ and
$I_{i,j}=0.00001$. The total population of infection seems to decay in a power
law in these numerical simulations. Figure 9(b) shows some snapshots of
$I_{i,j}$ at a section of $j=N/2$ when $\beta_{i,j}$ takes a uniform random
value between 0 and 3.2 at $t=50,100,\cdots$, 500. Localized clusters of
infection survive for long. The number of hot spots decreases with time, or a
coarsening occurs also in two dimensions.

Figure 9: (a) Time evolutions of $SI=\sum_{i,j}I_{i,j}$ for Eq. (23) at
$D_{S}=D_{I}=1,\gamma=2$ and $\beta_{m}=3.5$(solid line) and $\beta_{m}=3.2$
(dashed line). The dotted line is $SI\propto 1/t^{0.85}$. The system size is
$600\times 600$. (b) Snapshots of $I_{i,j}$ at a section of $j=N/2$ when
$\beta_{i,j}$ for $\beta_{m}=3.2$ at $t=50,100,\cdots$, 500.

The three-dimensional model is expressed as

$\displaystyle\frac{dS_{i,j,k}}{dt}$ $\displaystyle=$
$\displaystyle-\beta_{i,j,k}S_{i,j,k}I_{i,j,k}+D_{S}(S_{i+1,j,k}+S_{i-1,j,k}+S_{i,j+1,k}+S_{i,j-1,k}+S_{i,j,k+1}+S_{i,j,k-1}-6S_{i,j,k}),$
$\displaystyle\frac{dI_{i,j,k}}{dt}$ $\displaystyle=$
$\displaystyle\beta_{i,j,k}S_{i,j,k}I_{i,j,k}-\gamma
I_{i,j,k}+D_{I}(I_{i+1,j,k}+I_{i-1,j,k}+I_{i,j+1,k}+I_{i,j-1,k}+I_{i,j,k+1}+I_{i,j,k-1}-6I_{i,j,k}),$
(24)

Figure 10(a) shows time evolution of $I_{N/2+1,N/2+1,N/2}$ for Eq. (24) at
$D_{S}=D_{I}=0.2$, and $\gamma=1$ when $\beta_{i,j,k}=3$ for $i=N/2,N/2+1$,
$j=N/2,N/2+1$, $k=N/2,N/2+1$ and $\beta_{i,j,k}=0.9$ for the other sites. The
system size is $100\times 100\times 100$. The initial condition is
$S_{i,j,k}=1$ and $I_{i,j,k}=0.00001$. This is a three-dimensional simulation
similar to the case shown in Fig. 6(d). $I_{N/2+1,N/2+1,N/2}$ is almost
constant after $t>50$, that is, the localized infection is maintained for very
long. Figure 10(b) shows time evolutions of $SI=\sum_{i,j,k}I_{i,j,k}$ for Eq.
(24) at $D_{S}=D_{I}=0.2$, and $\gamma=5$ in random systems, where
$\beta_{i,j,k}$ takes a uniform random value between 0 and 10 (solid line) or
between 0 and 6 (dashed line). The dotted line is $SI\propto 1/t^{1.1}$. The
total population of infection seems to decay in a power law also in these
random three-dimensional systems.

Figure 10:  (a) Time evolution of $I_{N/2+1,N/2+1,N/2}$ for Eq. (24) at
$D_{S}=D_{I}=0.2$, and $\gamma=1$ when $\beta_{i,j,k}=3$ for $i=N/2,N/2+1$,
$j=N/2,N/2+1$, $k=N/2,N/2+1$ and $\beta_{i,j,k}=0.9$ for the other sites. The
system size is $100\times 100\times 100$. (b) Time evolutions of
$SI=\sum_{i,j,k}I_{i,j,k}$ for Eq. (24) at $D_{S}=D_{I}=0.2$, and $\gamma=5$
in random systems, where $\beta_{i,j,k}$ takes a uniform random value between
0 and 10 (solid line) or between 0 and 6 (dashed line). The dotted line is
$SI\propto 1/t^{1.1}$.

##  V Summary

We have found slow decay of infection in the Kermack-McKendrick model with
spatially inhomogeneous infection rate in some parameter range. First, we have
studied the Kermack-McKendrick model with one hot spot where the infection
rate is locally higher than the surrounding region. We have shown
theoretically a power-law decay of exponent $1/t^{1/2}$ in the one-dimensional
system with a spatially localized hot shot. The slow decay occurs as $1/\log
t$ in the two-dimensional system with one hot spot, and the decay seems to be
even slower in three dimensions. Next, we have studied the random Kermack-
McKendrick model, and found the power-law type slow decay in one, two, and
three dimensions. We found that the infection occurs locally at hot spots. The
uninfected persons in the surrounding area around the hot spots diffuse into
the hot spots and are infected at the hot spots. The number of hot spots
decreases in time and the lifetime of the survived hot spots become even
longer because the surrounding areas of the survived hot spots become larger.
The mechanism of the slow decay in our system is unique in that the diffusion,
coarsening, and quenched randomness are important, although there is some
similarity to the slows dynamics in the Griffiths phase in the contact process
Noest ; Cafiero  and the phase transition dynamics Gunton ; Onuki . However,
the slower decay in our two- and three- dimensional SIR models with one hot
spot, and the exponent of the power-law decay in one-, two-, and three-
dimensional random SIR models are not well understood. They are left to future
study. Our finding of the slow decay of the infected population suggests that
the infection is locally maintained for long time and the diseases are hardly
stamped out in some cases.

## References

  * (1) A. J. Lotka, J. Phys. Chem. 14, 271 (1910).
  * (2) R. M. May, Nature 261, 459 (1976).
  * (3) W. O. Kermack and A. G. McKendrick, Proc. Roy. Soc. A 115, 700 (1927).
  * (4) J. D. Murray Mathematical Biology (Springer-Verlag, Berlin Heidelberg, 2003).
  * (5) O. Diekmann and J. A. P. Heesterbeek, Mathematical Epidemiology of Infectious Diseases: Model Building, Analysis, and Interpretation (John Wiley and Sons, Chichester, 2000).
  * (6) V. Cappaso and G. Serio, Math. Bio. 42, 43 (1978).
  * (7) J. V. Nobel, Nature 250, 726 (1974).
  * (8) B. Shulgin, L. Stone, Z, Agur, Bul.Math. Biol.60, 1123 (1998).
  * (9) L. Lagzet and G. Turinici, Bul. Math. Biol. 77, 1955 (2015).
  * (10) e.g., L. Ferretti et al. Science 10.1126/science.abb6936 (2020). M. Chinazzim J.T. Davis, M. Ajelli et al. Science 10.1126/sciemce.aba9757.
  * (11) G. S. Costa, W. Cota, and S. C. Ferreira, Phys. Rev. Research 2, 043306 (2020).
  * (12) P. Grassberger and A. de la Torre, Ann. Phys. 122, 373 (1979).
  * (13) W. Kinzel, Z. Phys. B 58, 229 (1985).
  * (14) R. Dockman, Int. J. Mod. Phys.C 4, 271 (1993).
  * (15) M. Henkel, H. Hinrichsen, S. Lübeck, and M. Pleimling, Non-equilibrium phase transitions (Springer Dordrecht, 2008).
  * (16) T. Tomé and R. M. Ziff, Phys. Rev. E 82, 051921 (2010).
  * (17) R. Pastor-Satorras and A. Vespignani, Phys. Rev. Lett. 86, 3200 (2001).
  * (18) A. J. Noest, Phys. Rev. Lett. 57, 90 (1986).
  * (19) R. Cafiero, A. Gabrielli, and M. A. Muñoz, Phys. Rev. E 57, 5060 (1998).
  * (20) R. B. Griffiths, Phys. Rev. Lett. 23, 17 (1969).
  * (21) M. A. Muñoz, R. Juhász, C. Castellano, and G. Ódor, Phys. Rev. Lett. 105, 128701 (2010).
  * (22) P. Moretti and M. Mu$\tilde{\rm{n}}$oz, Nature Comm. 4, 1 (2011).
  * (23) W. Cota, S. C. Ferreira, and G. Ódor, Phys. Rev. E 93, 032322 (2016).
  * (24) W. Cota, G. Ódor, and S. C. Ferreira, Sci. Rep. 8, 9144 (2018).
  * (25) H. K. Lee, P.S. Shim. and J. D. Noh, Phys. Rev. E 87, 062812 (2013).
  * (26) W. Cota, S. C. Ferreira, and G. Ódor, Phys. Rev. E 93, 032322 (2-16).
  * (27) C. A. Angell, J. Phys. Chem. Solids, 49, 863 (1988).
  * (28) W. Gotze and L. Sjogren, Rep. Prog. Phys. 55, 241 (1992).
  * (29) G. Adams and J. H. Gibbs, J. Chem. Phys. 43, 139 (1965).
  * (30) J. D. Gunton, M. van Miguel, and P. S. Sahni, Phase Transitions and Critical Phenomena, Vol. 8 edited by C. Domb and J. Lebowitz (Academic Press 1983).
  * (31) A. Onuki, Phase Transition Dynamics (Cambrdge University Press, 2009).


# Semi-supervised Word Sense Disambiguation with Neural Models

Dayu Yuan Julian Richardson Ryan Doherty Colin Evans Eric Altendorf
Google, Mountain View CA, USA
<EMAIL_ADDRESS>
###### Abstract

Determining the intended sense of words in text – word sense disambiguation
(WSD) – is a long-standing problem in natural language processing. Recently,
researchers have shown promising results using word vectors extracted from a
neural network language model as features in WSD algorithms. However, a simple
average or concatenation of word vectors for each word in a text loses the
sequential and syntactic information of the text. In this paper, we study WSD
with a sequence learning neural net, LSTM, to better capture the sequential
and syntactic patterns of the text. To alleviate the lack of training data in
all-words WSD, we employ the same LSTM in a semi-supervised label propagation
classifier. We demonstrate state-of-the-art results, especially on verbs.

##  1 Introduction

00footnotetext:  This work is licensed under a Creative Commons Attribution
4.0 International License. License details:
http://creativecommons.org/licenses/by/4.0/

Word sense disambiguation (WSD) is a long-standing problem in natural language
processing (NLP) with broad applications. Supervised, unsupervised, and
knowledge-based approaches have been studied for WSD [Navigli, 2009]. However,
for all-words WSD, where all words in a corpus need to be annotated with word
senses, it has proven extremely challenging to beat the strong baseline, which
always assigns the most frequent sense of a word without considering the
context [Pradhan et al., 2007a, Navigli, 2009, Navigli et al., 2013, Moro and
Navigli, 2015]. Given the good performance of published supervised WSD systems
when provided with significant training data on specific words [Zhong and Ng,
2010], it appears lack of sufficient labeled training data for large
vocabularies is the central problem.

One way to leverage unlabeled data is to train a neural network language model
(NNLM) on the data. Word embeddings extracted from such a NNLM (often Word2Vec
[Mikolov et al., 2013]) can be incorporated as features into a WSD algorithm.
Iacobacci et al. [Iacobacci et al., 2016] show that this can substantially
improve WSD performance and indeed that competitive performance can be
attained using word embeddings alone.

In this paper, we describe two novel WSD algorithms. The first is based on a
Long Short Term Memory (LSTM) [Hochreiter and Schmidhuber, 1997]. Since this
model is able to take into account word order when classifying, it performs
significantly better than an algorithm based on a continuous bag of words
model (Word2vec) [Mikolov et al., 2013, Iacobacci et al., 2016], especially on
verbs.

We then present a semi-supervised algorithm which uses label propagation
[Talukdar and Crammer, 2009, Ravi and Diao, 2016] to label unlabeled sentences
based on their similarity to labeled ones. This allows us to better estimate
the distribution of word senses, obtaining more accurate decision boundaries
and higher classification accuracy.

The best performance was achieved by using an LSTM language model with label
propagation. Our algorithm achieves state-of-art performance on many SemEval
all-words tasks. It also outperforms the most-frequent-sense and Word2Vec
baselines by $10$% (see Section 5.2 for details).

Organization:  We review related work in Section 2. We introduce our
supervised WSD algorithm in Section 3, and the semi-supervised WSD algorithm
in Section 4. Experimental results are discussed in Section 5. We provide
further discussion and future work in Section 6.

##  2 Related Work

The development of large lexical resources, such as WordNet [Fellbaum, 1998]
and BabelNet [Navigli and Ponzetto, 2012], has enabled knowledge-based
algorithms which show promising results on all-words prediction tasks
[Ponzetto and Navigli, 2010, Navigli et al., 2013, Moro and Navigli, 2015].
WSD algorithms based on supervised learning are generally believed to perform
better than knowledge-based WSD algorithms, but they need large training sets
to perform well [Pradhan et al., 2007a, Navigli et al., 2007, Navigli, 2009,
Zhong and Ng, 2010]. Acquiring large training sets is costly. In this paper,
we show that a supervised WSD algorithm can perform well with $\sim 20$
training examples per sense.

In the past few years, much progress has been made on using neural networks to
learn word embeddings [Mikolov et al., 2013, Levy and Goldberg, 2014], to
construct language models [Mikolov et al., 2011], perform sentiment analysis
[Socher et al., 2013], machine translation [Sutskever et al., 2014] and many
other NLP applications.

A number of different ways have been studied for using word embeddings in WSD.
There are some common elements:

  * •

Context embeddings. Given a window of text $w_{n-k},...,w_{n},...,w_{n+k}$
surrounding a focus word $w_{n}$ (whose label is either known in the case of
example sentences or to be determined in the case of classification), an
embedding for the context is computed as a concatenation or weighted sum of
the embeddings of the words $w_{i},i\neq n$. Context embeddings of various
kinds are used in both [Chen et al., 2014] and [Iacobacci et al., 2016].

  * •

Sense embeddings. Embeddings are computed for each word sense in the word
sense inventory (e.g. WordNet). In [Rothe and Schütze, 2015], equations are
derived relating embeddings for word senses with embeddings for
undisambiguated words. The equations are solved to compute the sense
embeddings. In [Chen et al., 2014], sense embeddings are computed first as
weighted sums of the embeddings of words in the WordNet gloss for each sense.
These are used in an initial bootstrapping WSD phase, and then refined by a
neural network which is trained on this bootstrap data.

  * •

Embeddings as SVM features. Context embeddings [Iacobacci et al., 2016,
Taghipour and Ng, 2015b], or features computed by combining context embeddings
with sense embeddings [Rothe and Schütze, 2015], can be used as additional
features in a supervised WSD system e.g. the SVM-based IMS [Zhong and Ng,
2010]. Indeed Iacobacci et al. [Iacobacci et al., 2016] found that using
embeddings as the only features in IMS gave competitive WSD performance.

  * •

Nearest neighbor classifier. Another way to perform classification is to find
the word sense whose sense embedding is closest, as measured by cosine
similarity, to the embedding of the classification context. This is used, for
example, in the bootstrapping phase of [Chen et al., 2014].

  * •

Retraining embeddings. A feedforward neural network can be used to jointly
perform WSD and adjust embeddings [Chen et al., 2014, Taghipour and Ng,
2015b].

In our work, we start with a baseline classifier which uses $1000$-dimensional
embeddings trained on a $100$ billion word news corpus using Word2Vec [Mikolov
et al., 2013]. The vocabulary consists of the most frequent $1,000,000$ words,
without lemmatization or case normalization. Sense embeddings are computed by
averaging the context embeddings of sentences which have been labeled with
that sense. To classify a word in a context, we assign the word sense whose
embedding has maximal cosine similarity with the embedding of the context.
This classifier has similar performance to the best classifier in [Iacobacci
et al., 2016] when SemCor is used as a source of labeled sentences. The
Word2Vec embeddings are trained using a bag of words model, i.e. without
considering word order in the training context, and word order is also not
considered in the classification context. In Section 3 we show that using a
more expressive language model which takes account of word order yields
significant improvements.

Semi-supervised learning has previously been applied successfully to word
sense disambiguation. In [Yarowsky, 1995] bootstrapping was used to learn a
high precision WSD classifier. A low recall classifier was learned from a
small set of labeled examples, and the labeled set then extended with those
sentences from an unlabeled corpus which the classifier could label with high
confidence. The classifier was then retrained, and this iterative training
process continued to convergence. Additional heuristics helped to maintain the
stability of the bootstrapping process. The method was evaluated on a small
data set.

In [Niu et al., 2005], a label propagation algorithm was proposed for word
sense disambiguation and compared to bootstrapping and a SVM supervised
classifier. Label propagation can achieve better performance because it
assigns labels to optimize a global objective, whereas bootstrapping
propagates labels based on local similarity of examples.

In Section 4 we describe our use of label propagation to improve on nearest
neighbor for classification.

##  3 Supervised WSD with LSTM

Neural networks with long short-term memory (LSTM) units [Hochreiter and
Schmidhuber, 1997] make good language models which take into account word
order [Sundermeyer et al., 2012]. We train a LSTM language model to predict
the held-out word in a sentence. As shown in Figure 1, we first replace the
held-out word with a special symbol $, and then, after consuming the remaining
words in the sentence, project the $h$ dimensional hidden layer to a $p$
dimensional context layer, and finally predict the held out word with softmax.
By default, the LSTM model has $2048$ hidden units, $512$ dimensional context
layer and $512$ dimensional word embeddings. We also studied other settings,
see Section 5.2.2 for details. We train the LSTM on a news corpus of about
$100$ billion tokens, with a vocabulary of $1,000,000$ words. Words in the
vocabulary are neither lemmatized nor case normalized.

Our LSTM model is different from that of K geb ck and Salomonsson [Kågebäck
and Salomonsson, 2016]. We train a LSTM language model, which predicts a held-
out word given the surrounding context, with a large amount of unlabeled text
as training data. The huge training dataset allows us to train a high-capacity
model ($2048$ hidden units, $512$ dimensional embeddings), which achieves high
precision without overfitting. In our experiments, this directional LSTM model
was faster and easier to train than a bidirectional LSTM, especially given our
huge training dataset. K geb ck and Salomonsson’s LSTM directly predicts the
word senses and it is trained with a limited number of word sense-labeled
examples. Although regularization and dropout are used to avoid overfitting
the training data, the bidirectional LSTM is small with only $74+74$ neurons
and $100$ dimensional word embeddings [Kågebäck and Salomonsson, 2016].
Because our LSTM is generally applicable to any word, it achieves high
performance on all-words WSD tasks (see Section 5 for details), which is the
focus of this paper. K geb ck and Salomonsson’s LSTM is only evaluated on
lexical sample WSD tasks of SemEval 2 and 3 [Kågebäck and Salomonsson, 2016].

Figure 1: LSTM: Replace the focus word $w_{3}$ with a special symbol $ and
predict $w_{3}$ at the end of the sentence.

The behavior of the LSTM can be intuited by its predictions. Table 1 shows the
top $10$ words predicted by an LSTM language model for the word ‘stock’ in
sentences containing various senses of ‘stock’.

id |  sentence |  top 10 predictions from LSTM |  sense
---|---|---|---
1 |  Employee compensation is offered in the form of cash and/or stock. |  cash, stock, equity, shares, loans, bonus, benefits, awards, equivalents, deposits |  sense#1
2 |  The stock would be redeemed in five years, subject to terms of the company’s debt. |  bonds, debt, notes, shares, stock, balance, securities, rest, Notes, debentures |
3 |  These stores sell excess stock or factory overruns . |  inventory, goods, parts, sales, inventories, capacity, products, oil, items, fuel |  sense#2
4 |  Our soups are cooked with vegan stock and seasonal vegetables. |  foods, food, vegetables, meats, recipes, cheese, meat, chicken, pasta, milk |  sense#3
query |  In addition, they will receive stock in the reorganized company, which will be named Ranger Industries Inc. |  shares, positions, equity, jobs, awards, representation, stock, investments, roles, funds |  ?
Table 1: Top predictions of ‘stock’ in 5 sentences of different word senses

In our initial experiments, we computed similarity between two contexts by the
overlap between their bags of predicted words. For example (Table 1) the top
predictions for the query overlap most with the LSTM predictions for ‘sense#1’
—we predict that ‘sense#1’ is the correct sense. This bag of predictions,
while easily interpretable, is just a discrete approximation to the internal
state of the LSTM when predicting the held out word. We therefore directly use
the LSTM’s context layer from which the bag of predictions was computed as a
representation of the context (see Figure 1). Given context vectors extracted
from the LSTM, our supervised WSD algorithms classify a word in a context by
finding the sense vector which has maximum cosine similarity to the context
vector (Figure 2(a)). We find the sense vectors by averaging context vectors
of all training sentences of the same sense. We observed in a few cases that
the context vector is far from the held-out word’s embedding, especially when
the input sentence is not informative. For example, the LSTM language model
will predict “night” for the input sentence “I fell asleep at [work].” when we
hold out “work”. Currently, we treat the above cases as outliers. We would
like explore alternative solutions, e.g., forcing the model to predict words
that are close to one sense vector of the held-out word, in further work. As
can be seen in SemEval all-words tasks and Tables 6, this LSTM model has
significantly better performance than the Word2Vec models.

##  4 Semi-supervised WSD

(a) nearest neighbor

(b) semi-supervised WSD with label propagation

Figure 2: WSD classifiers. Filled nodes represent labeled sentences. Unfilled
nodes represent unlabeled sentences.

The non-parametric nearest neighbor algorithm described in Section 3 has the
following drawbacks:

  * •

It assumes a spherical shape for each sense cluster, being unable to
accurately model the decision boundaries given the limited number of examples.

  * •

It has no training data for, and does not model, the sense prior, omitting an
extremely powerful potential signal.

To overcome these drawbacks we present a semi-supervised method which augments
the labeled example sentences with a large number of unlabeled sentences from
the web. Sense labels are then propagated from the labeled to the unlabeled
sentences. Adding a large number of unlabeled sentences allows the decision
boundary between different senses to be better approximated.

A label-propagation graph consists of (a) vertices with a number of labeled
seed nodes and (b) undirected weighted edges. Label propagation (LP) [Talukdar
and Crammer, 2009] iteratively computes a distribution of labels on the
graph’s vertices to minimize a weighted combination of:

  * •

The discrepancy between seed labels and their computed labels distributions.

  * •

The disagreement between the label distributions of connected vertices.

  * •

A regularization term which penalizes distributions which differ from the
prior (by default, a uniform distribution).

We construct a graph for each lemma with labeled vertices for the labeled
example sentences, and unlabeled vertices for sentences containing the lemma,
drawn from some additional corpus. Vertices for sufficiently similar sentences
(based on criteria discussed below) are connected by an edge whose weight is
the cosine similarity between the respective context vectors, using the LSTM
language model. To classify an occurrence of the lemma, we create an
additional vertex for the new sentence and run LP to propagate the sense
labels from the seed vertices to the unlabeled vertices.

Figure 2(b) (b) illustrates the graph configuration. Spatial proximity
represents similarity of the sentences attached to each vertex and the shape
of each node represents the word sense. Filled nodes represent seed nodes with
known word senses. Unfilled nodes represent sentences with no word sense
label, and the ? represents the word we want to classify.

With too many edges, sense labels propagate too far, giving low precision.
With too few, sense labels do not propagate sufficiently, giving low recall.
We found that the graph has about the right density for common senses when we
ranked vertex pairs by similarity and connected the pairs above the $95$
percentile. This may still leave rare senses sparsely connected, so we
additionally added edges to ensure that every vertex is connected to at least
$10$ other vertices. Our experiments (Table 9) show that this setting achieves
good performance on WSD, and the performance is stable when the percentile
ranges between $85$ to $98$. Since it requires running LP for every
classification, the algorithm is slow compared to the nearest neighbor
algorithm.

##  5 Experiments

We evaluated the LSTM algorithm with and without label propagation on standard
SemEval all-words tasks using WordNet as the inventory. Our proposed
algorithms achieve state-of-art performance on many SemEval all-words WSD
tasks. In order to assess the effects of training corpus size and language
model capacity we also evaluate our algorithms using the New Oxford American
Dictionary (NOAD) inventory with SemCor [Miller et al., 1993] or MASC
111http://www.anc.org/MASC/About.html/.

###  5.1 SemEval Tasks

In this section, we study the performance of our classifiers on Senseval2
[Edmonds and Cotton, 2001], Senseval3 [Snyder and Palmer, 2004], SemEval-2007
[Pradhan et al., 2007b], SemEval-2013 Task 12 [Navigli et al., 2013] and
SemEval-2015 task 13 [Moro and Navigli, 2015] 222We mapped all senses to
WordNet3.0 by using maps in
https://wordnet.princeton.edu/wordnet/download/current-version/ and
http://web.eecs.umich.edu/ mihalcea/downloads.html. We focus the study on all-
words WSD tasks. For a fair comparison with related works, the classifiers are
evaluated on all words (both polysemous and monosemous).

Following related works, we use SemCor or OMSTI [Taghipour and Ng, 2015a] for
training. In our LP classifiers, unlabeled data for each lemma consists either
of $1000$ sentences which contain the lemma, randomly sampled from the web, or
all OMSTI sentences (without labels) which contain the lemma.

| Senseval2 | Senseval3 | SemEval7 | SemEval7-Coarse | SemEval13
---|---|---|---|---|---
model | all |  n. |  all |  n. |  all |  n. |  all |  n. |  n.
IMS + Word2Vec (T:SemCor) | 0.634 |  0.742 |  0.653 |  0.701 |  0.578 |  0.686 |  |  |
IMS + Word2Vec (T:OMSTI) | 0.683 |  0.777 |  0.682 |  0.741 |  0.591 |  0.715 |  |  |
Taghipour and Ng [Taghipour and Ng, 2015b] |  |  |  0.682 |  |  |  |  |  |
Chen et al. [Chen et al., 2014] |  |  |  |  |  |  |  0.826 |  0.853 |
Weissenborn et al. [Weissenborn et al., 2015] |  |  |  |  0.688 |  |  0.660 |  |  0.855 |  0.728
Word2Vec (T:SemCor) | 0.678 |  0.737 |  0.621 |  0.714 |  0.585 |  0.673 |  0.795 |  0.814 |  0.661
LSTM (T:SemCor) | 0.736 |  0.786 |  0.692 |  0.723 |  0.642 |  0.723 |  0.828 |  0.834 |  0.670
LSTM (T:OMSTI) | 0.724 |  0.777 |  0.643 |  0.680 |  0.607 |  0.673 |  0.811 |  0.820 |  0.673
LSTMLP (T:SemCor, U:OMSTI) | 0.739 |  0.797 |  0.711 |  0.748 |  0.637 |  0.704 |  0.843 |  0.834 |  0.679
LSTMLP (T:SemCor, U:1K) | 0.738 |  0.796 |  0.718 |  0.763 |  0.635 |  0.717 |  0.836 |  0.831 |  0.695
LSTMLP (T:OMSTI, U:1K) | 0.744 |  0.799 |  0.710 |  0.753 |  0.633 |  0.717 |  0.833 |  0.825 |  0.681
Table 2: F1 scores on SemEval all-words tasks. T:SemCor stands for models
trained with SemCor. U:OMSTI stands for using OMSTI as unlabeled sentences in
semi-supervised WSD. IMS + Word2Vec scores are from [Iacobacci et al., 2016]

Table 2 shows the Sem-Eval results. Our proposed algorithms achieve the
highest all-words F1 scores except for Sem-Eval 2013. Weissenborn et
al.[Weissenborn et al., 2015] only disambiguates nouns, and it outperforms our
algorithms on Sem-Eval 2013 by $4$%, but is ranked behind our algorithms on
Senseval-3 and SemEval-7 tasks with an F1 score more than $6$% lower than our
algorithms. Unified WSD [Chen et al., 2014] has the highest F1 score on Nouns
(Sem-Eval-7 Coarse), but our algorithms outperform Unified WSD on other part-
of-speech tags.

##### Settings

For a fair comparison of Word2Vec and LSTM, we do not use pre-trained word-
embeddings as in [Iacobacci et al., 2016], but instead train the Word2Vec and
LSTM models on a $100$ billion word news corpus 333The training corpus could
not be released, but we have plans to open source the well-trained models with
a vocabulary of the most frequent 1,000,000 words. Our self-trained word
embeddings have similar performance to the pre-trained embeddings, as shown in
Table 2. The Word2Vec word vectors are of dimension ${1024}$. The LSTM model
has $2048$ hidden units, and inputs are ${512}$-dimensional word vectors. We
train the LSTM model by minimizing sampled softmax loss [Jean et al., 2014]
with Adagrad [Duchi et al., 2011]. The learning rate is $0.1$. We experimented
with other learning rates, and observed no significant performance difference
after the training converges. We also downsample frequent terms in the same
way as [Mikolov et al., 2013].

##### Word2Vec vectors Vs. LSTM

To better compare LSTM with word vectors we also build a nearest neighbor
classifier using Word2Vec word embeddings and SemCor example sentences,
Word2Vec (T:SemCor). It performs similar to IMS + Word2Vec (T:SemCor), a SVM-
based classifier studied in [Iacobacci et al., 2016]. Table 2 shows that the
LSTM classifier outperforms the Word2Vec classifier across the board.

##### SemCor Vs. OMSTI

Contrary to the results observed in [Iacobacci et al., 2016], the LSTM
classifier trained with OMSTI performs worse than that trained with SemCor. It
seems that the larger size of the OMSTI training data set is more than offset
by noise present in its automatically generated labels. While the SVM
classifier studied in [Iacobacci et al., 2016] may be able to learn a model
which copes with this noise, our naive nearest neighbor classifiers do not
have a learned model and deal less well with noisy labels.

##### Label propagation

We use the implementation of DIST_EXPANDER [Ravi and Diao, 2016]. We test the
label propagation algorithm with SemCor or OMSTI as labeled data sets and
OMSTI or $1000$ random sentences from the web per lemma as unlabeled data. The
algorithm performs similarly on the different data sets.

Table 3 shows the results of Sem-Eval 2015. The LSTM LP classifier with an
LSTM language model achieves the highest scores on nouns and adverbs as well
as overall F1. The LSTM classifier has the highest F1 on verbs.

algorithm | all | n. | v. | adj. | adv.
---|---|---|---|---|---
LIMSI | 0.647 |  |  |  | 0.795
DFKI |  | 0.703 | 0.577 |  |
UNIBA |  |  |  | 0.790 |
BFS Baseline | 0.663 | 0.667 | 0.551 | 0.821 | 0.825
Word2Vec (T:SemCor) | 0.667 | 0.661 | 0.555 | 0.789 | 0.810
LSTM (T:SemCor) | 0.721 | 0.713 | 0.642 | 0.813 | 0.845
LSTMLP (T:SemCor, U:1K) | 0.726 | 0.728 | 0.622 | 0.813 | 0.857
Table 3: F1 Scores of SemEval-2015 English Dataset. The BFS baseline uses
BabelNet first sense.

###  5.2 NOAD Eval

Many dictionary lemmas and senses have no examples in SemCor or OSTMI, giving
rise to losses in all-words WSD when these corpora are used as training data.
The above SemEval scores do not distinguish errors caused by missing training
data for certain labels or inaccurate classifier. To better study the proposed
algorithms, we train the classifiers with the New Oxford American Dictionary
(NOAD) [Stevenson and Lindberg, 2010], in which there are example sentences
for each word sense.

####  5.2.1 Word Sense Inventory

The NOAD focuses on American English and is based on the Oxford Dictionary of
English (ODE) [Stevenson, 2010]. It distinguishes between coarse (core) and
fine-grained (sub) word senses in the same manner as ODE. Previous
investigations [Navigli, 2006, Navigli et al., 2007] using the ODE have shown
that coarse-grained word senses induced by the ODE inventory address problems
with WordNet’s fine-grained inventory, and that the inventory is useful for
word sense disambiguation.

For our experiments, we use NOAD’s core senses, and we also use lexicographer-
curated example sentences from the Semantic English Language Database
(SELD)444http://oxfordgls.com/our-content/english-language-content/, provided
by Oxford University Press. We manually annotated all words of the English
language SemCor corpus and MASC corpora with NOAD word senses in order to
evaluate performance 555https://research.google.com/research-
outreach.html#/research-outreach/research-datasets. Table 4 shows the total
number of polysemes (more than one core sense) and average number of senses
per polyseme in NOAD/SELD (hereafter, NOAD), SemCor and MASC. Both SemCor and
MASC individually cover around $45\%$ of NOAD polysemes and $62\%$ of senses
of those polysemes.

|  | noun | verb | adj. | adv.
---|---|---|---|---|---
Number of polysemous lemmas in dictionary / corpus | NOAD | 8097 | 2124 | 2126 | 266
SemCor | 2833 | 1362 | 911 | 193
MASC | 2738 | 1250 | 829 | 181
Sense count per polyseme | NOAD | 2.46 | 2.58 | 2.30 | 2.47
SemCor | 1.44 | 1.69 | 1.49 | 1.84
MASC | 1.48 | 1.66 | 1.48 | 2.01
Table 4: NOAD polysemous lemma in NOAD, SemCor and MASC

Table 5 gives the number of labeled sentences of these datasets. Note that
although NOAD has more labeled sentences than SemCor or MASC, the average
numbers of sentences per sense of these datasets are similar. This is because
NOAD has labeled sentences for each word sense, whereas SemCor (MASC) only
covers a subset of lemmas and senses (Table 4). The last column of Table 5
shows that each annotated word in SemCor and MASC has an average of more than
$4$ NOAD corse senses. Hence, a random guess will have a precision around
$1/4$.

| example count (in 1000’s) | example count per sense | number of candidate senses per example
---|---|---|---
dataset | all | n. | v. | adj. | adv
NOAD | 580 | 312 | 150 | 95 | 13 | 18.43 | 3.1
SemCor | 115 | 38 | 57 | 11.6 | 8.6 | 14.27 | 4.1
MASC | 133 | 50 | 57 | 12.7 | 13.6 | 17.38 | 4.2
Table 5: Number of examples in each dataset and the average sense count per
example.

In the default setting, we use NOAD example sentences as labeled training data
and evaluate on SemCor and MASC. We evaluate all polysemous words in the
evaluation corpus.

####  5.2.2 LSTM classifier

We compare our algorithms with two baseline algorithms:

  * •

Most frequent sense: Compute the sense frequency (from a labeled corpus) and
label word $w$ with $w$’s most frequent sense.

  * •

Word2Vec: a nearest-neighbor classifier with Word2Vec word embedding, which
has similar performance to cutting-edge algorithms studied in [Iacobacci et
al., 2016] on SemEval tasks.

eval data | SemCor | MASC |
---|---|---|---
model | train data |  all |  n. |  v. |  adj. |  adv. |  all |  n. |  v. |  adj. |  adv.
Frequent Sense | SemCor |  |  |  |  |  |  0.753 |  0.799 |  0.713 |  0.758 |  0.741
Frequent Sense | MASC |  0.752 |  0.751 |  0.749 |  0.737 |  0.789 |  |  |  |  |
Word2Vec | NOAD |  0.709 |  0.783 |  0.657 |  0.736 |  0.693 |  0.671 |  0.790 |  0.562 |  0.724 |  0.638
Word2Vec | SemCor |  |  |  |  |  |  0.692 |  0.806 |  0.592 |  0.754 |  0.635
Word2Vec | NOAD,SemCor |  |  |  |  |  |  0.678 |  0.808 |  0.565 |  0.753 |  0.604
Word2Vec | MASC |  0.698 |  0.785 |  0.619 |  0.766 |  0.744 |  |  |  |  |
Word2Vec | NOAD,MASC |  0.695 |  0.801 |  0.605 |  0.767 |  0.719 |  |  |  |  |
LSTM | NOAD |  0.786 |  0.796 |  0.782 |  0.781 |  0.784 |  0.786 |  0.805 |  0.772 |  0.776 |  0.786
LSTM | SemCor |  |  |  |  |  |  0.799 |  0.843 |  0.767 |  0.808 |  0.767
LSTM | NOAD,SemCor |  |  |  |  |  |  0.812 |  0.846 |  0.786 |  0.816 |  0.798
LSTM | MASC |  0.810 |  0.825 |  0.799 |  0.809 |  0.825 |  |  |  |  |
LSTM | NOAD,MASC |  0.821 |  0.834 |  0.814 |  0.818 |  0.821 |  |  |  |  |
Table 6: F1 scores of LSTM algorithm in comparison with baselines eval data | SemCor | MASC |
---|---|---|---
model | train data |  all |  n. |  v. |  adj. |  adv. |  all |  n. |  v. |  adj. |  adv.
LSTM | NOAD |  0.769 |  0.791 |  0.759 |  0.751 |  0.672 |  0.780 |  0.791 |  0.768 |  0.780 |  0.726
LSTM | SemCor |  |  |  |  |  |  0.656 |  0.663 |  0.668 |  0.643 |  0.581
LSTM | NOAD,SemCor |  |  |  |  |  |  0.796 |  0.805 |  0.790 |  0.794 |  0.742
LSTM | MASC |  0.631 |  0.653 |  0.606 |  0.617 |  0.600 |  |  |  |  |
LSTM | NOAD,MASC |  0.782 |  0.803 |  0.774 |  0.761 |  0.688 |  |  |  |  |
Table 7: Macro F1 scores of LSTM classifier

Table 6 compares the F1 scores of the LSTM and baseline algorithms. LSTM
outperforms Word2Vec by more than $10\%$ over all words, where most of the
gains are from verbs and adverbs. The results suggest that syntactic
information, which is well modeled by LSTM but ignored by Word2Vec, is
important to distinguishing word senses of verbs and adverbs.

##### Change of training data

By default, the WSD classifier uses the NOAD example sentences as training
data. We build a larger training dataset by adding labeled sentences from
SemCor and MASC, and study the change of F1 scores in Table 6. Across all part
of speech tags and datasets, F1 scores increase after adding more training
data. We further test our algorithm by using SemCor (or MASC) as training data
(without NOAD examples). The SemCor (or MASC) trained classifier is on a par
with the NOAD trained classifier on F1 score. However, the macro F1 score of
the former is much lower than the latter, as shown in Table 7, because of the
limited coverage of rare senses and words in SemCor and MASC.

##### Change of language model capacity

In this experiment, we change the LSTM model capacity by varying the number of
hidden units $h$ and the dimensions of the input embeddings $p$ and measuring
F1. Figure 3 shows strong positive correlation between F1 and the capacity of
the language model. However, larger models are slower to train and use more
memory. To balance the accuracy and resource usage, we use the second best
LSTM model ($h=2048$ and $p=512$) by default.

(a) eval on SemCor

(b) eval on MASC

Figure 3: F1 scores of LSTM models with different capacity: h is the number of
hidden units; p is the embedding dimension.

####  5.2.3 Semi-supervised WSD

We evaluate our semi-supervised WSD classifier in this subsection. We
construct the graph as described in Section 4 and run LP to propagate sense
labels from the seed vertices to the unlabeled vertices. We evaluate the
performance of the algorithm by comparing the predicted labels and the gold
labels on eval nodes. As can be observed from Table 8, LP did not yield clear
benefits when using the Word2Vec language model. We did see significant
improvements, 6.3% increase on SemCor and 7.3% increase on MASC, using LP with
the LSTM language model. We hypothesize that this is because LP is sensitive
to the quality of the graph distance metric.

eval data | SemCor | MASC |
---|---|---|---
model | train data |  all |  n. |  v. |  adj. |  adv. |  all |  n. |  v. |  adj. |  adv.
Word2Vec LP | NOAD |  0.642 |  0.733 |  0.554 |  0.705 |  0.725 |  0.643 |  0.752 |  0.524 |  0.726 |  0.664
LSTM LP | NOAD |  0.822 |  0.859 |  0.800 |  0.817 |  0.816 |  0.831 |  0.865 |  0.806 |  0.825 |  0.821
LSTM LP | NOAD,SemCor |  |  |  |  |  |  0.872 |  0.897 |  0.852 |  0.865 |  0.868
LSTM LP | NOAD,MASC |  0.873 |  0.883 |  0.870 |  0.858 |  0.874 |  |  |  |  |
Table 8: F1 scores of label propagation

##### Change of seed data:

As can be seen in Table 8, LP substantially improves classifier F1 when the
training datasets are SemCor+NOAD or MASC+NOAD. As discussed in Section 4, the
improvement may come from explicitly modeling the sense prior. We did not see
much performance lift by increasing the number of unlabeled sentences per
lemma.

##### Change of graph density:

By default, we construct the LP graph by connecting two nodes if their
affinity is above $95$% percentile. We also force each node to connect to at
least $10$ neighbors to prevent isolated nodes. Table 9 shows the performance
of the LP algorithm by changing the percentile threshold. The F1 scores are
relatively stable when the percentile ranges between $85$ to $98$, but
decrease when the percentile drops to $80$. Also, it takes longer to run the
LP algorithm on a denser graph. We pick the $95$ percentile in our default
setting to achieve both high F1 scores and short running time.

| SemCor | MASC
---|---|---
pos-tag |  98 |  95 |  90 |  85 |  80 |  70 |  98 |  95 |  90 |  85 |  80 |  70
overall |  0.823 |  0.822 |  0.823 |  0.818 |  0.813 |  0.800 |  0.827 |  0.831 |  0.835 |  0.830 |  0.824 |  0.806
n. |  0.848 |  0.859 |  0.852 |  0.846 |  0.840 |  0.828 |  0.863 |  0.865 |  0.868 |  0.865 |  0.861 |  0.847
v. |  0.810 |  0.800 |  0.803 |  0.797 |  0.792 |  0.778 |  0.800 |  0.806 |  0.806 |  0.799 |  0.794 |  0.769
Table 9: F1 scores of the LSTM LP trained on NOAD with varying graph density.

##  6 Conclusions and Future Work

In this paper, we have presented two WSD algorithms which combine (1) LSTM
neural network language models trained on a large unlabeled text corpus, with
(2) labeled data in the form of example sentences, and, optionally, (3)
unlabeled data in the form of additional sentences. Using an LSTM language
model gave better performance than one based on Word2Vec embeddings. The best
performance was achieved by our semi-supervised WSD algorithm which builds a
graph containing labeled example sentences augmented with a large number of
unlabeled sentences from the web, and classifies by propagating sense labels
through this graph.

Several unanswered questions suggest lines of future work. Since our general
approach is amenable to incorporating any language model, further developments
in NNLMs may permit increased performance. We would also like to better
understand the limitations of language modeling for this task: we expect there
are situations – e.g., in idiomatic phrases – where per-word predictions carry
little information.

We believe our model should generalize to languages other than English, but
have not yet explored this. Character-level LSTMs [Kim et al., 2015] may
provide robustness to morphology and diacritics and may prove useful even in
English for spelling errors and out of vocabulary words.

We would like to see whether our results can be improved by incorporating
global (document) context and multiple embeddings for polysemous words [Huang
et al., 2012].

Finally, many applications of WSD systems for nominal resolution require
integration with resolution systems for named entities, since surface forms
often overlap [Moro et al., 2014, Navigli and Ponzetto, 2012]. This will
require inventory alignment work and model reformulation, since we currently
use no document-level, topical, or knowledge-base coherence features.

We thank our colleagues and the anonymous reviewers for their insightful
comments on this paper.

## References

  * [Chen et al., 2014]  Xinxiong Chen, Zhiyuan Liu, and Maosong Sun.  2014\.  A unified model for word sense representation and disambiguation.  In EMNLP, pages 1025–1035. Citeseer.
  * [Duchi et al., 2011]  John Duchi, Elad Hazan, and Yoram Singer.  2011\.  Adaptive subgradient methods for online learning and stochastic optimization.  J. Mach. Learn. Res., 12:2121–2159, July.
  * [Edmonds and Cotton, 2001]  Philip Edmonds and Scott Cotton.  2001\.  SENSEVAL-2: Overview.  In Proceedings of SENSEVAL-2 Second International Workshop on Evaluating Word Sense Disambiguation Systems, pages 1–5, Toulouse, France, July. Association for Computational Linguistics.
  * [Fellbaum, 1998]  Christiane Fellbaum.  1998\.  WordNet.  Wiley Online Library.
  * [Hochreiter and Schmidhuber, 1997]  Sepp Hochreiter and Jürgen Schmidhuber.  1997\.  Long short-term memory.  Neural computation, 9(8):1735–1780.
  * [Huang et al., 2012]  Eric H Huang, Richard Socher, Christopher D Manning, and Andrew Y Ng.  2012\.  Improving word representations via global context and multiple word prototypes.  In Proceedings of the 50th Annual Meeting of the Association for Computational Linguistics: Long Papers-Volume 1, pages 873–882. Association for Computational Linguistics.
  * [Iacobacci et al., 2016]  Ignacio Iacobacci, Mohammad Taher Pilehvar, and Roberto Navigli.  2016\.  Embeddings for word sense disambiguation: An evaluation study.  In ACL. Citeseer.
  * [Jean et al., 2014]  S. Jean, K. Cho, R. Memisevic, and Y. Bengio.  2014\.  On Using Very Large Target Vocabulary for Neural Machine Translation.  ArXiv e-prints, December.
  * [Kågebäck and Salomonsson, 2016]  M. Kågebäck and H. Salomonsson.  2016\.  Word Sense Disambiguation using a Bidirectional LSTM.  ArXiv e-prints.
  * [Kim et al., 2015]  Yoon Kim, Yacine Jernite, David Sontag, and Alexander M Rush.  2015\.  Character-aware neural language models.  arXiv preprint arXiv:1508.06615.
  * [Levy and Goldberg, 2014]  Omer Levy and Yoav Goldberg.  2014\.  Neural word embedding as implicit matrix factorization.  In Advances in Neural Information Processing Systems, pages 2177–2185.
  * [Mikolov et al., 2011]  Tomáš Mikolov, Stefan Kombrink, Lukáš Burget, Jan Honza Černockỳ, and Sanjeev Khudanpur.  2011\.  Extensions of recurrent neural network language model.  In Acoustics, Speech and Signal Processing (ICASSP), 2011 IEEE International Conference on, pages 5528–5531. IEEE.
  * [Mikolov et al., 2013]  Tomas Mikolov, Ilya Sutskever, Kai Chen, Greg S Corrado, and Jeff Dean.  2013\.  Distributed representations of words and phrases and their compositionality.  In C. J. C. Burges, L. Bottou, M. Welling, Z. Ghahramani, and K. Q. Weinberger, editors, Advances in Neural Information Processing Systems 26, pages 3111–3119. Curran Associates, Inc.
  * [Miller et al., 1993]  George A Miller, Claudia Leacock, Randee Tengi, and Ross T Bunker.  1993\.  A semantic concordance.  In Proceedings of the workshop on Human Language Technology, pages 303–308. Association for Computational Linguistics.
  * [Moro and Navigli, 2015]  Andrea Moro and Roberto Navigli.  2015\.  Semeval-2015 task 13: multilingual all-words sense disambiguation and entity linking.  Proc. of SemEval, pages 288–297.
  * [Moro et al., 2014]  Andrea Moro, Alessandro Raganato, and Roberto Navigli.  2014\.  Entity linking meets word sense disambiguation: a unified approach.  Transactions of the Association for Computational Linguistics, 2:231–244.
  * [Navigli and Ponzetto, 2012]  Roberto Navigli and Simone Paolo Ponzetto.  2012\.  BabelNet: The automatic construction, evaluation and application of a wide-coverage multilingual semantic network.  Artificial Intelligence, 193:217–250.
  * [Navigli et al., 2007]  Roberto Navigli, Kenneth C Litkowski, and Orin Hargraves.  2007\.  Semeval-2007 task 07: Coarse-grained English all-words task.  In Proceedings of the 4th International Workshop on Semantic Evaluations, pages 30–35. Association for Computational Linguistics.
  * [Navigli et al., 2013]  Roberto Navigli, David Jurgens, and Daniele Vannella.  2013\.  Semeval-2013 task 12: Multilingual word sense disambiguation.  In Second Joint Conference on Lexical and Computational Semantics (* SEM), volume 2, pages 222–231.
  * [Navigli, 2006]  Roberto Navigli.  2006\.  Meaningful clustering of senses helps boost word sense disambiguation performance.  In Proceedings of the 21st International Conference on Computational Linguistics and the 44th Annual Meeting of the Association for Computational Linguistics, ACL-44, pages 105–112, Stroudsburg, PA, USA. Association for Computational Linguistics.
  * [Navigli, 2009]  Roberto Navigli.  2009\.  Word sense disambiguation: A survey.  ACM Computing Surveys (CSUR), 41(2):10.
  * [Niu et al., 2005]  Zheng-Yu Niu, Dong-Hong Ji, and Chew Lim Tan.  2005\.  Word sense disambiguation using label propagation based semi-supervised learning.  In Proceedings of the 43rd Annual Meeting on Association for Computational Linguistics, pages 395–402. Association for Computational Linguistics.
  * [Ponzetto and Navigli, 2010]  Simone Paolo Ponzetto and Roberto Navigli.  2010\.  Knowledge-rich word sense disambiguation rivaling supervised systems.  In Proceedings of the 48th Annual Meeting of the Association for Computational Linguistics, ACL ’10, pages 1522–1531, Stroudsburg, PA, USA. Association for Computational Linguistics.
  * [Pradhan et al., 2007a]  Sameer S Pradhan, Edward Loper, Dmitriy Dligach, and Martha Palmer.  2007a.  Semeval-2007 task 17: English lexical sample, SRL and all words.  In Proceedings of the 4th International Workshop on Semantic Evaluations, pages 87–92. Association for Computational Linguistics.
  * [Pradhan et al., 2007b]  Sameer S Pradhan, Edward Loper, Dmitriy Dligach, and Martha Palmer.  2007b.  Semeval-2007 task 17: English lexical sample, srl and all words.  In Proceedings of the 4th International Workshop on Semantic Evaluations, pages 87–92. Association for Computational Linguistics.
  * [Ravi and Diao, 2016]  Sujith Ravi and Qiming Diao.  2016\.  Large scale distributed semi-supervised learning using streaming approximation.  In AISTATS.
  * [Rothe and Schütze, 2015]  Sascha Rothe and Hinrich Schütze.  2015\.  Autoextend: Extending word embeddings to embeddings for synsets and lexemes.  In Proceedings of the 53rd Annual Meeting of the Association for Computational Linguistics and the 7th International Joint Conference on Natural Language Processing (Volume 1: Long Papers), pages 1793–1803, Beijing, China, July. Association for Computational Linguistics.
  * [Snyder and Palmer, 2004]  Benjamin Snyder and Martha Palmer.  2004\.  The english all-words task.  In Senseval-3: Third International Workshop on the Evaluation of Systems for the Semantic Analysis of Text, pages 41–43.
  * [Socher et al., 2013]  Richard Socher, Alex Perelygin, Jean Y Wu, Jason Chuang, Christopher D Manning, Andrew Y Ng, and Christopher Potts.  2013\.  Recursive deep models for semantic compositionality over a sentiment treebank.  In Proceedings of the conference on empirical methods in natural language processing (EMNLP), volume 1631, page 1642. Citeseer.
  * [Stevenson and Lindberg, 2010]  Angus Stevenson and Christine A. Lindberg.  2010\.  New Oxford American Dictionary.
  * [Stevenson, 2010]  Angus Stevenson.  2010\.  Oxford Dictionary of English.
  * [Sundermeyer et al., 2012]  Martin Sundermeyer, Ralf Schlüter, and Hermann Ney.  2012\.  LSTM neural networks for language modeling.  In INTERSPEECH, pages 194–197.
  * [Sutskever et al., 2014]  Ilya Sutskever, Oriol Vinyals, and Quoc V Le.  2014\.  Sequence to sequence learning with neural networks.  In Z. Ghahramani, M. Welling, C. Cortes, N. D. Lawrence, and K. Q. Weinberger, editors, Advances in Neural Information Processing Systems 27, pages 3104–3112. Curran Associates, Inc.
  * [Taghipour and Ng, 2015a]  Kaveh Taghipour and Hwee Tou Ng.  2015a.  One million sense-tagged instances for word sense disambiguation and induction.  CoNLL 2015, page 338.
  * [Taghipour and Ng, 2015b]  Kaveh Taghipour and Hwee Tou Ng.  2015b.  Semi-supervised word sense disambiguation using word embeddings in general and specific domains.  In Proceedings of the 2015 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, pages 314–323, Denver, Colorado, May–June. Association for Computational Linguistics.
  * [Talukdar and Crammer, 2009]  Partha Pratim Talukdar and Koby Crammer.  2009\.  New regularized algorithms for transductive learning.  In Proceedings of the European Conference on Machine Learning and Knowledge Discovery in Databases: Part II, ECML PKDD ’09, pages 442–457, Berlin, Heidelberg. Springer-Verlag.
  * [Weissenborn et al., 2015]  Dirk Weissenborn, Leonhard Hennig, Feiyu Xu, and Hans Uszkoreit.  2015\.  Multi-objective optimization for the joint disambiguation of nouns and named entities.  In Proceedings of the 53rd Annual Meeting of the Association for Computational Linguistics and the 7th International Joint Conference on Natural Language Processing of the Asian Federation of Natural Language Processing, pages 596–605.
  * [Yarowsky, 1995]  David Yarowsky.  1995\.  Unsupervised word sense disambiguation rivaling supervised methods.  In Proceedings of the 33rd annual meeting on Association for Computational Linguistics, pages 189–196. Association for Computational Linguistics.
  * [Zhong and Ng, 2010]  Zhi Zhong and Hwee Tou Ng.  2010\.  It makes sense: A wide-coverage word sense disambiguation system for free text.  In ACL, ACLDemos ’10.