# Efficient distributed representations beyond negative sampling

Lorenzo Dall’Amico<EMAIL_ADDRESS>
ISI Foundation
Via Chisola 5, 10126 Torino, Italy Enrico Maria Belliardo
<EMAIL_ADDRESS>
ISI Foundation
Via Chisola 5, 10126 Torino, Italy

###### Abstract

This article describes an efficient method to learn distributed
representations, also known as _embeddings_. This is accomplished minimizing
an objective function similar to the one introduced in the Word2Vec algorithm
and later adopted in several works. The optimization computational bottleneck
is the calculation of the _softmax_ normalization constants for which a number
of operations scaling quadratically with the sample size is required. This
complexity is unsuited for large datasets and _negative sampling_ is a popular
workaround, allowing one to obtain distributed representations in linear time
with respect to the sample size. Negative sampling consists, however, in a
change of the loss function and hence solves a different optimization problem
from the one originally proposed. Our contribution is to show that the sotfmax
normalization constants can be estimated in linear time, allowing us to design
an efficient optimization strategy to learn distributed representations. We
test our approximation on two popular applications related to word and node
embeddings. The results evidence competing performance in terms of accuracy
with respect to negative sampling with a remarkably lower computational time.

Keywords: Representation learning, softmax approximation, negative sampling,
Word2Vec, Node2Vec

##  1 Introduction

The ability of assessing the similarity of objects is one of the most crucial
roles of machine learning (Bellet et al., 2015). This task is particularly
hard to accomplish when dealing with complex entities – such as the words
appearing in a text – that do not have an easily defined encoding. Distributed
representations (or _embeddings_) allow one to map them into points of a high
dimensional space, providing an answer to this problem (Bengio et al., 2013).
In the aforementioned example, an effective representation maps _similar_
words into nearby vectors.

A very relevant contribution to the field is the Word2Vec algorithm proposed
by Mikolov _et al._ in (Mikolov et al., 2013) that owes its success to the
efficiency, scalability and flexibility. Since its proposal, in fact, several
works improved it (Bojanowski et al., 2017) or drew inspiration from it to
learn distributed representations of sentences and documents (Mikolov et al.,
2013), graphs (Perozzi et al., 2014; Grover and Leskovec, 2016; Gao et al.,
2019; Rozemberczki et al., 2019; Nickel and Kiela, 2017; Narayanan et al.,
2017), time (Kazemi et al., 2019), temporal contact sequences (Goyal et al.,
2020; Rahman et al., 2018; Nguyen et al., 2018; Sato et al., 2021; Torricelli
et al., 2020), biological entities (Du et al., 2019; Ng, 2017), tweets
(Dhingra et al., 2016) and higher order interactions (Billings et al., 2019)
among others. These algorithms build upon two steps: (1) for each
$i\in\mathcal{V}$ – the set of entities – a sample
$\mathcal{S}_{i}\subseteq\mathcal{V}$ is defined; (2) the set of embedding
vectors $\\{\bm{x}\\}$ in $\mathbb{R}^{d}$ is obtained maximizing the log-
likelihood $\mathcal{L}$,

$\displaystyle\mathcal{L}(\\{\bm{x}\\})$
$\displaystyle=\sum_{i\in\mathcal{V}}\sum_{j\in\mathcal{S}_{i}}{\rm
log}\leavevmode\nobreak\ p_{\bm{x}}(j|i)$ $\displaystyle p_{\bm{x}}(j|i)$
$\displaystyle=\frac{e^{\bm{x}_{i}^{T}\bm{x}_{j}}}{\sum_{k\in\mathcal{V}}e^{\bm{x}_{i}^{T}\bm{x}_{k}}}:=\frac{e^{\bm{x}_{i}^{T}\bm{x}_{j}}}{Z_{i}}.$
(1)

Intuitively, if $j$ is often sampled from $i$, then the probability
$p_{\bm{x}}(j|i)$ is large and the embedding vectors $\bm{x}_{i},\bm{x}_{j}$
become closer, obtaining a similar representation to related entities. The
known limitation of the optimization of this cost function resides in the
partition function $Z_{i}$, requiring $\mathcal{O}(n)$ operations to be
computed, where $n=|\mathcal{V}|$ is the number of items to be embedded. Since
$Z_{i}$ has to be computed for all $i\in\mathcal{V}$, optimizing $\mathcal{L}$
requires $\mathcal{O}(n^{2})$ operations, making it unsuited to approach large
datasets. One of the most commonly adopted workarounds is _negative sampling_
(Mikolov et al., 2013) consisting in replacing $\mathcal{L}$ with the cost
function $\tilde{\mathcal{L}}$ of a weighted logistic principal component
analysis (Landgraf and Bellay, 2017). Its minima should be close to those of
$\mathcal{L}$ but it has the advantage of being computed in a linear time with
respect to $n$. Negative sampling has been extensively adopted but, to date,
it is still only poorly understood from a theoretical standpoint and some
studies evidenced its potential weaknesses (Chen et al., 2018; Qin et al.,
2016; Landgraf and Bellay, 2017; Mimno and Thompson, 2017; Kamigaito and
Hayashi, 2021).

In this paper we propose an approximation method that allows one to estimate
all the $Z_{i}$’s in $\mathcal{O}(n)$ operations. Adding a regularization term
to the loss function of Equation (1), we formulate the problem in a very
general way so that our proposed algorithm can be easily adopted in a variety
of contexts. In particular, we suppose there exists a set of sparse sampling
probability vectors $\\{\bm{p}^{(i)}\\}_{i\in\mathcal{V}}$ in $\mathbb{R}^{n}$
encoding a measure of proximity between $i$ and the other elements of
$\mathcal{V}$. For instance, in the case of words, $p^{(i)}_{j}$ may indicate
the probability that $j$ is in the same sentence as $i$. The design of the set
$\\{\bm{p}^{(i)}\\}_{i\in\mathcal{V}}$ is problem-dependent and beyond the
scope of this article in which we describe an efficient framework to create a
distributed entity representation _given_ the probability set. Our main
contributions are:

  1. 1.

We provide a theoretical framework showing that _all_ the partition functions
$\\{Z_{i}\\}_{i\in\mathcal{V}}$ can be accurately estimated in
$\mathcal{O}(n)$ operations. Based on this result, we describe an efficient
strategy to obtain distributed representations.

  2. 2.

We show that the pure Python implementation of our algorithm is competitive or
outperforms in terms of speed an optimized C implementation of Word2Vec.111We
refer to the _Gensim_ implementation (Rehurek and Sojka, 2011), available at
radimrehurek.com/gensim/ run in parallel on $8$ cores.

  3. 3.

We display two use cases of our algorithm applied to node and words
embeddings, proving its efficacy also in terms of accuracy.

  4. 4.

We provide Python version of our algorithm at
github.com/lorenzodallamico/EDRep.

#### Notations

  1. $\bullet$

Sets are denoted with the font $\mathcal{A}$. The cardinality of a set is
denoted with $|\mathcal{A}|$.

  2. $\bullet$

Matrices are denoted with a standard capital letter, $A$. The identity matrix
of size $n$ is denoted with $I_{n}$.

  3. $\bullet$

The notation $a_{n}=O_{n}(b_{n})$ implies that there exists a finite $c$ so
that $\displaystyle\lim_{n\to\infty}\frac{a_{n}}{b_{n}}=c$. Similarly
$a_{n}=o_{n}(b_{n})$ implies that $\lim_{n\to\infty}\frac{a_{n}}{b_{n}}=0$.

  4. $\bullet$

Vectors are indicated in boldface (e.g. $\bm{a}$), while vector elements and
scalars in standard font (e.g. $a_{i}$). The all-ones vector of size $n$ is
denoted with $\bm{1}_{n}$.

##  2 Main result

###  2.1 Problem position

Consider a set of items $\mathcal{V}$ with $n=|\mathcal{V}|$. For all
$i\in\mathcal{V}$, we let $\bm{p}^{(i)}\in\mathbb{R}^{n}$ be a probability
distribution defined on $\mathcal{V}$, encoding similarity: a large
$p_{j}^{(i)}$ implies high affinity between $i$ and $j$. We further define
$\bm{p}_{0}\in\mathbb{R}^{n}$ as a baseline $\bm{p}^{(i)}$ should be compared
to. We aim at giving a _distributed representation_ of
$\\{\bm{p}^{(i)}\\}_{i\in\mathcal{V}}$, _i.e._ a vector embedding
$\\{\bm{x}_{i}\\}_{i\in\mathcal{V}}$ that preserves structural properties
encoded in the probability distributions. We equivalently denote the set
$\\{\bm{x}_{i}\\}_{i\in\mathcal{V}}$ as $X\in\mathbb{R}^{n\times d}$, a matrix
storing the embedding vectors in its rows. We define $X$ as the solution to
the following optimization problem:

$\displaystyle X$ $\displaystyle=\underset{Y\in\mathcal{C}_{n\times d}}{\rm
arg\leavevmode\nobreak\ min}\leavevmode\nobreak\
\sum_{i\in\mathcal{V}}\left[D_{\rm
KL}\left(\bm{p}^{(i)}\|\bm{q}^{(i)}_{Y}\right)-D_{\rm
KL}\left(\bm{p}_{0}\|\bm{q}^{(i)}_{Y}\right)\right]$
$\displaystyle\mathcal{C}_{n\times d}$
$\displaystyle=\\{Y\in\mathbb{R}^{n\times d}\leavevmode\nobreak\
:\leavevmode\nobreak\ \|\bm{y}_{i}\|=f_{i}\leavevmode\nobreak\
\forall\leavevmode\nobreak\ i\in\mathcal{V}\\}$ $\displaystyle q_{Y,j}^{(i)}$
$\displaystyle=\frac{1}{Z_{i}}e^{\bm{y}_{i}^{T}\bm{y}_{j}};\hskip
28.45274ptZ_{i}:=\sum_{k\in\mathcal{V}}e^{\bm{y}_{i}^{T}\bm{y}_{k}},$ (2)

where $D_{\rm KL}$ denotes the Kullback-Leibler divergence between two
distributions and $\\{f_{i}\\}_{i\in\mathcal{V}}$ are inputs of the
optimization problem. In words, we adopt a variational approach, approximating
$\bm{p}^{(i)}$ with the distribution $\bm{q}_{Y}^{(i)}$ and optimize over the
embeddings under the constraint that they have a predefined norm,
$f_{i}\in\mathbb{R}^{+}$. The term $D_{\rm KL}(\bm{p}_{0}\|\bm{q}_{Y}^{(i)})$
serves as a regularizer and “promotes” only the probability part exceeding the
baseline, thus preventing the embedding to encode a trivial piece of
information. This type of approach is – for instance to define the
_modularity_ (Newman, 2006) – to encode a null model.

######  Remark 1 (Extension to an asymmetric setting).

In Equation (2) and in the remainder, we only consider distributions
$\bm{p}^{(i)}$ that are defined and indexed over the same set, $\mathcal{V}$.
It is however possible to consider an asymmetric scenario in which the entries
$p^{(i)}_{a}$ are defined for $i\in\mathcal{V}$ and $a\in\mathcal{W}$. In this
case, one needs to define two embedding matrices
$X\in\mathbb{R}^{|\mathcal{V}|\times d}$ and
$Y\in\mathbb{R}^{|\mathcal{W}|\times d}$ and let

$\displaystyle
q_{X,Y,a}^{(i)}=\frac{1}{Z_{i}}e^{\bm{x}_{i}^{T}\bm{y}_{a}};\hskip
28.45274ptZ_{i}=\sum_{a\in\mathcal{W}}e^{\bm{x}_{i}^{T}\bm{y}_{a}}.$

All the framework described below has a straightforward generalization to this
setting.

Given Equation (2) our goal is to efficiently compute $Z_{i}$ for all
$i\in\mathcal{V}$ to obtain the embedding $X$ that solves the optimization
problem. We now proceed describing our main result on an efficient way to
estimate the normalization constants.

###  2.2 Efficient estimation of $Z_{i}$

In the following, we interpret the vectors
$\\{\bm{x}_{i}\\}_{i\in\mathcal{V}}$ as random variables. Under proper
assumptions, $Z_{i}/n$ is concentrated around its expected value.222A random
variable $Z$ is said to be concentrated around the pivot $a$ if
$\forall\leavevmode\nobreak\ t>0\leavevmode\nobreak\ :\leavevmode\nobreak\
\mathbb{P}(|Z-a|>t)\leq\alpha(t)$, where $\alpha\leavevmode\nobreak\
:\leavevmode\nobreak\ \mathbb{R}^{+}\to[0,1]$ is a non-increasing left
continuous function (Louart and Couillet, 2018). We then determine the
condition ensuring its concentration that only relies on the distribution of
the embedding vector norms.333Note that the $f_{i}$ values are fixed but they
can be seen as a realization of a random variable. To compute the expectation
of $Z_{i}$, an analytical expression thereof is needed and we empirically show
that a good approximation can be obtained assuming that
$\bm{x}_{i}^{T}\bm{x}_{j}$ is distributed as a mixture of $\kappa$ Gaussians,
where $\kappa$ is a input of our estimation algorithm. Finally, we describe a
practical way to efficiently compute $Z_{i}$ based on our results and detail
an algorithm to solve the optimization problem of Equation (2).

Let us first set sufficient conditions for the concentration of $Z_{i}/n$ to
its expected value.

######  Proposition 1 (Concentration of $Z_{i}$).

Let $X\in\mathbb{R}^{n\times d}$ be a random matrix so that
$\|\bm{x}_{i}\|=f_{i}$ where $f_{i}\in\mathbb{R}^{+}$ are random i.i.d.
variables distributed $f_{i}\sim p_{f}$. Suppose that the expected value of
$f_{i}=O_{n}(1)$ and that for all $\omega$ there is a $t^{*}>0$ so that for
$t>t^{*}$, $p_{f}(x>t)<2e^{-\omega t}$. Then, under these hypotheses, for all
large $n$, $Z_{i}/n$ is concentrated around its expected value.

The Proof of Proposition 1 is given in Appendix A. This result requires that
the distribution of the embedding vectors should not be too broad to guarantee
concentration. Making some mild assumptions on the distribution of the
$\\{\bm{x}_{i}\\}_{i\in\mathcal{V}}$ vectors, we now provide our main result
to approximate the $Z_{i}$ values.

######  Proposition 2 (Approximation of $Z_{i}$).

Let $\\{\bm{x}_{i}\\}_{i\in\mathcal{V}}$ be a set of random vectors in
$\mathbb{R}^{d}$ and let $\ell:\mathcal{V}\to[\kappa]$ be a labeling function.
Suppose that for all $j\in\mathcal{V}$ so that $\ell(j)=a$ the corresponding
embedding vectors are distributed as $\bm{x}_{i}\overset{\rm
i.i.d}{\sim}P_{a}(\bm{\mu}_{a},\Omega_{a})$, where
$\bm{\mu}_{a}\in\mathbb{R}^{d}$ and $\Omega_{a}\in\mathbb{R}^{d\times d}$ are
the mean and covariance matrix of $P_{a}$. Then, for $i$ fixed and
$\ell(j)=a$,

$\displaystyle\bm{x}_{i}^{T}\bm{x}_{j}\sim\mathcal{N}(\bm{x}_{i}^{T}\bm{\mu}_{a},\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}),$

where $\mathcal{N}$ is a Gaussian random variable. Letting $\pi_{a}$ be the
fraction of items in class $a$, we get

$\displaystyle\frac{Z_{i}}{n}\approx\sum_{a=1}^{\kappa}\pi_{a}e^{\bm{x}_{i}^{T}\bm{\mu}_{a}+\frac{1}{2},\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}}.$
(3)

Figure 1: Empirical validation of Proposition 2. _Top row_ : cumulative
density function (CDF) of the relative error $\epsilon$ between the true an
estimated value of $Z_{i}$ for $1000$ randomly selected indices $i$ and the
order of the mixture of Gaussian distribution $\kappa=1,\dots,8$, color-
encoded. _Bottom row_ : histogram of $\bm{x}_{i}^{T}\bm{x}_{j}$ for a fixed
$i$ in green and in black the mixture of Gaussian fit, obtained for
$\kappa=5$.

To best understand Proposition 2, let us comment both its hypotheses and
conclusion. Since we do not have an explicit expression of the embedding
vectors $\\{\bm{x}_{i}\\}_{i\in\mathcal{V}}$ distribution, we suppose that we
can write it as a mixture of distributions with finite norm and covariance. If
this conditions is verified, then $\bm{x}_{i}^{T}\bm{x}_{j}$ is distributed as
a Gaussian random variable, for which the expectation can be computed
analytically. Under general hypotheses, we did not prove the scalar product to
converge towards a Gaussian distribution but, in Appendix A, however, we
provide a theoretical support to this claim, further assuming that all
embedding dimensions are independent. The consequence of Proposition 2 is that
we can efficiently compute _all_ the $Z_{i}$ with $\mathcal{O}(n\kappa d^{2})$
operations: $\mathcal{O}(n\kappa d)$ operations are needed to compute a label
assignment with _k-means_ (Lloyd, 1982; Christopher D. Manning, 2010), for
instance, while the $Z_{i}$’s, given the parameters
$\\{\bm{\mu}_{a},\Omega_{a},\pi_{a}\\}_{a=1,\dots,\kappa}$ are obtained in
$O(n\kappa d^{2})$ operations.

We now proceed giving an empirical validation to Proposition 2. We considered
$4$ datasets taken from the NLPL word embeddings repository (Kutuzov et al.,
2017), representing word embeddings obtained with different algorithms:

  1. $0$.

_British National Corpus_ , trained with _Continuous Skip-Gram_ (Mikolov et
al., 2013), $n=163.473$, $d=300$;

  2. $7$.

_English Wikipedia Dump 02/2017_ , trained with _Global Vectors_ (Pennington
et al., 2014), $n=273.930$, $d=300$;

  3. $9$.

_English Wikipedia Dump of 02/2017_ , trained with _fastText Skipgram_ (Joulin
et al., 2017), $n=273.930$, $d=300$;

  4. $224$.

_Ukrainian CoNLL17 corpus_ , trained with _Continuous Bag-of-Words_ (Mikolov
et al., 2013), $n=99.884$, $d=200$.

For each of these datasets we (1) set the average of the embedding norm to
$1$, (2) sample at random $1000$ words and (3) compute the corresponding exact
and the estimated $Z_{i}$ values for different approximation orders $\kappa$.
Figure 1 shows the results of this procedure: in the top row the cumulative
density function (CDF) of the relative error is displayed for several values
of $\kappa$. As expected, the precision is an increasing function of $\kappa$,
but very good approximations are obtained on all datasets for $\kappa=1$
already. Note that in this case, the clustering step can be omitted and is,
therefore, of particular interest. The bottom row shows the empirical and the
predicted distribution of the $\bm{x}_{i}^{T}\bm{x}_{j}$ for a randomly chosen
fixed $i$ and confirms the goodness of the approximation.

As a remark, the computational complexity of our proposed method scales
quadratically with $d$ because we need to compute the covariance matrix
$\Omega_{a}$ for all $a$. If one assumes that all dimensions are independent,
then the set of $\Omega_{a}$ becomes diagonal and it can be computed in
$\mathcal{O}(n\kappa d)$ operations. We empirically verified that under this
assumption results are worse but, for the $4$ considered datasets, we still
obtain good approximations, with a relative error of the order of $1\%$ or
smaller for most items. Let us further point out that in the _random matrix
regime_ in which $n$ and $d$ are comparable, the empirical covariance matrix
is not a good estimator and more sophisticated techniques should be adopted
(Tiomoko et al., 2019). In the code we provide, however, we only implement the
sample covariance matrix to estimate the set of $\Omega_{a}$ since in no case
$n$ and $d$ are comparable in our settings. On the other hand, we let the user
choose whether to compute the full covariance matrix or only its diagonal.

###  2.3 EDRep: efficient distributed representations

Let us now detail the main steps needed to translate the result of Proposition
2 into a practical algorithm to produce efficient distributed representations.

#### Optimization strategy

According to Equation (2), the embedding matrix $X\in\mathbb{R}^{n\times d}$
is obtained by optimizing a cost function $\mathcal{L}(X)$ that, after some
trivial calculations, can be rewritten in the following form:

$\displaystyle\mathcal{L}(X)=-\sum_{i\in\mathcal{V}}\left[\sum_{j\in\mathcal{V}}\left(p^{(i)}_{j}-p_{0,j}\right)\bm{x}_{i}^{T}\bm{x}_{j}+{\rm
log}(Z_{i})\right].$

We now exploit the result of Proposition 2 and make the notation more compact.
Let $P\in\mathbb{R}^{n\times n}$ be a matrix so that $P_{ij}=p^{(i)}_{j}$.
Then, $\mathcal{L}(X)$ can be rewritten as:

$\displaystyle\mathcal{L}(X)$ $\displaystyle=-{\rm
tr}\left[X^{T}\left(P-\bm{1}_{n}\bm{p}_{0}^{T}\right)X\right]+\sum_{i\in\mathcal{V}}{\rm
log}\leavevmode\nobreak\ \sum_{a=1}^{\kappa}\underbrace{\pi_{a}{\rm
exp}\left\\{\bm{x}_{i}^{T}\bm{\mu}_{a}+\frac{1}{2}\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}\right\\}}_{Z_{ia}}$
$\displaystyle:=-{\rm
tr}\left[X^{T}\left(P-\bm{1}_{n}\bm{p}_{0}^{T}\right)X\right]+\sum_{i\in\mathcal{V}}{\rm
log}\leavevmode\nobreak\ (Z\bm{1}_{\kappa})_{i}$ (4)

where we also introduced the matrix $Z\in\mathbb{R}^{n\times\kappa}$ with
elements $Z_{ia}$ for convenience. The proposed strategy to optimize
$\mathcal{L}$ is summarized in Algorithm 1 and consists of alternating a
gradient descent step with a forgetting factor and the normalization. Letting
$M\in\mathbb{R}^{d\times\kappa}$ have the
$\\{\bm{\mu}_{a}\\}_{a=1,\dots,\kappa}$ values in its rows, the
$i\in\mathcal{V}$ and $\alpha\in\\{1,\dots,d\\}$ gradient component reads

Input: $\\{\bm{p}^{(i)}\\}_{i\in\mathcal{V}}$, set of probabilities;
$\bm{p}_{0}\in\mathbb{R}^{n}$ null model probability; $d$, embedding
dimension; $\bm{f}$ vector of embedding norms;
$\bm{\ell}\in\\{1,\dots,\kappa\\}^{n}$ node label vector; $\eta$, largest
learning rate; n_epochs, number of training epochs

Output: $X^{n\times d}$, embedding matrix

begin

$X\leftarrow$ initialize the embedding matrix;

for  _$1\leq i\leq$ n_epochs_ do

$\bm{\mu}_{a=1,\dots,\kappa}\leavevmode\nobreak\
,\Omega_{a=1,\dots,\kappa}\leftarrow$ update the parameters;

$G\leftarrow$ gradient matrix as in Equation (2.3);

$\eta_{\rm c}\leftarrow$ largest learing rate so that
$\forall\leavevmode\nobreak\ i,\bm{x}^{T}_{i}(t+1)\bm{x}_{i}(t)>0$ (Equation
6);

$X\leftarrow(1-\eta_{\rm c})X-\eta_{\rm c}G$; gradient descent step;

$\|\bm{x}_{i}\|\leftarrow f_{i}\leavevmode\nobreak\
\forall\leavevmode\nobreak\ i$; normalize the embedding vectors

end for

end

Algorithm 1 EDRep optimization

$\displaystyle
G_{i\alpha}=-[(P+P^{T})X]_{i\alpha}+[\bm{1}_{n}(\bm{p}_{0}^{T}X)+\bm{p}_{0}(\bm{1}_{n}^{T}X)]_{i\alpha}+\frac{1}{Z_{i}}\left([ZM]_{i\alpha}+\sum_{a=1}^{\kappa}Z_{ia}[X\Omega_{a}]_{i\alpha}\right).$
(5)

At each iteration, the embedding learns only from the gradient component which
is orthogonal to the current embedding. As the minimum is approached, this
component decreases in magnitude, hence determining the algorithm convergence.
To improve the stability, we impose for all $i$ and all $t$ that
$\bm{x}^{T}_{i}(t+1)\bm{x}_{i}(t)>0$, where $\bm{x}_{i}(t)$ indicates the
embedding vector at the $t$-th iteration. From a straightforward calculation,
one can see that, provided $f_{i}^{2}+\bm{g}_{it}\bm{x}_{it}>0$, the learning
rate $\eta_{\rm c}$ must satisfy the following condition for all $i$:

$\displaystyle\eta_{\rm
c}<\frac{f_{i}^{2}}{(f_{i}^{2}+\bm{g}_{i,t}\bm{x}_{i,t})}.$ (6)

If $f_{i}^{2}+\bm{g}_{it}\bm{x}_{it}<0$, any positive $\eta_{\rm c}$ can be
used, and we preset the value $\eta$ to be adopted.

#### Computational complexity of the algorithm

As an important remark, we must mention that Algorithm 1 requires the node
labeling vector $\bm{\ell}$ as an input which is, generally, not known in
advance. To cope with this problem, we run the optimization with $\kappa=1$
(for which $\bm{\ell}=\bm{1}_{n}$), obtain an embedding $X$ and compute a new
$\bm{\ell}$ as the result $\kappa$-class clustering on $X$.

The cost of this optimization procedure requires $\mathcal{O}(n\kappa d)$
operations to compute $\bm{\ell}$, $\mathcal{O}(n\kappa d^{2})$ for the “Z
part” of the gradient, and $\mathcal{O}(Ed)$ to compute the product $PX$,
where $E$ denotes the number of non-zero entries of $P$. This complexity
becomes prohibitive when $E\gg n$, but this is not the case in typical
settings and the product can be performed very efficiently. Nonetheless, even
when $E$ is very large, if $P$ can be written as the product (or sum of
products) of sparse matrices, $PX$ can still be computed very efficiently. In
fact, let $P=P_{m}\cdot P_{m-1}\dots P_{1}$ for some positive $m$, then $PX$
can be obtained without computing $P$, simply taking the products from right
to left:

$\displaystyle PX=(\overset{{\Longleftarrow}}{P_{m}\cdot
P_{m-1}\cdot\dots\cdot P_{1}}X),$

thus speeding up the computational bottleneck of our algorithm. In our
implementation we explicitly consider this representation of $P$ as an input.
When a non-factorized dense matrix $P$ is provided, one could envision to
adopt a method such as the one presented in (Le Magoarou and Gribonval, 2016)
to approximate a dense matrix $P$ with the product of sparse matrices to
considerably speed up the algorithm.

#### The role of the vector of norms

The vector of norms $\bm{f}$ is used as an input of our algorithm and it is
necessary to preserve its numerical stability. An effective choice for this
parameter depends on how the embedding is used in the successive steps. After
convergence is reached, in fact, we may write

$\displaystyle P_{ij}\approx e^{f_{i}f_{j}{\rm
cosSim}(\bm{x}_{i},\bm{x}_{j})},$

where ${\rm
cosSim}(\bm{x},\bm{y})=\frac{\bm{x}^{T}\bm{y}}{\|\bm{x}\|\|\bm{y}\|}$ denotes
the cosine similarity between two vectors. $P_{ij}$, that is the probability
to sample $j$ from $i$, may be large because of their affinity as well as
because of the popularity.444For instance, popularity can be the number of
times a word appears in a text. If $f_{i}$ is chosen to be an increasing
function of $i$-th popularity, one allows instances in which $P_{ij}$ may be
large for popular items with a low affinity. On the opposite, if it is a
decreasing function of popularity, one favors the alignment with popular
items. By default, in our implementation, we leave $\bm{f}=\bm{1}_{n}$.

We now proceed showing two use cases of our algorithm, applied to the
generation of node and word distributed representations.

##  3 Two use cases

In this section we analyze two use cases for two popular applications of
distributed representations: node and words embeddings. In the former case,
the goal is to determine a distributed representation of the vertices of a
graph that preserves its structural properties (Cai et al., 2018; Goyal and
Ferrara, 2018). On the other hand, distributed word representations are one of
the most popular applications of this class of algorithms that allow one to
encode words semantic similarity (Almeida and Xexéo, 2019). In both cases,
pairs of _similar_ items (nodes or words) should correspond to embedding
vectors with a high cosine similarity.

######  Remark 2.

To test our algorithm on practical applications, we must define the set
$\\{\bm{p}^{(i)}\\}_{i\in\mathcal{V}}$. In the remainder, we adopt simple
strategies to obtain this set and show, nonetheless, that the EDRep approach
achieves competitive results in terms of performance. We would like to
underline, however, that the sampling probabilities choice is a hard and
problem-dependent task and addressing it optimally is beyond the scope of this
article. Moreover, we formulated EDRep as an _offline_ algorithm while
Word2Vec is _online_. This difference influences the algorithm performances
but not the essence our main result – Proposition 2 – that provides an
efficient estimation of the partition functions _regardless_ of the details of
the implementation strategy.

###  3.1 Node embeddings

A graph is a tuple $\mathcal{G}(\mathcal{V},\mathcal{E})$ where $\mathcal{V}$
is the set of _nodes_ and $\mathcal{E}$ the set of edges connecting nodes
(Newman, 2003; Borgatti and Halgin, 2011). Graphs are a versatile tool to
model interacting systems and have a range spectrum of application, including
social networks (Scott, 1988; Wasserman and Faust, 1994), technological
networks (Amaral et al., 2000), spatial networks (Barthélemy, 2011),
biological networks (Alm and Arkin, 2003) and many others. A graph
$\mathcal{G}$ can be easily represented in matrix form with the _adjacency
matrix_ $A\in\mathbb{R}^{n\times n}$, so that $A_{ij}=1$ if
$(ij)\in\mathcal{E}$ and equals zero otherwise. Here and in the remainder we
only focus on undirected ($A_{ij}=A_{ji}$) and unweighted graphs
($A_{ij}\in[0,1]$). Our algorithm can, however, easily deal with the
challenging directed and weighted cases, provided that an appropriate sampling
strategy is defined.

One of the most relevant problems in graph learning is _community detection_
that can be informally defined as the task of determining a non-overlapping
node partition, unveiling more densely connected groups of nodes (Fortunato
and Hric, 2016). In the following – similarly to what done in _spectral
clustering_ (Von Luxburg, 2007) – we adopt the EDRep algorithm to produce a
node embedding and then perform clustering in the embedded space. Our choice
of $P$ in this case is the _random-walk Laplacian_ , obtained from $A$ and
imposing that its rows sum up to $1$. The entry $P_{ij}$ of this matrix is the
limiting probability that a random walker on $\mathcal{G}$ goes from node $i$
to $j$ in one step. We further set $\bm{p}_{0}=\bm{1}_{n}/n$, that is
equivalent to comparing $\bm{p}^{(i)}$ with a uniform prior.

With these choices, we generate a distributed node representation $X$ and then
cluster the nodes of $\mathcal{G}$ into groups. To evaluate the accuracy we
consider synthetic graphs generated from the _degree-corrected stochastic
block model_ in which a ground truth label assignment is known and real-world
graphs with a class annotation. The accuracy is measured in terms of
_normalized mutual information_ (NMI), a score that ranges between $0$ (random
label assignment) and $1$ (perfect assignment).

#### Synthetic graphs

The _degree-corrected stochastic block model_ (DCSBM) (Karrer and Newman,
2011) is a random graph generative model capable to create graphs with a
community structure and a arbitrary degree distribution.555The degree of a
node is its number of connections. This model is convenient to study community
detection as it provides a clear definition of communities but it also allows
one to model the broad degree distribution, typically observed in real-world
graphs (Barabási and Albert, 1999).

######  Definition 1 (DCSBM).

Let $\omega:\mathcal{V}\to\\{1,\dots,q\\}$ be a class labeling function, where
$q$ is the number of classes. Let $\mathbb{P}(\omega_{i}=a)=q^{-1}$ and
consider two positive integers satisfying $c_{\rm in}>c_{\rm out}\geq 0$.
Further Let $\theta\sim p_{\theta}$ be a random variable that encodes the
intrinsic node connectivity, so that $\mathbb{E}[\theta]=1$ and
$\mathbb{E}[\theta^{2}]=O_{n}(1)$. For all $i\in\mathcal{V},\theta_{i}$ is
drawn independently at random from $p_{\theta}$. The entries of the graph
adjacency matrix are generated independently (up to symmetry) at random with
probability

$\displaystyle\mathbb{P}(A_{ij}=1)=\frac{\theta_{i}\theta_{j}}{n}\cdot\begin{cases}c_{\rm
in}\quad&{\rm if}\leavevmode\nobreak\ \omega(i)=\omega(j)\\\ c_{\rm
out}\quad&{\rm else}\end{cases}$

Figure 2: Community detection on synthetic graphs. _Left panel_ : NMI score as
a function of $\alpha$ (footnote 6) for a DCSBM graph (Definition 1) with
$n=30\leavevmode\nobreak\ 000$, $c=10$, $q=4,\kappa=1$ and varying values of
$c_{\rm out}$. The three lines correspond to three different embedding
techniques. The purple dots are the EDRep algorithm with $d=32$, the cyan
diamonds is DeepWalk with $d=32$, while the orange narrow diamonds are the
spectral clustering algorithm of (Dall’Amico et al., 2021). _Right panel_ :
corresponding computational time in logarithmic scale. _Both plots_ : averages
are over $10$ samples and the error bar width equals the standard deviation.

In words, nodes in the same community ($\omega(i)=\omega(j)$) are connected
with a higher probability than nodes in different communities. From a
straightforward calculation, the expected degree is
$\mathbb{E}[c_{i}]\propto\theta_{i}$. For this model, the goal is to infer the
node label assignment from a realization of $A$. It was shown that this is
theoretically feasible (in the large $n$ regime) if and only if a function of
the model parameters666Letting $c$ be the expected average degree
$c=\frac{c_{\rm in}+(k-1)c_{\rm out}}{k}$, $\alpha=(c-c_{\rm
out})\sqrt{\frac{\mathbb{E}[\theta^{2}]}{c}}$.  $\alpha>1$ (Gulikers et al.,
2015, 2017). In left panel of Figure 2 we show the NMI score obtained by
different algorithms as a function of $\alpha$, while the right panel shows
the computation time in logarithmic scale. We compare:

  1. $\bullet$

the spectral algorithm of (Dall’Amico et al., 2019, 2021) that was shown to be
nearly Bayes-optimal for this task,

  2. $\bullet$

an algorithm combining the EDRep optimization approach with $\kappa=1$ and _k-
means_. The vector of norms is set to $f_{i}\propto\sqrt{d_{i}}g(d_{i})$,
where $g(d_{i})$ is a Gaussian cut-off,

  3. $\bullet$

an algorithm combining DeepWalk (Perozzi et al., 2014) with _k-means_. We used
the C++ implementation available at github.com/thibaudmartinez/node2vec. The
algorithm is deployed with its default values.

The results show that the EDRep-based algorithm performs almost as well as the
optimal algorithm of (Dall’Amico et al., 2021) and a slight mismatch is only
observed for hard problems, _i.e._ for $\alpha$ approaching $1$. This is a
particularly challenging setting in which few algorithms are capable to
retrieve the community structure. The method of (Dall’Amico et al., 2021) was
explicitly designed to work well in this regime in which, instead, the
spectral algorithm based on the _random walk Laplacian_ (Shi and Malik, 2000)
(our choice of $P$) performs poorly. When compared to the DeepWalk approach,
our method generally yields better results for all $\alpha$.

On top of accuracy, the computational time gives the EDRep approach a decisive
advantage with respect to the C++ efficient implementation of DeepWalk being
more than $10$ times faster with a improvement observed also with respect to
the spectral clustering algorithm. The main advantage in this case is however
related to the algorithm’s computational complexity. The considered spectral
clustering algorithm runs in $\mathcal{O}(nq^{3})$ operations, while the
complexity of EDRep is independent of $q$.777A relevant setting is the one in
which $\kappa=q$, _i.e._ the order of the GMM approximation equals the number
of communities. It this case, EDRep’s complexity scales as
$\mathcal{O}(q^{2})$ or $\mathcal{O}(q)$ when considering diagonal matrices
$\Omega$. It has to be said, however, that for very large $q$ it seems
pointless to set $\kappa=q$ since $\kappa\approx 8$ gave very good estimates
in all cases we tested. Actually, while the estimation of $Z_{i}$ is more
accurate for large $\kappa$, we did not see any significant gain in the
overall performance when choosing $\kappa>1$. Consequently, for large $q$,
spectral clustering has to pay a significant price in terms of computational
efficiency and it becomes unpractical.

Name | $n$ | $q$ |  EDRep – NMI | NS – NMI |  EDRep – $t$ | NS – $t$
---|---|---|---|---|---|---
_amazon_ | 9 317 | 232 | 0.92 | 0.93 | 4 s | 17 s
_youtube_ | 26 931 | 608 | 0.60 | 0.65 | 12 s | 72 s
_dblp_ | 60 530 | 252 | 0.52 | 0.51 | 26 s | 180 s
_livejournal_ | 64 504 | 1662 | 0.92 | 0.92 | 34 s | 196 s
Table 1: Community detection on real graphs. For each graph we report: 1) the
name; 2) the size $n$; 3) the number of communities $q$; 4,5) the NMI score
for EDRep and negative sampling (NS) and 6,7) the corresponding computational
time in seconds. The results are obtained averaging over $5$ samples.

#### Real-world graphs

We further test our algorithm on more realistic data available from the
Stanford Network Analysis Platform (Leskovec and Krevl, 2014) that come with a
community label annotation. In most cases, for the considered datasets, each
node may be associated to more than one community. For simplicity, in order to
keep the same framework as the one detailed above, each node was associated to
the largest community it belongs to. In Table 1 we report the results
summarizing the network properties, the obtained NMI and the computational
time. In this case, spectral clustering is not considered since its complexity
does not allow us to obtain a partition in a reasonable time. Similar to what
already mentioned for DCSBM-generated graphs, the EDRep algorithm obtains
competing performances with respect to negative sampling at a much lower
computational time.

###  3.2 Word embeddings

Word embeddings are one of the most prominent applications for the class of
algorithms EDRep belongs to. The goal is to encode in a vector-representation
the words semantic similarity as inferred from a text (Selva Birunda and
Kanniga Devi, 2021). For instance, the words _flower_ and _petal_ are closely
related and should be represented by vectors with a high cosine similarity.
This representation is mathematically convenient and opens up to a great span
of applications such as those discussed in (Wang et al., 2018; Khattak et al.,
2019) to provide some examples.

Unlike graphs, for word representations the matrix $P$ is not given and it has
to be built from the text. To do so we define a co-occurrence matrix
$W\in\mathbb{R}^{n\times n}$, where $n$ denotes the vocabulary size and
$W_{ij}$ is the number of times the word $j$ appeared in the _context_ of word
$i$, _i.e._ within a maximum distance from $i$. For example, in Albus
Dumbledore’s sentence

_It takes a great deal of bravery to stand up to our enemies,
but just as much to stand up to our friends_

the words _stand, up, to, our, but, just, as, much_ appear in a context of
size $4$ around _enemies_. Given $W$ and normalizing its rows we obtain our
choice for the matrix $P$ in which large entries correspond to pairs of words
frequently appearing close to one-another. As one can observe in the example
above, however, some context-independent words tend to appear with high
frequency without carrying much information. This is the case of “ _to_ ” that
has four occurrences in the short sentence. Prior to building the co-
occurrence matrix, we first proceed with the standard cleaning procedure
introduced in (Mikolov et al., 2013) in which each word is dropped with a
probability $1-\sqrt{e^{-5}/\omega_{i}}$, where $\omega_{i}$ denotes the
number of appearances of word $i$ in the text. For the $\bm{p}_{0}$ vector, we
drew inspiration from negative sampling and used
$p_{0,i}\propto\omega_{i}^{\gamma}$ with $\gamma=0.75$. The embedding norms
are chosen to be a slow increasing function of the word frequency as well and
in particular $f_{i}\propto{\rm log}(\omega_{i}+1)$.

We trained the EDRep algorithm on three datasets: million, text8 and text9.
These files are obtained from a dump of English Wikipedia performed on March
3rd, 2006 (mattmahoney.net/dc/textdata.html) and represent the first $10^{6}$
words, $10^{8}$ and $10^{9}$ bytes of this dump, respectively. Surprisingly,
we observed that a context window of size $1$ gives best performances for our
algorithm, but in our implementation we leave it as an input parameter. Figure
3a shows the histogram of the non-zero entries of the matrix $P$. As expected
(Zipf, 2013), most of them correspond to very rare occurrences. The
computational complexity of EDRep scales linearly with the number of non-zero
entries and very small values only marginally change the result, but majorly
affect efficiency. For this reason we adopted a sparsification strategy,
keeping only the $100$ largest entries per row in $P$ and observed small
losses in the performance.

Figure 3: Test on word embeddings. (a) _distribution of the non-zero entries
of the $P$ matrix_: the plot is shown for _text8_ and is plotted in
logarithmic scale on the $y$ axis. (b, c, d): performance achieved on
different benchmarks indicated on the $x$ axis for the million, text8, text9
datasets. The results of the EDRep algorithm (orange dot) are compared with
Word2Vec (green cross). The plotted quantity is average Spearman correlation
between the benchmark similarity scores and those obtained with the two
embedding methods. Averages are over $20$ samples for million and text8 and
over $10$ samples for text9. Reading the plots from left to right (b, c, d)
the _online_ negative sampling approach profits from longer training unlike
the _offline_ EDRep method, as detailed in Section 3.2 and motivated in Remark
2.

A standard way to evaluate word embeddings is to check if the cosine
similarity between pairs of words relates to its semantic affinity. To obtain
a quantitative performance measure, we used the datasets made available by
(Jastrzebski et al., 2017) in which a list of pairs of words $(a,b)$ are
associated to a score $s_{ab}$. For each of these pairs, we compute the cosine
similarity obtained from the embeddings and use the Spearman correlation
coefficient to quantify the agreement between the two measures. Higher
correlation values indicate a good match with manually annotated data. The
results are compared in Figure 3 with the Word2Vec algorithm using $d=128$ as
embedding dimension. In particular, Figure 3b evidences a neat advantage of
the EDRep approach when dealing with a small text, as it allows us to learn
reasonable word representation from a very short text already. Increasing the
input file size, in text8 (Figure 3c) we observe very similar performances for
the two algorithms, while on text9 (Figure 3d), higher performances are
obtained by Word2Vec on all benchmarks and EDRep performs essentially as well
as on text8.

These results can be easily interpreted when looking more closely to the two
algorithms. EDRep must first learn a matrix text representation and then a
distributed representation thereof. Consequently, the algorithm training time
is independent of the text length. On the opposite, Word2Vec, learns the
embeddings in an _online_ fashion and longer texts imply more training epochs.
These two algorithms are indeed different and an online training method may
effectively be more powerful than an offline one for this specific task. With
this perspective, given that our main result of Proposition 2 would still
hold, one could envision to adapt EDRep to different training strategies.
Moreover, EDRep’s ability to obtain significant representation from sparse
data could be exploited to pre-train the Word2Vec embeddings on a short text
to speed up convergence.

Coming to the computation time, we confirm also in this case that our approach
is faster on all three texts by approximately a factor $2$, even though
Word2Vec is run in parallel on all available cores. One can notice that the
gain in computational efficiency is smaller than the one obtained in Section
3.1. This is because in that case a graph readily comes with a matrix
representation and Node2Vec has to generate a “text” from the graph.
Contrariwise, here the matrix must be obtained from the text and, with our
implementation, it requires approximately $50\%$ of the total time.

As a final remark, we would like to stress that the scores reported in Figure
3 are an only marginally objective way to evaluate a word embedding. In fact,
we tried other design choices of the matrix $P$ and in many cases we ended up
obtaining larger scores on some benchmarks but lower on others. For this
reason, we comfortably state that, once a proper evaluation metric is defined,
there is much room for improvement in the design of the $P$ matrix and these
results should be interpreted as an algorithm deployment example.

##  4 Conclusion

In this article we showed that, under proper conditions, an approximation of
the normalization constants of _softmax_ can be obtained in linear time. This
allowed us to optimize a loss function similar to the one originally proposed
in (Mikolov et al., 2013) to efficiently obtain distributed representations.
We detailed a theoretical framework justifying our claim that was validated on
real data and then translated our result into a practical algorithm. This
method relies on the definition of a probability matrix $P$ capable to encode
the similarity between the entities to be represented. To foster the impact of
our work, we showed two simple examples applied to node and word embeddings.

In both cases, our proposed approach proved to be faster than an optimized
negative sampling C implementation. This is a striking result, but we are
firmly convinced that there is a lot of room for improvement. In fact, our
implementation is fully coded in Python – that is known to be slower than C –,
it is not parallelized and it can reasonably gain a considerable advantage if
carefully optimized. Moreover, our main result fits in a larger frame than the
one described in this paper and can potentially be exploited with other
optimization strategies, as we already suggested along the lines.

Negative sampling appeared to be the natural benchmark of comparison of our
work, but we would like to stress that it consists in optimizing a different
cost function from the one optimized in EDRep. We noticed a remarkable
difference between the case in which an affinity matrix representation is
given (graphs) from the one in which it is not (texts). In the earlier case
all the information related to the corpus (the graph) is represented in the
form of a sparse matrix, hence the _offline_ implementation of EDRep works
optimally and significantly outperformed Node2Vec both in terms of speed and
of accuracy. In the latter case, when encoding the text into a matrix, all the
complex relations between the words should be captured by the matrix
structure. This is a non-trivial task that might penalize infrequent co-
occurrences between pairs of words. We believe, however, that a crucial
advantage of our approach lies in its high interpretability. For a given
practical deployment of our algorithm, a practitioner only needs to define a
sampling strategy that meaningfully encodes proximity for the problem at hand.
The simplicity of our framework makes it also easier to study than negative
sampling with the techniques described in (Jaffe et al., 2020), for instance.
Moreover, our main result holds in a wider range of settings with respect to
the one explored and it can be easily generalized to non symmetric, or non
normalized $P$ matrices.

Acknowledgments and Disclosure of Funding

LD acknowledges support from the the Lagrange Project of the ISI Foundation
funded by CRT Foundation and from the European Union’s Horizon 2020 research
and innovation programme under grant agreement No. 101016233 (PERISCOPE). EMB
acknowledges support from CRT Lagrange Fellowships in Data Science for Social
Impact of the ISI Foundation. The authors thank Ciro Cattuto for fruitful
discussions during the project conception.

## References

  * Alm and Arkin (2003) E. Alm and A. P. Arkin.  Biological networks.  _Current opinion in structural biology_ , 13(2):193–202, 2003.
  * Almeida and Xexéo (2019) F. Almeida and G. Xexéo.  Word embeddings: A survey.  _arXiv preprint arXiv:1901.09069_ , 2019.
  * Amaral et al. (2000) L. A. N. Amaral, A. Scala, M. Barthelemy, and H. E. Stanley.  Classes of small-world networks.  _Proceedings of the national academy of sciences_ , 97(21):11149–11152, 2000.
  * Barabási and Albert (1999) A.-L. Barabási and R. Albert.  Emergence of scaling in random networks.  _science_ , 286(5439):509–512, 1999.
  * Barthélemy (2011) M. Barthélemy.  Spatial networks.  _Physics reports_ , 499(1-3):1–101, 2011.
  * Bellet et al. (2015) A. Bellet, A. Habrard, and M. Sebban.  Metric learning.  _Synthesis lectures on artificial intelligence and machine learning_ , 9(1):1–151, 2015.
  * Bengio et al. (2013) Y. Bengio, A. Courville, and P. Vincent.  Representation learning: A review and new perspectives.  _IEEE transactions on pattern analysis and machine intelligence_ , 35(8):1798–1828, 2013.
  * Billings et al. (2019) J. C. W. Billings, M. Hu, G. Lerda, A. N. Medvedev, F. Mottes, A. Onicas, A. Santoro, and G. Petri.  Simplex2vec embeddings for community detection in simplicial complexes.  _arXiv preprint arXiv:1906.09068_ , 2019.
  * Bojanowski et al. (2017) P. Bojanowski, E. Grave, A. Joulin, and T. Mikolov.  Enriching word vectors with subword information.  _Transactions of the association for computational linguistics_ , 5:135–146, 2017.
  * Borgatti and Halgin (2011) S. P. Borgatti and D. S. Halgin.  On network theory.  _Organization science_ , 22(5):1168–1181, 2011\.
  * Cai et al. (2018) H. Cai, V. W. Zheng, and K. C.-C. Chang.  A comprehensive survey of graph embedding: Problems, techniques, and applications.  _IEEE transactions on knowledge and data engineering_ , 30(9):1616–1637, 2018.
  * Chen et al. (2018) L. Chen, F. Yuan, J. M. Jose, and W. Zhang.  Improving negative sampling for word representation using self-embedded features.  In _Proceedings of the Eleventh ACM International Conference on Web Search and Data Mining_ , pages 99–107, 2018.
  * Christopher D. Manning (2010) H. S. Christopher D. Manning, Prabhakar Raghavan.  Introduction to information retrieval.  _Natural Language Engineering_ , 16(1):100–103, 2010.
  * Dall’Amico et al. (2019) L. Dall’Amico, R. Couillet, and N. Tremblay.  Revisiting the bethe-hessian: improved community detection in sparse heterogeneous graphs.  In _Advances in Neural Information Processing Systems_ , pages 4039–4049, 2019.
  * Dall’Amico et al. (2021) L. Dall’Amico, R. Couillet, and N. Tremblay.  A unified framework for spectral clustering in sparse graphs.  _Journal of Machine Learning Research_ , 22(217):1–56, 2021.
  * Dhingra et al. (2016) B. Dhingra, Z. Zhou, D. Fitzpatrick, M. Muehl, and W. W. Cohen.  Tweet2vec: Character-based distributed representations for social media.  page 269 – 274, 2016.  doi: 10.18653/v1/p16-2044.  Cited by: 76; All Open Access, Green Open Access, Hybrid Gold Open Access.
  * Du et al. (2019) J. Du, P. Jia, Y. Dai, C. Tao, Z. Zhao, and D. Zhi.  Gene2vec: distributed representation of genes based on co-expression.  _BMC genomics_ , 20(1):7–15, 2019.
  * Feller (1968) W. Feller.  An introduction to probability theory and its applications.  Technical report, 1968.
  * Fortunato and Hric (2016) S. Fortunato and D. Hric.  Community detection in networks: A user guide.  _Physics reports_ , 659:1–44, 2016.
  * Gao et al. (2019) Z. Gao, G. Fu, C. Ouyang, S. Tsutsui, X. Liu, J. Yang, C. Gessner, B. Foote, D. Wild, Y. Ding, et al.  edge2vec: Representation learning using edge semantics for biomedical knowledge discovery.  _BMC bioinformatics_ , 20(1):1–15, 2019.
  * Goyal and Ferrara (2018) P. Goyal and E. Ferrara.  Graph embedding techniques, applications, and performance: A survey.  _Knowledge-Based Systems_ , 151:78–94, 2018.
  * Goyal et al. (2020) P. Goyal, S. R. Chhetri, and A. Canedo.  dyngraph2vec: Capturing network dynamics using dynamic graph representation learning.  _Knowledge-Based Systems_ , 187:104816, 2020.
  * Grover and Leskovec (2016) A. Grover and J. Leskovec.  node2vec: Scalable feature learning for networks.  In _Proceedings of the 22nd ACM SIGKDD international conference on Knowledge discovery and data mining_ , pages 855–864, 2016.
  * Gulikers et al. (2015) L. Gulikers, M. Lelarge, and L. Massoulié.  An impossibility result for reconstruction in a degree-corrected planted-partition model.  _arXiv preprint arXiv:1511.00546_ , 2015.
  * Gulikers et al. (2017) L. Gulikers, M. Lelarge, and L. Massoulié.  A spectral method for community detection in moderately sparse degree-corrected stochastic block models.  _Advances in Applied Probability_ , 49(3):686–721, 2017.
  * Jaffe et al. (2020) A. Jaffe, Y. Kluger, O. Lindenbaum, J. Patsenker, E. Peterfreund, and S. Steinerberger.  The spectral underpinning of word2vec.  _Frontiers in applied mathematics and statistics_ , 6:593406, 2020.
  * Jastrzebski et al. (2017) S. Jastrzebski, D. Leśniak, and W. M. Czarnecki.  How to evaluate word embeddings? on importance of data efficiency and simple supervised tasks.  _arXiv preprint arXiv:1702.02170_ , 2017.
  * Joulin et al. (2017) A. Joulin, E. Grave, P. Bojanowski, and T. Mikolov.  Bag of tricks for efficient text classification.  In _Proceedings of the 15th Conference of the European Chapter of the Association for Computational Linguistics: Volume 2, Short Papers_ , pages 427–431. Association for Computational Linguistics, April 2017.
  * Kamigaito and Hayashi (2021) H. Kamigaito and K. Hayashi.  Why does negative sampling not work well? analysis of convexity in negative sampling.  2021\.
  * Karrer and Newman (2011) B. Karrer and M. E. Newman.  Stochastic blockmodels and community structure in networks.  _Physical review E_ , 83(1):016107, 2011.
  * Kazemi et al. (2019) S. M. Kazemi, R. Goel, S. Eghbali, J. Ramanan, J. Sahota, S. Thakur, S. Wu, C. Smyth, P. Poupart, and M. Brubaker.  Time2vec: Learning a vector representation of time.  _arXiv preprint arXiv:1907.05321_ , 2019.
  * Khattak et al. (2019) F. K. Khattak, S. Jeblee, C. Pou-Prom, M. Abdalla, C. Meaney, and F. Rudzicz.  A survey of word embeddings for clinical text.  _Journal of Biomedical Informatics_ , 100:100057, 2019.
  * Kutuzov et al. (2017) A. Kutuzov, M. Fares, S. Oepen, and E. Velldal.  Word vectors, reuse, and replicability: Towards a community repository of large-text resources.  In _Proceedings of the 58th Conference on Simulation and Modelling_ , pages 271–276. Linköping University Electronic Press, 2017.
  * Landgraf and Bellay (2017) A. J. Landgraf and J. Bellay.  Word2vec skip-gram with negative sampling is a weighted logistic pca.  _arXiv preprint arXiv:1705.09755_ , 2017.
  * Le Magoarou and Gribonval (2016) L. Le Magoarou and R. Gribonval.  Flexible multilayer sparse approximations of matrices and applications.  _IEEE Journal of Selected Topics in Signal Processing_ , 10(4):688–700, 2016.
  * Leskovec and Krevl (2014) J. Leskovec and A. Krevl.  SNAP Datasets: Stanford large network dataset collection.  http://snap.stanford.edu/data, June 2014.
  * Lloyd (1982) S. Lloyd.  Least squares quantization in pcm.  _IEEE transactions on information theory_ , 28(2):129–137, 1982.
  * Louart and Couillet (2018) C. Louart and R. Couillet.  Concentration of measure and large random matrices with an application to sample covariance matrices.  _arXiv preprint arXiv:1805.08295_ , 2018.
  * Mikolov et al. (2013) T. Mikolov, I. Sutskever, K. Chen, G. S. Corrado, and J. Dean.  Distributed representations of words and phrases and their compositionality.  _Advances in neural information processing systems_ , 26, 2013.
  * Mimno and Thompson (2017) D. Mimno and L. Thompson.  The strange geometry of skip-gram with negative sampling.  In _Empirical Methods in Natural Language Processing_ , 2017.
  * Narayanan et al. (2017) A. Narayanan, M. Chandramohan, R. Venkatesan, L. Chen, Y. Liu, and S. Jaiswal.  graph2vec: Learning distributed representations of graphs.  _arXiv preprint arXiv:1707.05005_ , 2017.
  * Newman (2003) M. E. Newman.  The structure and function of complex networks.  _SIAM review_ , 45(2):167–256, 2003.
  * Newman (2006) M. E. Newman.  Modularity and community structure in networks.  _Proceedings of the national academy of sciences_ , 103(23):8577–8582, 2006.
  * Ng (2017) P. Ng.  dna2vec: Consistent vector representations of variable-length k-mers.  _arXiv preprint arXiv:1701.06279_ , 2017.
  * Nguyen et al. (2018) G. H. Nguyen, J. B. Lee, R. A. Rossi, N. K. Ahmed, E. Koh, and S. Kim.  Continuous-time dynamic network embeddings.  In _Companion proceedings of the the web conference 2018_ , pages 969–976, 2018.
  * Nickel and Kiela (2017) M. Nickel and D. Kiela.  Poincaré embeddings for learning hierarchical representations.  _Advances in neural information processing systems_ , 30, 2017.
  * Pennington et al. (2014) J. Pennington, R. Socher, and C. D. Manning.  Glove: Global vectors for word representation.  In _Proceedings of the 2014 conference on empirical methods in natural language processing (EMNLP)_ , pages 1532–1543, 2014.
  * Perozzi et al. (2014) B. Perozzi, R. Al-Rfou, and S. Skiena.  Deepwalk: Online learning of social representations.  In _Proceedings of the 20th ACM SIGKDD international conference on Knowledge discovery and data mining_ , pages 701–710, 2014.
  * Qin et al. (2016) P. Qin, W. Xu, and J. Guo.  A novel negative sampling based on tfidf for learning word representation.  _Neurocomputing_ , 177:257–265, 2016.
  * Rahman et al. (2018) M. Rahman, T. K. Saha, M. A. Hasan, K. S. Xu, and C. K. Reddy.  Dylink2vec: Effective feature representation for link prediction in dynamic networks.  _arXiv preprint arXiv:1804.05755_ , 2018.
  * Rehurek and Sojka (2011) R. Rehurek and P. Sojka.  Gensim–python framework for vector space modelling.  _NLP Centre, Faculty of Informatics, Masaryk University, Brno, Czech Republic_ , 3(2), 2011.
  * Rozemberczki et al. (2019) B. Rozemberczki, R. Davies, R. Sarkar, and C. Sutton.  Gemsec: Graph embedding with self clustering.  In _Proceedings of the 2019 IEEE/ACM international conference on advances in social networks analysis and mining_ , pages 65–72, 2019.
  * Sato et al. (2021) K. Sato, M. Oka, A. Barrat, and C. Cattuto.  Predicting partially observed processes on temporal networks by dynamics-aware node embeddings (dyane).  _EPJ Data Science_ , 10(1), 2021.  doi: 10.1140/epjds/s13688-021-00277-8.  Cited by: 2; All Open Access, Gold Open Access, Green Open Access.
  * Scott (1988) J. Scott.  Trend report social network analysis.  _Sociology_ , pages 109–127, 1988.
  * Selva Birunda and Kanniga Devi (2021) S. Selva Birunda and R. Kanniga Devi.  A review on word embedding techniques for text classification.  _Innovative Data Communication Technologies and Application: Proceedings of ICIDCA 2020_ , pages 267–281, 2021.
  * Shi and Malik (2000) J. Shi and J. Malik.  Normalized cuts and image segmentation.  _IEEE Transactions on pattern analysis and machine intelligence_ , 22(8):888–905, 2000.
  * Tiomoko et al. (2019) M. Tiomoko, R. Couillet, F. Bouchard, and G. Ginolhac.  Random matrix improved covariance estimation for a large class of metrics.  In _International Conference on Machine Learning_ , pages 6254–6263. PMLR, 2019.
  * Torricelli et al. (2020) M. Torricelli, M. Karsai, and L. Gauvin.  weg2vec: Event embedding for temporal networks.  _Scientific Reports_ , 10(1):1–11, 2020.
  * Von Luxburg (2007) U. Von Luxburg.  A tutorial on spectral clustering.  _Statistics and computing_ , 17:395–416, 2007.
  * Wang et al. (2018) Y. Wang, S. Liu, N. Afzal, M. Rastegar-Mojarad, L. Wang, F. Shen, P. Kingsbury, and H. Liu.  A comparison of word embeddings for the biomedical natural language processing.  _Journal of biomedical informatics_ , 87:12–20, 2018.
  * Wasserman and Faust (1994) S. Wasserman and K. Faust.  Social network analysis: Methods and applications.  _Cambridge university press_ , 1994.
  * Zipf (2013) G. K. Zipf.  _The psycho-biology of language: An introduction to dynamic philology_.  Routledge, 2013.

##  A Appendix

In this section we provide the proofs of the statements appearing in the main
article. Proposition 1 states the sufficient condition to guarantee that
$Z_{i}/n$ is concentrated around its expectation value.

######  Proposition (1).

Let $X\in\mathbb{R}^{n\times d}$ be a random matrix so that
$\|\bm{x}_{i}\|=f_{i}$ where $f_{i}\in\mathbb{R}^{+}$ are random i.i.d.
variables distributed $f_{i}\sim p_{f}$. Suppose that the expected value of
$f_{i}=O_{n}(1)$ and that for all $\omega$ there is a $t^{*}>0$ so that for
$t>t^{*}$, $p_{f}(x>t)<2e^{-\omega t}$. Then, under these hypotheses, for all
large $n$, $Z_{i}/n$ is concentrated around its expected value.

######  Proof.

Let $y_{k}^{(i)}={\rm exp}\\{\|\bm{x}_{i}\|\|\bm{x}_{k}\|\alpha_{ik}\\}$,
where $\alpha_{ik}$ is the cosine-similarity between $\bm{x}_{i}$ and
$\bm{x}_{k}$. With this notation, $Z_{i}$ can be written as

$\displaystyle Z_{i}=\sum_{k\in\mathcal{V}}y_{k}^{(i)}.$

Let us now compute the expectation value of $y_{k}^{(i)}$. We assume
$\bm{x}_{k}$ to be a random variable and $\bm{x}_{i}$ to be fixed and that the
$\alpha_{ik}$ is independent from $f_{k}$. We then write:

$\displaystyle\mathbb{E}[y_{k}^{(i)}]$
$\displaystyle=\int_{-1}^{1}d\alpha\leavevmode\nobreak\
p_{X}(\alpha)\underbrace{\int_{0}^{\infty}d\theta\leavevmode\nobreak\
p_{f}(\theta)\leavevmode\nobreak\
e^{\|\bm{x}_{i}\|\theta\alpha}}_{h(\alpha)<\infty}=\int_{-1}^{1}d\alpha\leavevmode\nobreak\
p_{X}(\alpha)\leavevmode\nobreak\ h(\alpha)<\infty,$

where we exploited the fact that $p_{f}$ goes to zero faster than an
exponential distribution. With a similar calculation one sees that, under the
same assumption, $\mathbb{E}[(y_{k}^{(i)})^{2}]$ is also finite. Chebyshev
inequality (Feller, 1968) can then be applied to show that $Z_{i}/n$ is
concentrated around its expected value, in fact

$\displaystyle\mathbb{P}\left(\left|\frac{Z_{i}}{n}-\frac{\mathbb{E}[Z_{i}]}{n}\right|>t\right)\leq\frac{\mathbb{V}[Z_{i}]}{(nt)^{2}},$

where $\mathbb{V}$ denotes the variance. The condition that the mean of
$\mathbb{E}[f_{i}]=O_{n}(1)$ ensures that $\mathbb{V}[Z_{i}]\propto n$ and
hence, that the right hand-side is smaller than $1$ for all large $n$.

Our main result is Proposition 2 that we first recall here.

######  Proposition (2).

Let $\\{\bm{x}_{i}\\}_{i\in\mathcal{V}}$ be a set of random vectors in
$\mathbb{R}^{d}$ and let $\ell:\mathcal{V}\to[\kappa]$ be a labeling function.
Suppose that for all $j\in\mathcal{V}$ so that $\ell(j)=a$ the corresponding
embedding vectors are distributed as $\bm{x}_{i}\overset{\rm
i.i.d}{\sim}P_{a}(\bm{\mu}_{a},\Omega_{a})$, where
$\bm{\mu}_{a}\in\mathbb{R}^{d}$ and $\Omega_{a}\in\mathbb{R}^{d\times d}$ are
the mean and covariance matrix of $P_{a}$. Then, for $i$ fixed and
$\ell(j)=a$,

$\displaystyle\bm{x}_{i}^{T}\bm{x}_{j}\sim\mathcal{N}(\bm{x}_{i}^{T}\bm{\mu}_{a},\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}),$

where $\mathcal{N}$ is a Gaussian random variable. Letting $\pi_{a}$ be the
fraction of items in class $a$, we get

$\displaystyle\frac{Z_{i}}{n}\approx\sum_{a=1}^{\kappa}\pi_{a}e^{\bm{x}_{i}^{T}\bm{\mu}_{a}+\frac{1}{2},\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}}.$

We now proceed with a theoretical justification of our main result.

######  Proof.

The first step to prove Proposition 2 requires to show that
$\bm{x}_{i}^{T}\bm{x}_{j}$ is distributed according to a Gaussian random
variable. The intuition is that the scalar product is a sum of random
variables with a finite mean and variance. For simplicity, we assume all
entries of $\bm{x}_{j}$ to be independent, each one a random variable
$R(m_{i},v_{i})$ with mean $m_{a}$ and variance $v_{a}$. Under this
assumption, we may write

$\displaystyle\bm{x}_{i}^{T}\bm{x}_{j}=\sum_{a=1}^{d}x_{ia}R(m_{a},v_{a})=\sum_{a=1}^{d}R({x_{ia}m_{a},x_{ia}^{2}v_{a}}).$

This is a sum of independent random variables. To show the convergence towards
a Gaussian distribution we use Lyapunov central limit theorem. To apply it, we
must first check that the Lyapunov condition holds, as a consequence of the
hypotheses of Proposition 1:

$\displaystyle\sum_{j\in\mathcal{V}}|\bm{x}_{i}^{T}\bm{x}_{j}-\mathbb{E}[\bm{x}_{i}^{T}\bm{x}_{j}]|^{3}\leq
n\leavevmode\nobreak\ \underset{j\in\mathcal{V}}{\rm max}\leavevmode\nobreak\
2f_{j}^{6}\sim 2n[{\rm log}(n)]^{6}=o_{n}(n^{3/2}),$

where we exploited the fact the maximum of an exponential distribution grows
as ${\rm log}(n)$. From this relation the Lyapunov central limit theorem can
be applied to obtain the result of the claim. To finally obtain the estimated
value of $Z_{i}$, we solve the Gaussian integral:

$\displaystyle\mathbb{E}\left[\frac{Z_{i}}{n}\right]$
$\displaystyle=\frac{1}{n}\sum_{j\in\mathcal{V}}\mathbb{E}\left[e^{\bm{x}_{i}^{T}\bm{x}_{j}}\right]$
$\displaystyle=\frac{1}{n}\sum_{a=1}^{k}\sum_{j\in\mathcal{V}_{a}}\mathbb{E}\left[e^{\bm{x}_{i}^{T}\bm{x}_{j}}\right]$
$\displaystyle=\frac{1}{n}\sum_{a=1}^{\kappa}\sum_{j\in\mathcal{V}_{a}}\int_{-\infty}^{\infty}dt\leavevmode\nobreak\
\frac{1}{\sqrt{2\pi\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}}}\cdot{\rm
exp}\left\\{-\frac{(\bm{x}_{i}^{T}\bm{\mu}_{a}-t)^{2}}{2\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}}+t\right\\}$
$\displaystyle=\frac{1}{n}\sum_{a=1}^{\kappa}\sum_{j\in\mathcal{V}_{a}}e^{\bm{x}_{i}^{T}\bm{\mu}_{a}+\frac{1}{2}\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}}$
$\displaystyle=\sum_{a=1}^{\kappa}\pi_{a}e^{\bm{x}_{i}^{T}\bm{\mu}_{a}+\frac{1}{2}\bm{x}_{i}^{T}\Omega_{a}\bm{x}_{i}}.$

Given Proposition 1, $Z_{i}/n$ is concentrated around its expected value.