the underlying interaction mechanisms to remain the same over this period.
Slicing this data into time intervals would provide deeper insights into the
interaction mechanisms at stake and their evolution over time.

###  .7 Dynamic interactions

####  .7.1 Introduction

Dynamic networks are powerful tools to visualize and model interactions
between different entities that can evolve over time. The network’s nodes
represent the interacting entities, and ties between these nodes represent an
interaction. In many real-world situations, the strength of the ties can vary
over time –on music streaming websites for instance, users’ affinity with
various musical genres can vary greatly over time [73, 122]. Such network is
said dynamic.

Now, every interaction does not have the same importance. A music listener can
like both Rock and Jazz, but might prefer one over the other. This person’s
tie to each musical genre does not have the same intensity; each tie is
associated with a number, representing the strength of the interaction. The
network is said to be dynamic and valued.

However, in many real-world situations, valued networks are not enough to
fully represent a given situation. The same music listener as before can have
different opinions on musical genres; she can like it, dislike it, be bored of
it, like to listen to it only in the morning, or at night, etc. Each of these
relations can be represented by its own tie in the network, each having its
own value. The network is said to be dynamic, valued and labelled.

######  Inferring dynamic, valued and labelled networks

Networks are high-dimensional objects, whose exact inference is a difficult
problem –as stated in the previous section. Several ways to achieve this task
have been proposed, the Stochastic Block Models (SBM) family being one of the
most popular approaches, see Chapter - and [60, 50, 27]. The underlying
assumption is that certain sets of nodes behave similarly. They can be
described using a single mathematical object instead of individual ones: the
so-called clusters. Instead of modelling every edge for every node in a
network, only the edges between these sets are modelled, which makes the task
much more tractable. Each cluster is then associated with a labelled edge, and
each node is associated with a cluster. A variant of SBM that allows more
expressive power is called the Mixed-Membership SBM (MMSBM), where each node
can belong to several clusters in different proportions at the same time [144,
4, 42, 120]. A major advantage of this model’s family is that it yields
readily interpretable results (interpretable clusters), unlike most existing
neural-network-based modelling of labelled networks [33].

######  Overview of the proposed approach

The goal of this section is therefore to provide a way to infer networks that
are dynamic, valued and labelled by assuming a block structure – by using a
mixed-membership SBM. After a careful review of dynamic network inference
literature, it seems that no prior work tackles such a task. Although some
previous works handle similar problems, their adaptation to the case at hand
is not trivial. Besides, the solution we propose here is conceptually much
simpler than those present in the literature. Last but not least, our method
is readily pluggable into most existing MMSBM for labelled and valued networks
(such as SIMSBM and its special cases, comprising the IMMSBM and each of [4,
42, 120]) as their temporal extension.

We develop an EM optimization algorithm that scales linearly with the size of
the dataset, demonstrate the effectiveness of our approach on both synthetic
and real-world datasets, and further detail a possible application scenario.

Figure 13: Users’ attachment to groups can vary over time — A music listener
could cyclically prefer Rock, Jazz, or Pop music (top), or listen to either of
these without any specific pattern (bottom). For 200 epochs each containing
only 5 observations, our approach (in red) infers any smooth dynamic
membership pattern and does it more accurately than static models (in green
[42]) and models that consider each time slice independently (in blue [120]).

####  .7.2 State of the art and limitations

#####  .7.2.a Notations

We consider a network of $I$ nodes and $O$ labels. All the clustering models
discussed in this section can be represented using a matrix
$\theta\in\mathbb{R}^{I\times K}$ accounting for memberships over a set of $K$
possible clusters for each of $I$ nodes, and a block-interaction matrix
$p\in\mathbb{R}^{K\times K\times O}$ linking each of $K$ clusters to every of
$O$ labels. Both $\theta$ and $p$ can vary over time. The network is said to
be unlabelled when $O=1$, and binary (as opposed to valued) if an edge can
only exist or not exist.

#####  .7.2.b Dynamic unlabelled networks - Single-membership

Single-membership SBM considers a membership matrix such as
$\theta\in\\{0;1\\}^{I\times K}$: each membership vector equals $1$ for one
cluster, and $0$ everywhere else. Literature also speaks of “hard” clustering.
The authors in [136, 135] proposed to model a binary unlabelled dynamic
network using a label-switching inference method along with a Sequential Monte
Carlo algorithm [66]. Both the membership and the interaction matrices can
vary over time, thus supposing two independent underlying Markov processes
[66].

In [138, 119], the authors propose to model a binary dynamic and unlabelled
network. The cluster interaction matrix $p$ can vary over time while keeping
the memberships $\theta$ fixed over time. The entries of $p$ are drawn from a
Dirichlet distribution and expressed as a Chinese Restaurant Process. This
process converges to a Dirichlet distribution and allows to infer a
potentially infinite number of clusters. This model is therefore non-
parametric and inferred using an MCMC algorithm.

A novel way of tackling the problem has been proposed in [82, 81], where the
authors propose to model the cluster interaction and membership matrices as
Poisson processes, that explicitly model the temporal dependency without
slicing the time dimension into episodes. The method allows to infer varying
membership and interaction matrices for dynamic binary or valued networks, but
their results have shown that allowing both to vary simultaneously leads to
identifiability and label switching issues [38]. This conclusion seems
reasonable, given none of these SBM algorithms can reach a global optimum in
the likelihood function. During optimization, a model with both membership
$\theta$ and interaction $p$ matrices is all the more likely to get stuck into
a local maximum if both can vary over time.

Finally, we can mention the existence of SBM variants that account for dynamic
degree-correction [131] or that enforce a scale-free characteristic [132].
However, all these methods consider unlabelled networks and consider a hard
clustering which does not allow for as much expressive power as the Mixed-
Membership approaches.

#####  .7.2.c Dynamic unlabelled networks - Mixed-membership

Mixed-membership SBM considers a membership matrix such as
$\theta\in\mathbb{R}^{I\times K}$, where each membership vector is normalized
to $1$. Literature also refers to it as “soft” clustering.

Similar to [138, 119], a method for inferring dynamical binary unlabelled
networks has been proposed in [34]. The membership vector of each piece of
information is drawn from a Chinese Restaurant Process (CRP) according to the
number of times a node has already been associated with each cluster before.
The resulting process yields a distribution over an infinity of available
clusters. The formulation as a CRP arises naturally because the prior on
membership vectors is typically a Dirichlet distribution; CRP naturally
converges to this distribution. The block-interaction matrix $p$ does not
evolve with time. The article shows a complexity analysis that suggests the
methods run with a complexity of $\mathcal{O}(N^{2})$ which makes it unfit for
large-scale real-world applications.

The work the most closely related to ours is [133]. This seminal work proposed
the dMMSB as a way to model dynamic binary unlabelled networks using a
variational algorithm [76]. To do so, the authors modify the original MMSBM
[4] to consider a logistic normal distribution as a prior on the membership
vectors $\vec{\theta_{i}}$. This choice allows to model correlations between
membership vectors’ evolution [3]. The membership vectors are then embedded in
a state-space model, which is a space where we can define a linear transition
between two time points for a given variable. The authors define such a
trajectory for the membership vectors as a linear function of the previous
time point. The trajectory is estimated and smoothed using a Kalman Filter.
This approach is the most closely related to ours, as it proposes to consider
the temporal dependency directly in a prior distribution over memberships,
noted $P(\theta)$.

However, this model is not fit for the task at hand. It considers unlabelled
and binary networks [76], and extension to labelled and valued networks is not
trivial. The proposed optimization algorithm requires a loop until convergence
at each EM iteration, making it unable to handle large datasets. It is not
designed to let the clusters interaction matrix $p$ depend on time, which we
alleviate here. And most importantly, it assumes a linear transition between
time steps in the state space, while we do not assume any kernel function in
our proposed approach. [133] has been extended to consider $P(\theta)$ as a
logistic normal mixture prior [58], which improves its expressive power.
However, it does not address the aforementioned points.

#####  .7.2.d Static labelled networks - Mixed-membership

Recent years saw a rise of Bayesian methods for inferring static valued
labelled networks using MMSBM variants [42, 120, 102]. Note that in [120], the
authors consider a temporal slicing of the data and consider each slice as
independent from the others; a time slice is considered as a node in a
tripartite network. We will compare our approach to this modelling later. We
do not present the details of these works here, as we will develop their in-
depth functioning in Section .7.3.

The method we propose here uses these works as a base. This section focuses on
making the prior probability of both $\theta$ and $p$ time-dependent to model
these parameters’ dynamics. We provide a ready-to-use temporal plug-in for
each of the works we have just presented in this section and in Section .6. It
applies to dynamical, valued, and labelled networks in a mixed-membership
context, and inference is conducted with a scalable variational EM algorithm.

####  .7.3 SDSBM – Simple Dynamic labelled MMSBM

#####  .7.3.a Base model

In this section, we present the simplest form of a labelled MMSBM, or
SIMSBM(1), for demonstration purposes. We illustrate the structure of this
model in Fig. 14. Its derivation trivially extends to [42, 120, 102], the
IMMSBM, and any of the SIMSBM iterations discussed Section .6.2 (we detail our
work’s inclusion in these more complex models in Appendix, Section F).

Figure 14: Illustration of the SIMSBM(1), which is the base model coupled to
the SDSBM prior. Nodes are associated with clusters, and clusters are
associated with labels. Ties between each entity represent a membership and
can take values between 0 and 1.

We consider a set of $I$ nodes that can be associated with $O$ possible labels
on a discrete time interval, or epoch, written $t$. We assume that nodes can
be efficiently represented as a mixture of $K$ available clusters at each time
step, each of which is in turn linked to the labels. The membership of each of
$I$ nodes into each of the $K$ possible clusters at time $t$ is encoded in the
membership matrix $\theta^{(t)}\in\mathbb{R}^{I\times K}$. One vector
$\theta_{i}^{(t)}$ represents the probability that $i$ belongs to any of the
$K$ clusters at time $t$, and is normalized as:

$\sum_{k\in K}\theta_{i,k}^{(t)}=1\ \forall i,t$ (23)

The probability of each cluster to be associated with each label at time $t$
is encoded in the matrix $p^{(t)}\in\mathbb{R}^{K\times O}$. An entry
$p_{k}^{(t)}(o)$ represents the probability that cluster $k$ is associated
with label $o$ at time $t$, and thus is normalized as:

$\sum_{o\in O}p_{k}^{(t)}(o)=1\ \forall k,^{(t)}$ (24)

Finally, the probability that a node $i$ is associated with label $o$ at time
$t$ (this can be seen as the probability an edge between $i$ and $o$ exists at
time $t$) is written:

$P^{(t)}(i\rightarrow o)=\sum_{k\in K}\theta_{i,k}^{(t)}p_{k}^{(t)}(o)$ (25)

Given a set $R^{\circ}$ of observed triplets $(i,o,t)$, the model’s posterior
distribution can be written [42, 120]:

$\displaystyle P(\theta,p|R^{\circ})$ $\displaystyle\propto
P(R^{\circ}|\theta,p)\prod_{t}P(\theta^{(t)})P(p^{(t)})$ (26)
$\displaystyle=\prod_{(i,o,t)\in R^{\circ}}$ $\displaystyle\sum_{k\in
K}\theta_{i,k}^{(t)}p_{k}^{(t)}(o)\prod_{t}\left(\prod_{i}P(\theta_{i}^{(t)})\prod_{k}P(p_{k}^{(t)})\right)$

Now, before we describe the optimization procedure, we must choose the priors
$P(\theta^{(t)})$ and $P(p^{(t)})$.

#####  .7.3.b Simple Dynamic prior

We formulate the prior distribution over $\theta^{(t)}$ and $p^{(t)}$
following a simple assumption: the parameters at a given time are unlikely to
vary abruptly at small time scales — the apple does not fall far from the
tree. It means an entry $\theta^{(t_{1})}_{ik}$ should not differ so much from
$\theta^{(t_{2})}_{ik}$ for every $t_{2}$ close enough to $t_{1}$. Such
entries close to a reference time are called temporal neighbours.

Our a priori knowledge on each entry $\theta_{i}^{(t)}$ and $p_{k}^{(t)}$ is
that they should not differ significantly from their temporal neighbours. This
is a fundamental difference with [133], where the next parameters values are
estimated using a Kalman Filter that only considers the previous time step.
Moreover, the authors assume a linear transition function, while we do not
make such a hypothesis. An illustration of the proposed approach is given in
Fig. 15, where the prior probability of a membership vector depends on its
temporal neighbours (in white).

######  Dirichlet distribution

Since each entry $\theta_{i}^{(t)}$ and $p_{k}^{(t)}$ is normalized to 1, we
consider a Dirichlet distribution as a prior, which naturally yields
normalized vectors such that $\sum_{n}x_{n}=1$. It reads:

$Dir(x|\alpha)=\frac{1}{B(\alpha)}\prod_{n}x_{n}^{\alpha_{n}-1}$ (27)

where $B(\cdot)$ is the multivariate beta function. In Eq. 27, the vector
$\alpha$ is called the concentration parameter and must be provided to the
model. Importantly to the model introduced in this section, it allows to
control the mode and the variance of the Dirichlet distribution.

We consider a concentration parameter such as $\alpha=1+\beta\alpha_{0}$, so
that for $\beta=0$ we recover a uniform prior over the simplex, and $\alpha>1$
so that the prior has a unique mode. The most frequent value drawn from Eq. 27
(or mode) is
$\mathbb{M}(x_{n})=\frac{\alpha_{n}-1}{\sum_{n}^{\prime}(\alpha_{n}^{\prime}-1)}=\frac{\alpha_{0,n}}{\sum_{n}^{\prime}\alpha_{0,n}}$.
We recover a uniform prior for $\beta=0$; the variance vanishes with $\beta\gg
1$ as $\frac{1}{\beta}$. The effect of various values of $\beta$ on the prior
distribution is illustrated in Fig. 15.

Figure 15: Prior probability on a membership vector for various values of
$\beta$ according to temporal neighbourhood — Darker means higher probability.
Projected on a simplex tri-space (each of 3 axes ranges from 0 to 1). The
white dots represent the temporal neighbours of the considered 3D vector.
Their average is given as $\langle\theta\rangle$ using a uniform weight
function $\kappa(t,t^{\prime})$ for illustration purpose. $\beta$ controls the
variable’s prior variance around its neighbours.

######  Prior’s mode

Our main assumption states that $\theta_{i}^{(t)}$ and $p_{k}^{(t)}$ do not
vary abruptly. To enforce this, we define their prior probability mode with
respect to their close temporal neighbours. The hyper-parameter $\beta$
controls the variance of the prior –that is, how much it should impact the
inference procedure. We express the Simple Dynamic prior parameters for
$\theta_{i}^{(t)}$ as:

$\alpha_{i,k}^{(t,\theta)}=1+\beta\underbrace{\left(\frac{\sum_{t^{\prime}\neq
t}\kappa(t,t^{\prime})\theta_{i,k}^{(t^{\prime})}}{\sum_{t^{\prime}\neq
t}\kappa(t,t^{\prime})}\right)}_{\langle\theta_{i,k}^{(t)}\rangle}$ (28)

where $\kappa(t,t^{\prime})$ is a weight function, and $\alpha^{(t,\theta)}$
corresponds to the concentration parameter for $\theta$ at time $t$. In
following experiments, we define the weight function as
$\kappa(t,t^{\prime})=\frac{N_{t^{\prime}}}{|t-t^{\prime}|}$, where
$N_{t^{\prime}}$ is the number of observations made at time $t^{\prime}$, so
that temporal neighbours’ influence decrease as the inverse of temporal
distance. We illustrate the influence of this particular kernel function on
the prior probability on membership at all times in Fig. 16. In particular we
see that with this expression, the prior probability variance collapses to 0
when the considered time is very close to a temporal neighbour.

Figure 16: Prior probability’s variance on memberships at all times according
to the temporal neighbourhood — Variance of the prior over a membership entry
(filled curves, we represented 3 such entries as illustration) as a function
of time, given some temporal neighbours (black dots). This illustration
considers an averaging kernel as
$\kappa(t,t^{\prime})=\frac{1}{|t-t^{\prime}|}$. When inferring a parameter
$x^{(t)}$ at a time $t$, the variance of its prior probability $P(x^{(t)})$
depends on $t$ relative to the temporal neighbours. Here for instance, the
variance is null at $t=2$ because $\kappa(t,t^{\prime})$ diverges, and so does
$\alpha^{(t)}$, hence the variance collapsing to 0.

The mode of the prior over a variable is then the average value of its the
temporal neighbours weighted by $\kappa(t,t^{\prime})$, noted
$\langle\theta_{i,k}^{(t)}\rangle$. Note that this holds because
${\sum_{k}\langle\theta_{i,k}^{(t)}\rangle=1\ \forall i,t}$. Besides, the
prior variance is a decreasing function of $\beta$; when $\beta=0$ the prior
is uniform over the simplex, and when $\beta\rightarrow\infty$ the variance
goes to $0$, as illustrated Fig. 15. The same reasoning holds for
$p_{k}^{(t)}$, with prior parameters $\alpha_{k,o}^{(t,p)}=1+\beta\langle
p_{k}^{(t)}(o)\rangle$.

######  Priors expression

Finally, we give the final log-priors on $\theta_{i}^{(t)}$ and $p_{k}^{(t)}$:

$\displaystyle
P(\theta_{i}^{(t)}|\\{\theta_{i,k}^{(t^{\prime})}\\}_{t^{\prime}\neq t})$
$\displaystyle\propto\prod_{k}{\theta_{i,k}^{(t)}}^{\beta\langle\theta_{i,k}^{(t)}\rangle}$
(29) $\displaystyle
P(p_{k}^{(t)}(o)|\\{p_{k}^{(t^{\prime})}(o)\\}_{t^{\prime}\neq t})$
$\displaystyle\propto\prod_{o}{p_{k}^{(t)}(o)}^{\beta\langle
p_{k}^{(t)}(o)\rangle}$

We omitted the normalisation factor for clarity –it does not influence the
inference procedure. Note that the dependence of each entry of
$\theta_{i}^{(t)}$ and $p_{k}^{(t)}$ raises additional questions from a
statistical perspective. It is unsure whether the specified conditional
distributions would allow to retrieve a global, joint probability distribution
for all entries of those vectors. Therefore, our approach works as a
penalization of the likelihood function, but does not imply the existence of a
statistical law on $\theta_{i}^{(t)}$ and $p_{k}^{(t)}$.

#####  .7.3.c Inference

######  E step

We develop an EM inference procedure for maximizing the log-posterior
distribution defined in Eq. 26. The expectation step computes the expected
probability of a latent variable (here a cluster $k$) being chosen given each
entry of $R^{\circ}$. Since such latent variables do not appear in the priors
expressions, the expectation step remains unchanged by the introduction of the
Simple Dynamic Priors; in general, prior distributions do not intervene in the
computation of the expectation step [16]. The E step for such labelled
networks has already been derived on many occasions. Therefore, we give the
expectation step equation without explicit derivation (full derivation is
however provided in Appendix, Section D, and in Section .6.2.b.1 and Section
.6.3.b.1).

The expectation of the latent variable $k$ given an observation $(i,o,t)\in
R^{\circ}$, written $\omega_{i,o}^{(t)}(k)$, is defined as:

$\omega_{i,o}^{(t)}(k)=\frac{\theta_{i,k}^{(t)}p_{k}^{(t)}(o)}{\sum_{k^{\prime}}\theta_{i,k^{\prime}}^{(t)}p_{k^{\prime}}^{(t)}(o)}$
(30)

Using this expression, we can rewrite the log-likelihood $\log
P(R^{\circ}|\theta,p)$ as [16, 42, 120]:

$\log P(R^{\circ}|\theta,p)=\sum_{(i,o,t)\in R^{\circ}}\sum_{k\in
K}\omega_{i,o}^{(t)}(k)\log\frac{\theta_{i,k}^{(t)}p_{k}^{(t)}(o)}{\omega_{i,o}^{(t)}(k)}$
(31)

######  M step

Taking back the first line of Eq. 26 and substituting with Eq. 29 and Eq. 31,
we get an unconstrained expression of the posterior distribution. We introduce
Lagrange multipliers to account for the constraints of Eq. 23
($\phi_{i}^{(t)}$) and Eq. 24 ($\psi_{i}^{(t)}$), and finally compute the
maximization equations with respect to the model’s parameters. Starting with
the membership matrix entries $\theta_{i,k}^{(t)}$:

$\displaystyle\frac{\partial\left(\log
P(\theta,p|R^{\circ})-\sum_{i^{\prime},t^{\prime}}\phi_{i^{\prime}}^{(t^{\prime})}(\sum_{k^{\prime}}\theta_{i^{\prime},k^{\prime}}^{(t^{\prime})}-1)\right)}{\partial\theta_{i,k}^{(t)}}=0$
$\displaystyle\Leftrightarrow\sum_{o\in\partial(i,t)}\frac{\omega_{i,o}^{(t)}(k)}{\theta_{i,k}^{(t)}}+\frac{\beta\langle\theta_{i,k}^{(t)}\rangle}{\theta_{i,k}^{(t)}}-\phi_{i}^{(t)}=0$
$\displaystyle\Leftrightarrow\sum_{o\in\partial(i,t)}\omega_{i,o}^{(t)}(k)+\beta\langle\theta_{i,k}^{(t)}\rangle=\phi_{i}^{(t)}\theta_{i,k}^{(t)}$
$\displaystyle\Leftrightarrow\sum_{o\in\partial(i,t)}\underbrace{\sum_{k}\omega_{i,o}^{(t)}(k)}_{=1\text{
(Eq.\leavevmode\nobreak\ \ref{eq-
omega})}}+\beta\underbrace{\sum_{k}\langle\theta_{i,k}^{(t)}\rangle}_{=1\text{
(Eq.\leavevmode\nobreak\ \ref{eq-normtheta})}}=\phi_{i}^{(t)}$
$\displaystyle\Leftrightarrow\frac{\sum_{o\in\partial(i,t)}\omega_{i,o}^{(t)}(k)+\beta\langle\theta_{i,k}^{(t)}\rangle}{N_{i,t}+\beta}=\theta_{i,k}^{(t)}$
(32)

where $\partial(i,t)=\\{o|(i,\cdot,t)\in R^{\circ}\\}$ is the subset of labels
associated with both $i$ and $t$, and $N_{i,t}=|\partial(i,t)|$ is the size of
this set. Note that for $\beta=0$ we recover the M-step of standard static
MMSBM models, see [42, 120] and Section .6.

The derivation of the M-step for the entries $p_{k}^{(t)}(o)$ is identical and
yields (see Appendix, Section E, for details):

$\displaystyle p_{k}^{(t)}(o)=\frac{\sum_{(i,t)\in\partial
o}\omega_{i,o}^{(t)}(k)+\beta\langle p_{k}^{(t)}(o)\rangle}{\sum_{(i,o,t)\in
R^{\circ}}\omega_{i,o}^{(t)}(k)+\beta}$ (33)

#####  .7.3.d Discussion

######  Easy to use

We briefly review some key points of the Simple Dynamic prior. Its
introduction induces minor changes in the existing works on MMSBM for labelled
networks. Compared to [42, 120] and Section .6, its introduction boils down to
simply adding a term $\beta\langle x\rangle$ to the numerator of maximization
equations, and the corresponding normalizing term $\beta$ to the denominator.
This way, our approach is ready-to-use to make these models, and variants
built on them, dynamic –explicit derivations are provided in Appendix, Section
F.

######  Flexible dynamic modelling

The prior allows us to consider that some parameters are dynamic and that
others are not. For instance, when several membership matrices are involved,
as in [42, 120]), setting $\beta=0$ for some makes them time-invariant (or
universal). The SD prior also allows choosing whether the block-interaction
tensor $p$ is dynamic.

In general, $\beta$ does not have to be identical for every membership matrix,
or even every entry $i$ of each of them. Moreover, it is not mandatory for
$\beta$ to be constant over time. A dynamic parameter $\beta(t)$ is especially
relevant when epochs are not evenly spaced over time; $\beta(t)$ would
typically lower the variance (by increasing) when temporal neighbours are
closer (right plot in Fig. 15).

To summarize, $\beta$ allows controlling the time scale over which variables
may vary. This allows greater modelling flexibility, allowing to jointly model
universal variables ($\beta=0$) and dynamical ones ($\beta\neq 0$).

######  Tuneable temporal dependence

Finally, the choice of the averaging kernel function $\kappa(t,t^{\prime})$ is
important. It allows choosing the range over which the inference of a variable
should rely on its temporal neighbours. A formulation as the inverse of time
difference seems relevant: the weight of a neighbour appearing at a time
$\delta t$ later should diverge as $\delta t\rightarrow 0$, so that continuity
is ensured. Besides, it allows controlling the smoothness of the curve with
respect to time by tuning the weight function as
$\kappa(t,t^{\prime})=\frac{N_{t^{\prime}}}{|t-t^{\prime}|^{a}}$ where
$a=1,2,...$ for instance, where $N_{t^{\prime}}$ is the number of observations
in the time slice $t^{\prime}$.

Overall, the Simple Dynamic prior works by inferring the variables using both
microscopic and mesoscopic temporal scales. If a time slice $t$ has very few
observations but some of its neighbours have a greater number of them for
instance; learning the parameters at $t$ is helped mostly by the population of
its closest ($\frac{1}{|\Delta t|}$) and most populated ($N_{t^{\prime}}$)
epochs, and less influenced by further and less populated epochs. This is what
is illustrated in Fig. 16.

#####  .7.3.e Experiments

######  Synthetic data

In this section, we address several situations in which our method
(abbreviated SDSBM for Simple Dynamic MMSBM) could be useful. Experiments are
run for $I=100$, $K=3$ and $O=3$, which are standard testing parameters in the
literature of dynamic networks inference [34, 82]. We choose to infer a
dynamic membership matrix $\theta^{(t)}$ and to provide a universal block-
interaction matrix $p\ \forall t$. Note that the model yields good
performances when $p$ also has to be inferred, but due to identifiability and
label-switching issues raised in [82], there is no unbiased way to assess the
correctness of the inferred memberships values. An experiment where $p$ is
inferred jointly to $\theta$ is provided and discussed in Appendix, Section G.
The expression we use for this matrix $p$ is given Eq. 34. We systematically
test two variation patterns for $\theta$: a sinusoidal pattern (Fig. 13-top)
and a broken-line pattern (Fig. 13-bottom). Each pattern is generated with a
different coefficient for each item; the memberships still sum to 1 at all
times.

$p=\begin{bmatrix}1-s&s&0\\\ 0&1-s&s\\\ s&0&1-s\end{bmatrix}$ (34)

To the best of our knowledge, the only attempt to model dynamic parameters in
labelled valued and dynamic networks using a MMSBM is [120]. In this work,
each epoch is modelled independently from the others. We refer to this
baseline as the “No coupling” or “NC” baseline. For reference, we also compare
to a baseline that does not consider the temporal dimension and infers a
single universal value for each variable (classical MMSBM) [42] and the models
Section .6.

We systematically perform a 5-fold cross-validation. The model is trained on
80% of the data, $\beta$ is tuned using 10% as a validation set, and the model
is evaluated on the 10% left. We choose as metrics the AUCROC and RMSE on the
real values of $\theta$ (black line in Fig. 13). The procedure is repeated 5
times; the error bars reported in the experimental results represent the
standard error over these folds.

Figure 17: Results on synthetic data — (a.) SDSBM retrieves the correct
dynamic memberships and is little influenced by the data slicing. (b.) SDSBM
works well on tiny datasets. (c.) SDSBM retrieves correct dynamic memberships
in challenging situations.

######  SDSBM unveils complex temporal patterns

In Fig. 17a., we consider 1,000 observations for each item $i\in I$ and vary
the number of epochs from 10 to 1,000. In the expression of $p$, $s$ is set to
0.05. For both the sinusoidal and line-broken memberships, the model shows
better predictive performances (in terms of AUCROC) than the proposed
baselines. Interestingly, the SDSBM performances remain stable as the number
of epochs increases unlike the NC baseline, which means it alleviates a bias
of the temporal modelling proposed in [120]. The RMSE with respect to the true
parameters remains low over the whole range of the tested number of epochs.
The RMSE increases as the number of epochs grows because the number of
parameters to estimate increases with it; this makes the inference more
subject to local variations, which in turn mechanically increases the RMSE.
Overall, SDSBM recovers dynamic variations of the membership vectors with
superior performance; a sample of the inferred dynamic memberships is shown in
Fig. 13.

######  SDSBM works with little data

A major problem that arises when considering temporal data is the scarcity of
observations, because slicing implies reducing the number of observations in
each slice. This concern mostly arises in social sciences, where data
retrieval cannot be automated and requires long-lasting human labour. Here, we
demonstrate that our method works in challenging conditions when data is
scarce. In Fig. 17b., we vary the number of observations available for each
item from 100 to 10.000, distributed over 100 epochs. Thus, in the most
challenging situation, there is only one observation per epoch used to
determine $I$ dynamic memberships over 3 clusters. In the expression of $p$,
$s$ is set to 0.05. We see Fig. 17b. that for both patterns, the predictive
power of SDSBM remains high in such conditions. Moreover, the RMSE on the true
dynamic memberships in this case is fairly low and decreases rapidly as the
number of observations increases. When the number of observations is high, the
“no coupling” baseline [120] reaches the performances of SDSBM. This is
because as the number of observations in each slice goes to infinity, the
models need to rely less on temporal neighbours. However, even for 10.000
observations per item (100 observations per epoch), SDSBM yields better
results than the proposed baselines. As an illustration, the results presented
in Fig. 13 have been inferred using only 5 observations per epoch.

######  SDSBM handles highly stochastic interaction patterns

Finally, we control the deterministic character of the block-interaction
matrix $p$ by varying $s$. We express such character as the mean entropy of
$p$ $\langle S(p)\rangle$ with respect to its possible outputs: $\langle
S(p)\rangle=\frac{1}{K}\sum_{k\in K}\sum_{o}p_{k}(o)\log p_{k}(o)$. The
maximum entropy for the proposed expression of $p$ is reached for $s=0.5$. We
consider 1,000 observations spread over 100 epochs. We show in Fig. 17c. that
the predictive performance of all three methods drops as the entropy
increases. This is expected, as the observations are generated from the true
model with a higher variance; each observation becomes less informative about
the underlying generative structure as $s$ grows. However, the RMSE on the
real parameters inferred using SDSBM remains low even at the maximum entropy,
meaning the model recovers the correct membership parameters.

Table 7: Numerical results on real-world datasets — Metrics abbreviations
stand for the area under the ROC curve (ROC), the Average Precision (AP), the
Normalized Coverage Error (NCE). Metrics for models stand for Simple Dynamic
SDM (SDSBM), No Coupling baseline (NC) and the classical static mixed
membership SBM (MMSBM). Overall, our approach allows for a higher predictive
power. The standard error over the folds is given in standard notation –
$0.123(12)\Leftrightarrow 0.123\pm 0.012$. The models presented in this
chapter are underlined.

1 c|  ROC AP NCE Epi SDSBM $0.9025(11)$ $0.3700(17)$ $0.1151(11)$ NC
$0.8420(22)$ $0.3435(36)$ $0.1582(19)$ MMSBM $0.8597(12)$ $0.2141(16)$
$0.1451(13)$ Lastfm SDSBM $0.8942(8)$ $0.0168(1)$ $0.1284(11)$ NC $0.8393(5)$
$0.0157(2)$ $0.1785(7)$ MMSBM $0.8647(5)$ $0.0115(2)$ $0.1493(4)$ Wiki SDSBM
$0.9759(2)$ $0.0659(9)$ $0.0459(3)$ NC $0.9092(7)$ $0.0608(10)$ $0.1195(8)$
MMSBM $0.9576(7)$ $0.0622(4)$ $0.0565(8)$ Reddit SDSBM $0.9803(3)$
$0.4295(54)$ $0.0312(3)$ NC $0.8508(5)$ $0.3598(17)$ $0.1846(7)$ MMSBM
$0.9798(2)$ $0.4269(40)$ $0.0322(3)$

######  Real-world data

Finally, we demonstrate the validity of our approach on real-world data to
argue for its usefulness and scalability. SDSBM builds on previous works on
labelled MMSBM and shares the same linear complexity
$\mathcal{O}(|R^{\circ}|)$ with $|R^{\circ}|$ the size of the dataset [42].
For our experiments, we consider the recent and documented datasets from [73],
namely the Reddit dataset (10.000 users, 984 labels, $\sim$670k observations),
the LastFm dataset (980 users, 1000 labels, $\sim$1.3M observations) and the
Wikipedia (Wiki) dataset (8227 users, 1000 labels, $\sim$157k observations).
The Reddit and Wikipedia datasets contain 1 month of data; we slice them in 1
day long temporal intervals. The LastFm dataset spans over approximately 5
years; we slice it into periods of 3 days each. In addition, we build an
additional dataset (Epi) about historical epigraphy data [26]. The dataset is
made of 117.000 Latin inscriptions comprising one or several of 18 social
statuses (slave, soldier, senator, etc.) and its location as one of 62
possible regions, along with an estimated dating spanning from 100BC to 400AD.
The goal is to guess the region where a status has been found, with respect to
time. The goal is to recover statuses diffusion in territories newly conquered
by the Roman Empire. We slice this dataset in epochs of one year each.

Evaluation is again conducted using a 5-fold cross-validation with 80% of
training data, 10% of validation data and 10% of testing data for each fold.
For each pair $(i,o_{true})$ in the test set, we query the probability for
every output $o$ given $i$ and build the confusion matrix by comparing them to
$o_{true}$. In Table .7.3.e.4, we present the results of our method compared
to the proposed baselines for various metrics: Area under the ROC curve
(AUCROC), Average Precision (AP) and Normalized Coverage Error (NCE). The
first two metrics evaluate how well models assign probabilities to
observations, and the latter evaluates the order in which possible outputs are
ranked. Overall, we see that our method exhibits a greater predictive power,
except for the Reddit dataset where the static SBM performs as well as SDSBM.
We explain this by the lack of significant temporal variation over the
considered interval. This could be expected, since the dataset comprises
roughly 80% of repeated actions [73], meaning that users do not significantly
explore new communities over a month. This result shows that SDSBM also works
well in the static case. On the other datasets, SDSBM performs better often by
a large margin, especially for the AUCROC, meaning that SDSBM is efficient at
distinguishing classes from each other. We recall that the model used here is
deliberately simplistic for demonstration purposes; low metrics do not mean
the Simple Dynamic prior does not work, but instead that it should be coupled
to a more complex model.

Figure 18: Geographic evolution of status distribution from Latin graves
(100BC - 500AD) — We applied the SDSBM to the Epigraphy dataset. We recall
that our goal is to predict a roman region (e.g., Illyria, Hispania, etc.)
given a status (e.g., Slave, Senator, etc.) and a year. We plot the temporal
evolution of statuses membership to the five manually labelled clusters (in
grey). For clarity, we removed small membership transfers from the data, which
explains why the total cluster’s population may vary from one time to another.
This plot allows us to visualize some global historical trends about the
evolution of the Roman Empire (e.g., 3rd-century crisis, the spread of
military presence in Europe, Italy’s demilitarization, etc.).

As an illustration of what SDSBM has to offer, we plot in Fig. 18 a possible
visualization of the membership ’s evolution over time for the epigraphy
dataset. On the left and the right, we show the items that are considered in
the visualization. The time goes from left (100BC) to right (500AD), and the
flows represent the membership transfers between epochs. The grey bars
represent the clusters. We manually annotated them by looking at their
composition –explicit clusters composition is given in Appendix, Section H.
From this figure, we can recover several historical facts: military presence
in Rome was scarce for most of the times considered; Italy concentrates less
military presence as time goes (due to its spread over the now extended
empire), until the 3rd-century crisis that led to its re-militarization; most
of the slaves that have been accorded an inscription are located in Italia
throughout time; the religious functions (Augusti) are evenly spread on the
territory at all times; the libertii (freed slaves) inscriptions are
essentially present in Rome and Italy, etc. Obviously, dedicated works are
needed to support these illustrative claims, and we believe SDSBM can provide
such extended comprehension of the processes at stake.

#####  .7.3.f Conclusion

We introduced a simple way to model time in dynamic valued and labelled
networks by assuming a dynamic Mixed-Membership SBM. Our method defines the
Simple Dynamic prior, ready to plug into any iteration of SIMSBM–see Section
.6.2.a. Time is considered under the single assumption that a network’s ties
do not vary abruptly over time.

We assessed the performance of the proposed method by defining the SDSBM and
evaluating it in several controlled situations on synthetic datasets. In
particular, we show that our prior shows stable performances with respect to
the dataset slicing, and that it works well under challenging conditions
(small amounts of data or high entropy blocks interaction matrix). We also
evaluated SDSBM on large scale real-world datasets and showed how accounting
for time increases yields better results than the two proposed baselines.
Finally, we illustrate an application interest on a dataset of Latin
inscriptions that indirectly narrates the evolution of the Roman Empire.

In the discussion section, we argued that our temporal prior offers great
modelling flexibility: uneven slicing of observations over time, heterogeneous
dynamic time scales for items (or clusters), time-dependent blocks-interaction
matrix, informativeness of the prior, and temporal neighbours’ dependence with
respect to the averaging function. Future works exploring these directions on
real-world data may help retrieve meaningful clusters for useful applications.
On a further note, we believe that a key interest of our approach is the
amount of data needed to get satisfactory predictive performances. As
mentioned above, this point is fundamental to several social sciences, and we
believe our approach could ease the incorporation of automated learning
methods in these fields.

###  .8 Conclusions

######  Global SIMSBM framework

In this section, we first developed a global framework, SIMSBM, that
generalizes several models from the literature as particular cases, such as
MMSBM, Bi-MMSBM, IMMSBM and T-MBM.

This results in a highly flexible model that can be applied to a broad range
of problems. In particular, we cited throughout the text several experimental
studies conducted in medicine, social behaviour and recommendation using
special cases of our model.

Our framework answers the two challenges raised in Section .6.1.b: it extends
MMSBMs to any number of entity types that can have higher-order interactions.

######  Case study with SIMSBM(2)

We then proposed to study interactions in information spread by considering a
special case of the SIMSBM: SIMSBM(2), or IMMSBM. We demonstrate that
SIMSBM(2) allows us to investigate the detail of interactions strength in
several datasets.

Our conclusions on interactions in information spread show that their effect
is not trivial and that taking them into account increases predictive
performance. However, these interactions are sparse: only 5% of significantly
interacting pairs of clusters. In most cases, non-interacting models seem to
provide an accurate enough description of output predictions – typically a
spreading action.

######  Modelling dynamic memberships of interacting entities

The previous conclusions have been made by supposing the underlying
interaction mechanisms to remain stable over one month. While this assumption
holds on some datasets, most of the time there is a need to model their
temporal evolution. In this perspective, we introduced a simple temporal
prior, ready to plug into any of the SIMSBM iterations. The time is considered
under the single assumption that the network’s ties do not vary abruptly over
time.

The performance of the proposed method is evaluated on several real-world
datasets and shows that taking the time into account improves the results on
each of them over the equivalent static model. On a further note, a key
interest of this approach is the small amount of data needed to get
satisfactory performances. This point is fundamental to several social
sciences and extends beyond the scope of interactions modelling.

######  Interactions are sparse

The overall conclusion of this section is the following: interactions are
sparse. Significant interactions take place only between a limited fraction of
cluster pairs, and between an even smaller fraction of entity pairs. This
underlines the necessity of considering clusters of entities to efficiently
model such sparse interactions.

######  Time range of interactions?

However, all the models discussed consider static interactions. The dynamic
version of SIMSBM proposed in Section .7 allows us to consider static
interactions whose expression evolves with time, but the interactions
themselves are not dynamic. To illustrate, a reasonable assumption is that the
strength of interactions is not constant over time; a tweet does not have the
same influence on a user if she saw it ten minutes or ten days earlier. The
SDSBM approach accounts for the long-term changes in the evolution of the
interacting mechanisms, but not for their immediate temporal evolution. A
Twitter user might be influenced by tweets about sports at some point, and
about politics at some other point, but being always influenced in the same
way by either of them; the membership might change but not the interaction
dynamics of each tweet. Modelling the temporal evolution of interaction
strength using MMSBMs might be possible. However, other methods that consider
time as a continuous variable instead of slicing it into episodes seem more
relevant to the task. In the next section, we propose to investigate how
entities ’ interaction strength evolves with time.

##  Chapter Temporal diffusion networks – Interactions are brief

We saw in the previous chapter that despite being sparse, interactions between
pieces of information (the entities) globally play a substantial role in
individuals’ actions: the adoption of a product, the spread of news, a
strategy choice, etc. However, the sparsity of significant interactions could
be due to temporal considerations; interactions may fade with time, making
them harder to spot for static models.

Section .9, we discuss and justify the role of time in interactions modelling.

Section .10, we show that this aspect of the underlying interaction mechanisms
remains unexplored in the literature.

Section .11, we introduce an efficient method to infer both the entities ’
interaction network and its evolution according to the temporal distance
separating interacting entities. We develop a convex model using multi-kernel
inference, named InterRate. Here, time is considered as a continuous variable,
unlike what we proposed in Section .7. The temporal evolution of interaction
intensity is what we call the interaction profile.

Section .12, we consider a timestamped sequence of exposures to entities
(e.g., URL, ads, situations) and the actions a user exerts on them (e.g.,
share, click, decision). We study how users exhibit different behaviours
according to combinations of entities they have been exposed to in the past.
We show that the joint effect of two exposures on a user is more than the
disjoint sum of their individual effect –there is an interaction. InterRate
allows for non-parametric convex optimization and can be solved in parallel.

Section .13, we show our method recovers state-of-the-art conclusions on
interaction processes on three real-world datasets. It outperforms the
proposed baselines in the inference of the underlying data generation
mechanisms. Finally, we show that interaction profiles can be visualized
intuitively, making the interpretation of the model easier.

The overall conclusion of this section is that interactions are brief. In
real-world diffusion processes, the interaction between two entities is
significant only when the two entities are close to each other in time. For
instance, a tweet does not have the same influence on the last tweet a user
saw if the first appeared ten minutes or ten days earlier in her news feed.
Typically, the intensity of an interaction decreases exponentially with the
time separating the interacting entities. Our conclusion emphasizes the
necessity of modelling the temporal aspect of interactions in social networks.

Published works: [97]

###  .9 Introduction

####  .9.1 Temporal evolution of interactions

In Chapter -, we showed that interactions play a significant role in
information spread. However, the models introduced there are all static: an
interaction between several entities has a constant intensity over time. Many
everyday examples contradict this assumption. The writing of a scientific
article is more influenced by last year’s publications than by last century’s
ones; a user on social media is more likely to answer a recent post than a
years-old one; a Spotify user is more likely to listen to a band if she
listened to it five minutes ago than if she listened to it five years ago. In
general, interactions between entities last for a given time; eventually
entities fade to a non-interacting “ground-state” as time goes forward. In
Chapter -, we called this ground-state virality, which is the probability of
an outcome in the absence of interactions.

The histograms presented in Fig.19 illustrate this assumption: the probability
of an action on an entity (like, share, comment, etc.) varies according to the
temporal distance separating any two given entities. We refer to such figures
as interaction profiles. They represent an action’s probability evolution
given the time separating the interacting entities.

The study of this quantity is a novel perspective: interactions between pieces
of information have been little explored in the literature, and no previous
work unveils trends in the information interaction mechanisms.

Figure 19: Interaction profiles between pairs of entities — Examples of
interaction profiles on Twitter; here is shown the effect of URL shortening
services migre.me (left), bit.ly (right-top), tinyurl (right-middle) and t.co
(right-bottom) on the probability of tweeting a t.co URL and its evolution
over time. This interaction profile shows, for instance, that there is an
increased probability of retweeting for a t.co URL when it appears shortly
after a migre.me one (interaction). This increase fades when the time
separation grows (no more interaction). In blue, the interaction profile
inferred by our model.

The study of interactions between entities has several applications in real-
world systems. We can mention the fields of recommender systems (the
probability of adoption is influenced by what a user saw shortly before it),
news propagation and control (when to expose users to an entity to maximize
its spreading probability) [123], advertising (same reasons as before) [22],
choice behaviour (what influences a choice and how) [27].

####  .9.2 Proposed approach

In this chapter, we propose to unveil the temporal mechanisms at stake within
those interacting processes; we infer their interaction profiles. Imagine, for
instance, that an internet user is exposed to a tweet at time $t_{1}$ and to
another at time $t_{2}>t_{1}$. We suppose that the exposure to the first one
influences the user’s sensitivity (likeliness of a retweet) to the second one
a time $t_{2}-t_{1}$ later. Modelling this process involves quantifying the
influence a tweet exerts on the other and how this influence varies with the
time separating the exposures. The representation of this probability
evolution is what we call the interaction profile –illustrated Fig.19). It
represents the influence the exposure to an entity exerts on an outcome
(click, buy, retweet, etc.) for another exposure to an entity a given time
later.

To perform this task, we introduce an efficient method to infer both the
entities ’ interaction network and its evolution according to the temporal
distance separating the interacting entities. In the proposed InterRate model,
nodes are entities, and edges between them represent the intensity of their
interaction; this intensity is a continuous-time function that depends on the
time separating two exposures to the entities (or nodes). This intensity
function is the interaction profile.

####  .9.3 Workflow

First, we review the state of the art in temporal interaction between entities
and underline the open challenges they rise in Section .10. Then, we answer
them with InterRate in Section .11, a model for inferring all the interaction
profiles between every node pair in a given set. This is performed in a
continuous-time setup using multi-kernel inference methods [28]. We show in
Section .12 that the inference of the parameters boils down to a convex
optimization problem for specific kernel families. Moreover, the problem can
be subdivided into as many subproblems as entities, which can be solved in
parallel. The convexity of the problem guarantees convergence to the
likelihood’s global optimum for each subproblem and, therefore, to the
problem’s optimal likelihood. In Section .12.2, we use InterRate to
investigate the role of interaction profiles on synthetic data and in various
corpora from different fields of research: advertisement (the exposure to an
ad influences the adoption of other ads [22]), social dilemmas (the previous
actions of one influences an other’s actions [27]) and information spread on
Twitter (the last tweets read influence what a user retweets [87]). Finally,
in Section .13, we provide analysis leads and show that our method recovers
state-of-the-art results on interaction processes on each of the three
considered datasets.

####  .9.4 Contributions

The main contributions developed in this chapter are the following:

  * •

We introduce the interaction profile. It represents the evolution of
interaction intensities as the time separating interacting entities increases.
The interaction profile is a powerful tool to understand how interactions take
place in a given corpus (see Fig.22) and it has not been developed in the
literature up to now. Its introduction is the main contribution of the present
section.

  * •

We design a convex non-parametric algorithm that can be solved in parallel,
baptized InterRate. InterRate automatically infers the interaction profile for
every node pair in a network. The collection of interaction profiles can be
interpreted as a temporal interaction network.

  * •

We show that InterRate yields better results than non-interacting or non-
temporal baseline models on several real-world datasets. Furthermore, our
model can recover several conclusions about the datasets from state-of-the-art
works.

  * •

We discuss the temporal aspect of entities interactions. Specifically, we show
that interactions are brief in real-world diffusion processes. We also recover
the interactions sparsity that has been observed in Chapter -.

###  .10 State of the art on temporal interaction network inference

####  .10.1 Temporal interactions in general

Previous efforts in investigating the role of interactions in information
diffusion have shown their importance in the underlying spreading processes.
Several works study the interaction of information with users’ attention
[129], closely linked to information overload concepts [46], but not the
interaction between the pieces of information themselves. On the other hand,
whereas most of the modelling of spreading processes is based on either no
competition [110, 96] or perfect competition [103] assumption, intermediate
competitions lead to a better description of their spread [14] –with the
example of Firefox and Chrome web browsers, whose respective popularity are
correlated but not perfectly as in [103].

####  .10.2 Modelling interactions

A significant effort has been put in elaborating complex processes to simulate
interaction on real-world networks [103, 147]. However, fewer works have been
developed to tackle interaction in information spread from a machine learning
point of view. The correlated cascade mode [145] infers an interacting
spreading process’s latent diffusion network. In this work, the interaction is
modelled by a hyper-parameter $\beta$. It tunes the intensity of interactions
according to an exponentially decaying kernel. In their conclusion, the
authors formulate the open problem of learning several kernels and the
interaction intensity parameter $\beta$, which we address in the present
chapter.

To our knowledge, the attempt the closest to our task to model the interaction
intensity parameter $\beta$ is Clash of the contagions [87]. It predicts
retweets on Twitter based on tweets seen by a user. This model estimates the
probability of retweeting a piece of information, given the last tweets a user
has been exposed to, according to their relative position in the Twitter feed.
The method suffers various flaws (scalability, non-convexity). It also defines
interactions based on an arguable hypothesis made on the prior probability of
a retweet (in the absence of interactions) that makes its conclusions about
interactions sloppy (see Section .6.1). It is worth noting that in [87], the
authors outline the problem of the inference of the interaction profile but do
so without searching for continuous trends such as the one shown in Fig.19.

In Chapter -, we addressed the various flaws observed in [87] and suggested a
more general approach to the estimation of the interaction intensity
parameters. However, we neglected the temporal aspect of interactions. To take
back the Twitter case study, it implies that in the case of a retweet at time
$t$, a tweet appearing at $t_{1}\ll t$ in the news feed has the same influence
on the retweet as a tweet that appeared at $t_{2}\approx t$; the interaction
profile is constant over time.

####  .10.3 Temporal network inference

For several years, temporal network inference has been a subject of interest.
Significant advances have been made using survival theory modelling applied to
partial observations of independent cascades of contagions [45, 48]. In this
context, an infected node tries to contaminate every other node at a rate that
is tuned by $\beta$. While this work is not directly linked to ours, it has
been a strong influence on the interaction profile inference problem we
develop here; the problems are different, but the methodology they introduce
greatly helped building InterRate (development and convexity of the problem,
analogy between interaction profile and hazard rate). Moreover, advances in
network inference based on the same works propose a multi-kernel network
inference method that we adapted to the problem we tackle here [28]. Inspired
by these works, we develop a flexible approach that allows for the inference
of the best interaction profile by using several candidate kernels.

###  .11 InterRate – Interaction dynamics

This work has been published, see [97]

####  .11.1 Problem definition

We illustrate the process to model in Fig. 20. It runs as follows: a user is
first exposed to a piece of information at time $t_{0}$. The user then chooses
whether to act on it at time $t_{0}+t_{s}$ (an act can be a retweet, a buy, a
booking, etc.); $t_{s}$ can be interpreted as the “reaction time” of the user
to the exposure, assumed constant. The user is then exposed to the next piece
of information a time $\delta t$ later, at $t_{1}=t_{0}+\delta t$ and decides
whether to act on it a time $t_{s}$ later, at $t_{1}+t_{s}$, and so on. Here,
$\delta t$ is the time separating two consecutive exposures, and $t_{s}$ is
the reaction time, separating the exposure from the possible contagion. In the
remaining of this chapter, we refer to the user’s action on an exposure (tweet
appearing in the feed, exposure to an ad, etc.) as a contagion (retweet or the
tweet, click on an ad, etc.).

Figure 20: Illustration of the interacting process — Light orange squares
represent the exposures, dark orange squares represent the exposures that are
followed by contagions and empty squares represent the exposures to the
information we do not consider in the datasets (they only play a role in the
distance between exposures when we consider the order of appearance as a time
feature). A contagion occurs at a time $t_{s}$ after the corresponding
exposure. Each new exposure arrives at a time $\delta$t after the previous
one. Contagion takes place with a probability conditioned by all previous
exposures. In the example, the contagion by A at time $t+t_{s}$ depends on the
effect of the exposure to A at times $t$ and $t-3\delta t$, and to B at times
$t-\delta t$ and $t-2\delta t$.

This choice of modelling comes with several hypotheses. First, the pieces of
information a user is exposed to appear independently from each other. It is
the main difference between our work and survival analysis literature: the
pseudo-survival of an entity is conditioned by the random arrival of pieces of
information. Therefore, users’ actions cannot be modelled as a survival
process. This assumption holds in our experiments on real-world datasets,
where users have no influence on what information they are exposed to. Second
hypothesis, the user is contaminated solely based on the previous exposures in
the feed [87, 145]. Third, the reaction time separating the exposure to a
piece of information from its possible contagion, $t_{s}$, is constant (i.e.,
the time between a read and a retweet in the case of Twitter). Importantly,
this hypothesis is a deliberate simplification of the model for clarity
purposes; relaxing this hypothesis is straightforward by extending the kernel
family, which preserves convexity and time complexity. Note that this
simplification does not always hold, as shown in recent works concluding that
response time can have complex time-dependent mechanisms [142].

####  .11.2 Likelihood

We now define the likelihood of the model whose process is described in Fig.
20. Let $t_{i}^{(x)}$ be the exposure to $x$ at time $t_{i}$, and
$t_{i}^{(x)}+t_{s}$ the time of its possible contagion. Consider now the
instantaneous probability of contagion (hazard function)
${H(t_{i}^{(x)}+t_{s}|t_{j}^{(y)},\beta_{xy})}$, that is the probability that
a user exposed to the piece of information $x$ at time $t_{i}$ is contaminated
by $x$ at $t_{i}+t_{s}$ given an exposure to $y$ at time $t_{j}\leq t_{i}$.
The matrix of parameters $\beta_{ij}$ is what the model infers. $\beta_{ij}$
is used to characterize the interaction profile between entities. We define
the set of exposures preceding the exposure to $x$ at time $t_{i}$ (or history
of $t_{i}^{(x)}$) as $\mathcal{H}_{i}^{(x)}:=\\{t_{j}^{(y)}\leq
t_{i}^{(x)}\\}_{j,y}$. Let $\mathcal{D}$ be the whole dataset such as
$\mathcal{D}:=\\{(\mathcal{H}_{i}^{(x)},t_{i}^{(x)},c_{t_{i}}^{(x)})\\}_{i,x}$.
Here, c is a binary variable that account for the contagion
($c_{t_{i}}^{(x)}=1$) or non-contagion ($c_{t_{i}}^{(x)}=0$) of $x$ at time
$t_{i}+t_{s}$. The likelihood for one exposure in the sequence given
$t_{j}^{(y)}$ is:

$\begin{split}&L(\beta_{xy}|\mathcal{D},t_{s})=P(\mathcal{D}|\beta_{xy},t_{s})=\\\
&\underbrace{H(t_{i}^{(x)}+t_{s}|t_{j}^{(y)},\beta_{xy})^{c_{t_{i}}^{(x)}}}_{\text{contagion
at $t_{i}^{(x)}+t_{s}$ due to
$t_{j}^{(y)}$}}\cdot\underbrace{(1-H(t_{i}^{(x)}+t_{s}|t_{j}^{(y)},\beta_{xy}))^{(1-c_{t_{i}}^{(x)})}}_{\text{Survival
at $t_{i}^{(x)}+t_{s}$ due to $t_{j}^{(y)}$}}\end{split}$

The likelihood of a sequence (as defined in Fig.20) is then the product of the
previous expression over all the exposures that happened before the contagion
event $t_{i}^{(x)}+t_{s}$ e.g. for all $t_{j}^{(y)}\in\mathcal{H}_{i}^{(x)}$.
Finally, the likelihood of the whole dataset $\mathcal{D}$ is the product of
$L(\beta_{x}|\mathcal{D},t_{s})$ over all the observed exposures
$t_{i}^{(x)}$. Taking the logarithm of the resulting likelihood, we get the
final log-likelihood to maximize:

$\begin{split}\ell(\beta|\mathcal{D},t_{s})=\ \ \ \ \ \ \ \ \ \ \ &\\\
\sum_{\mathcal{D}}\
\sum_{t_{j}^{(y)}\in\mathcal{H}_{i}^{(x)}}&c_{t_{i}}^{(x)}\log\left(H(t_{i}^{(x)}+t_{s}|t_{j}^{(y)},\beta_{xy})\right)\\\
+\
(1-c_{t_{j}}^{(y)})&\log\left(1-H(t_{i}^{(x)}+t_{s}|t_{j}^{(y)},\beta_{xy})\right)\end{split}$
(35)

####  .11.3 Proof of convexity

The convexity of a problem guarantees to retrieve its optimal solution and
allows using dedicated fast optimization algorithms.

######  Proposition 1.

The inference problem $\min_{\beta}-\ell(\beta|\mathcal{D},t_{s})\
\forall\beta\geq 0$, is convex in all of the entries of $\beta$ for any hazard
function that obeys the following conditions:

$\begin{cases}&H^{\prime 2}\geq H^{\prime\prime}H\\\ &H^{\prime
2}\geq-H^{\prime\prime}(1-H)\\\ &H\in\left]0;1\right[\end{cases}$ (36)

where ′ and ′′ denote the first and second derivative with respect to $\beta$,
and H is the shorthand notation for
$H(t_{i}^{(x)}+t_{s}|t_{j}^{(y)},\beta_{xy})\ \forall i,j,x,y$.

###### Proof.

The negative log-likelihood as defined in Eq.35 is a summation of $-\log H$
and $-\log(1-H)$; therefore $H\in\left]0;1\right[$. The second derivative of
these expressions according to any entry $\beta_{mn}$ (noted ′′) reads:

$\begin{cases}&\left(-\log
H\right)^{\prime\prime}=\left(\frac{-H^{\prime}}{H}\right)^{\prime}=\frac{H^{\prime
2}-H^{\prime\prime}H}{H^{2}}\\\
&\left(-\log(1-H)\right)^{\prime\prime}=\left(\frac{H^{\prime}}{1-H}\right)^{\prime}=\frac{H^{\prime
2}+H^{\prime\prime}(1-H)}{(1-H)^{2}}\end{cases}$ (37)

The convexity according to a single variable holds when the second derivative
is positive, which leads to Eq.36. The convexity of the problem then follows
from the composition rules of convexity. ∎

Several functions obey the conditions of Eq.36, such as the exponential
($e^{-\beta t}$), Rayleigh ($e^{-\frac{\beta}{2}t^{2}}$), power-law
($e^{-\beta\log t}$) functions, and any log-linear combination of those [28].
These functions are standard in survival theory literature [47].

The final convex problem can then be written $\min_{\beta\geq
0}-\ell(\beta|\mathcal{D},t_{s})$. An interesting feature of the proposed
method is that the problem can be subdivided into N convex subproblems that
can be solved independently (one for each piece of information). To solve the
subproblem of the piece of information $x$, that is to find the vector
$\beta_{x}$, we need to consider only the subset of $\mathcal{D}$ where $x$
appears. Explicitly, each subproblem consists in maximizing Eq.35 over the set
of observations
$\mathcal{D}^{(x)}:=\\{(\mathcal{H}_{i}^{(x)},t_{i}^{(x)},c_{t_{i}}^{(x)})\\}_{i}$.

###  .12 Experiments

####  .12.1 Experimental setup

#####  .12.1.a Kernel choice

######  Gaussian RBF kernel family (IR-RBF)

Based on [28], we consider a log-linear combination of Gaussian radial basis
function (RBF) kernels as hazard function. We also consider the time-
independent kernel needed to infer the base probability of contagion discussed
in the section “Background noise in the data” below. The resulting hazard
function is then:

$\log
H(t_{i}^{(x)}+t_{s}|t_{j}^{(y)},\beta_{ij})=-\beta_{ij}^{(bg)}-\sum_{s=0}^{S}\frac{\beta^{(s)}_{ij}}{2}(t_{i}+t_{s}-t_{j}-s)^{2}$

The parameters $\beta^{(s)}$ of Rayleigh kernels are the amplitude of a
Gaussian distribution centered on time $s$. The parameter S represents the
maximum time shift we consider. In our setup, we set S=20. We think it is
reasonable to assume that an exposition does not significantly affect a
possible contagion 20 steps later. The parameter $\beta_{ij}^{(bg)}$
corresponds to the time-independent kernel –base probability of contagion by
i, or virality. The formulation allows the model to infer complex
distributions from a reduced set of parameters whose interpretation is
straightforward.

######  Exponentially decaying kernel (IR-EXP)

We also consider an exponentially decaying kernel. We consider the following
form for the hazard function and refer to this modelling as IR-EXP:

$\log
H(t_{i}^{(x)}+t_{s}|t_{j}^{(y)},\beta_{xy})=-\beta_{ij}^{(bg)}-\beta_{ij}(t_{i}+t_{s}-t_{j})$

where $\beta_{ij}^{(bg)}$ once again accounts for the background noise in the
data discussed further in this section.

#####  .12.1.b Parameters learning

Datasets are made of sequences of exposures and contagions, as shown in
Fig.20. To assess the robustness of the proposed model, we apply a 5-fold
cross-validation method. After shuffling the dataset, we use 80% of the
sequences as a training set and the 20% left as a test set. We repeat this
slicing five times, taking care that an interval cannot be part of the test
set more than once. The optimization is made in parallel for each piece of
information via the convex optimization module for Python CVXPY.

We also set the time separating two exposures $\delta t$ as constant. It means
that we consider only the order of arrival of exposures instead of their
absolute arrival time. The hypothesis that the order of exposures matters more
than the absolute exposure times has already been used with success in the
literature [87]. Besides, in some situations, the exact exposure time cannot
be collected, while the exposures’ order is known. For instance, in a Twitter
corpus, we only know in what order a user reads her feed, unlike the exact
time she read each of the posts. However, from its definition, our model works
the same with non-integer and non-constant $\delta t$ in datasets where
absolute time matters more than the order of appearance.

#####  .12.1.c Background noise in the data

Figure 21: Underlying generation process vs observed data — The red curve
represents the underlying probability of contagion by C given an exposure
observed $\Delta t$ steps before C. The orange bars represent the observed
probability of such events. We see that there is a noise $P_{0}(C)$ in the
observed data. The underlying generation process can then only be observed in
the dataset when its effect is larger than some threshold $P_{0}(C)$.

Because the dataset is built looking at all exposure-contagion correlations in
a sequence, there is inherent noise in the resulting data. To illustrate this,
we look at the illustrated example Fig.20 and consider the exposure to C
leading to a contagion happening at time $t$. We assume that in the underlying
interaction process, the contagion by C at time $t+t_{s}$ only took place
because C appeared at time $t$. However, when building the dataset, the
contagion by C is also attributed to A appearing at times $t-\delta t$,
$t-3\delta t$ and $t-6\delta t$, and to B appearing at times $t-4\delta t$ and
$t-5\delta t$. It induces noise in the data. In general, for any contagion in
the dataset, several observations (pair exposure-contagion) come from the
random presence of entities unrelated to this contagion.

We now illustrate how this problem introduces noise in the data. In Fig.21, we
see that the actual underlying data generation process (probability of a
contagion by C given an exposure present $\Delta t$ step earlier) does not
exactly fit the collected resulting data: the data gathering process induces a
constant noise whose value is noted $P_{0}(C)$ –that is the average
probability of contagion by C. Thus, the interaction effect can only be
observed when its associated probability of contagion is larger than
$P_{0}(C)$. Consequently, the performance improvement of a model that accounts
for interactions may seem small compared with a baseline that only infers
$P_{0}(C)$. That is what we observe in the experimental section. However, in
this context, a small improvement in performance shows an extended
comprehension of the underlying interacting processes at stake (see Fig.21,
where the red line explains the data better than a constant baseline). Our
method efficiently infers $P_{0}(C)$ via a time-independent kernel function
$\beta_{i,j}^{P_{0}(i)}$. Inferring $P_{0}(i)$ is in line with the discussion
of Section .6.1 on the necessity to infer the virality along with interaction
terms.

#####  .12.1.d Evaluation criteria

The main difficulty in evaluating these models is that interactions might
occur between a small number of entities only. It is the case here, where many
pairs of entities have little to no interaction (see the Discussion section).
This makes it difficult to evaluate how good a model is at capturing them. To
this end, our principal metric is the residual sum of squares (RSS). The RSS
is the sum of the squared difference between the observed and the expected
frequency of an outcome. This metric is particularly relevant in our case,
where interactions may occur between a small number of entities: any deviation
from the observed frequency of contagion is accounted for, which is what we
aim at predicting here. We also consider the Jensen-Shannon (JS) divergence;
the JS divergence is a symmetric version of the Kullback–Leibler divergence,
which makes it usable as a metric [91].

We finally consider the best-case F1-score (BCF1) of the models, that is, the
F1-score of the best scenario of evaluation. It is not the standard F1 metric
(that poorly distinguishes the models since few interactions occur), although
its computation is similar. Explicitly, it generalizes F1-score for comparing
probabilities instead of comparing classifications; the closer to 1, the
closer the inferred and observed probabilities. It is derived from the best-
case confusion matrix, whose building process is as follows: we consider the
set of every information that appeared before information i at time $t_{i}$ in
the interval, that we denote $\mathcal{H}_{i}$. We then compute the contagion
probability of i at time $t_{i}+t_{s}$ to every exposure event
$t_{j}^{(y)}\in\mathcal{H}_{i}$. Confronting this probability with the
observed frequency f of contagions of i at time $t_{i}+t_{s}$ given
$t_{j}^{(y)}$ among N observations, we can build the best-case confusion
matrix. In the best-case scenario, if out of N observations the observed
frequency is f and the predicted frequency is p, the number of True Positives
is $N\times\min\\{p,f\\}$, the number of False Positives is
$N\times\max\\{p-f,0\\}$, the number of True Negatives is
$N\times\min\\{1-p,1-f\\}$, the number of False Positives is
$N\times\max\\{f-p,0\\}$.

Finally, when synthetic data is considered, we also compute the mean squared
error of the $\beta$ matrix inferred according to the $\beta$ matrix used to
generate the observations, that we note MSE $\beta$.

We purposely ignore evaluation in prediction because, as we show later,
interactions influence quickly fades over time: probabilities of contagion at
large times are mainly governed by the background noise discussed in previous
sections. Therefore, it would be irrelevant to evaluate our approach’s
predictive power on the whole range of times where it does not bring any
improvement over a naive baseline (see Fig.19). A way to alleviate this
problem would be to make predictions only when interaction effects are
above/below a certain threshold (at short times, for instance). However, such
an evaluation process would be debatable. Here, we choose to focus on the
descriptive aspect of InterRate.

#####  .12.1.e Baselines

######  Naive baseline

For a given piece of information i, the contagion probability is defined as
the number of times a user acts on it divided by the number of its
occurrences.

######  Clash of the contagions

We use the work presented in [87] as a baseline. In this work, the authors
model the probability of a retweet given the presence of a tweet in a user’s
feed. This model does not look for trends in the way interactions take place
(it does not infer an interaction profile), considers discrete time steps
(while our model works in a continuous-time framework), and is optimized via a
non-convex SGD algorithm (which does not guarantee convergence towards the
optimal model). More details on implementation are provided in Appendix,
Section I.

######  IMMSBM

The Interactive Mixed-Membership Stochastic Block Model (IMMSBM) is a model
that takes interactions between pieces of information into account to compute
the probability of a (non-)contagion –see Section .6.3.a. Note that this
baseline does not take the position of the interacting pieces of information
into account (time-independent) and assumes that interactions are symmetric
(the effect of A on B is the same as B on A).

######  ICIR

The Independent Cascade InterRate (ICIR) is a reduction of our main IR-RBF
model to the case where interactions are not considered. We consider the same
dataset, enforcing the constraint that off-diagonal terms of $\beta$ are null.
The (non-)contagion of a piece of information i is then determined solely by
the previous exposures to i itself.

####  .12.2 Results

#####  .12.2.a Synthetic data

Table 8: Experimental results on synthetic data — Our model outperforms all of
the baselines in almost every dataset for every evaluation metric. The
standard deviations of the 5 folds cross-validation are negligible. The models
presented in this chapter are underlined.

RSS JS div. BCF1 MSE $\beta$ Synth-20 IR-RBF $18.415\,145$ $0.002\,284\,2$
$0.918\,832$ $0.000\,503\,27$ ICIR $139.5926$ $0.009\,980\,1$ $0.827\,047\,7$
$0.015\,881\,1$ 2-6 Naive $145.5132$ $0.010\,378\,5$ $0.822\,139$ CoC
$123.0583$ $0.009\,383\,8$ $0.822\,015\,7$ IMMSBM $222.055\,495$
$0.017\,287\,5$ $0.726\,541\,3$ Synth-5 IR-RBF $0.116\,926\,7$ $0.000\,217\,4$
$0.974\,213\,3$ $0.005\,297\,7$ ICIR $8.266\,074\,4$ $0.008\,117\,4$
$0.849\,863\,7$ $0.019\,206\,1$ 2-6 Naive $10.026\,425\,4$ $0.009\,955\,6$
$0.821\,412\,8$ CoC $0.115\,429\,7$ $0.000\,197\,4$ $0.976\,287\,2$ IMMSBM
$11.693\,641\,5$ $0.013\,622\,3$ $0.769\,258\,0$

######  Data generation

We generate synthetic data according to the process described in Fig.20 for a
given $\beta$ matrix using the RBF kernel family. First, we generate a random
matrix $\beta$, whose entries are between 0 and 1. A piece of information is
then drawn with uniform probability and can result in a contagion according to
$\beta$, the RBF kernel family and its history. We simulate the outcome by
drawing a random number and finally increment the clock. The process then
keeps on by randomly drawing a new exposure and adding it to the sequence. We
set the maximum length of intervals to 50 steps and generate datasets of
20,000 sequences.

######  Numerical results

We present in Tab. 8 the results of the various models with generated
interactions between 20 (Synth-20) and 5 (Synth-5) entities. The interactions
are generated using the RBF kernel, hence the fact we are not evaluating the
IR-EXP model –its use would be irrelevant. The InterRate model outperforms the
proposed baselines for every metric considered. It is worth noting that
performances of non-interacting and/or non-temporal baselines are good on the
JS divergence and F1-score metrics due to the constant background noise
$P_{0}$. For cases where interactions do not play a significant role, IMMSBM
and Naive models perform well by fitting only the background noise. By
contrast, the RSS metric distinguishes very well the models that are better at
modelling interactions.

Note that while the baseline [87] yields good results when few interactions
are simulated (Synth-5), it performs as bad as the naive baseline when this
number increases (Synth-20). This is due to the non-convexity of the proposed
model, which struggles to reach a global maximum of the likelihood even after
100 runs (see Appendix, Section I for implementation details).

#####  .12.2.b Real data

We consider 3 real-world datasets. For each dataset, we select a subset of
entities that are likely to interact with each other. For instance, it has
been shown that the interaction between the various URL shortening services on
Twitter is non-trivial [145].

We provide details on the way datasets have been built from raw data. For each
of the real-world datasets, we choose to consider only the order of the
various entities’ apparition instead of their absolute appearance times. It
implies setting the time separating two successive exposures as constant, that
we note $\delta t$. This choice is supported by state-of-the-art works [87],
and we observed in our experiments that it is more relevant than considering
absolute times. Besides, we do not consider the first 10 pieces of information
of any sequence to avoid boundary effects (the first 5 steps for the PD
dataset): the history of exposures is incomplete in this case and could lead
to biased results. For each dataset entities list, the number before the
entity name is the key used in Fig. 22. The entities subsets have been chosen
by computing the co-occurrence matrix of all the entities and then selecting
the ones that are part of a cluster using a K-means algorithm. The datasets
are:

######  Datasets

  * •

Twitter dataset [59]: a collection of all the tweets containing URLs that have
been posted on Twitter during October 2010, with the associated followers’
network. A tweet read by a user in her feed is an exposition, and its possible
retweet is a contagion. We consider only the URLs associated with the
following URL shortening websites, the same as in [145]: {0: migre.me, 1:
bit.ly, 2: tinyurl, 3: t.co}. The final dataset is made of 104,349 sequences
of average length 53.5 steps (1 step = $t_{s}$), for 1,276,670,965 observed
interactions.

  * •

Prisoner’s dilemma dataset (PD): contains ordered sequences of repeated
Prisoner’s dilemma games between two players. From the dataset introduced in
[89], we consider the sub-dataset noted BR-risk 0 (first entry of Tab.2 in the
reference); we choose this subset to have decisions made in a homogeneous
context, where players struggle in a dilemma that is hard to solve (which
depends on the combination of the parameters T, R, S and P discussed further).
Within each round, the two players can defect or cooperate. Each duel is made
of 10 rounds. If both cooperate, the reward R is high; if both defect, the
reward P is low; if one player cooperates while the other defects, this one
gets a penalty S, while the other gets a reward T. To make the game a Prisoner
dilemma, the variables have to obey T$>$R$>$P$>$S. We refer to the combination
of players’ actions ("the user cooperated, and the opponent defected at time
t") as exposures and to the defect actions of the player in the following
round as a contagion. We defined the action of cooperating as a non-contagion.
We therefore have 4 possible situations ({0: Player cooperated, and opponent
defected, 1: Both players defected, 2: Both players cooperated, 3: Player
defected and opponent cooperated}) and 2 possible outcomes (Player cooperates
or defects). The final dataset is made of 2,337 sequences of average length of
10.0 steps, for 189,297 observed interactions.

  * •

Taobao dataset (Ads): contains all ads exposures for 1,140,000 randomly
sampled users from the website of Taobao for 8 days (5/6/2017-5/13/2017) [22].
Taobao is one of the largest e-commerce websites and is owned by Alibaba. Each
exposure is associated with the corresponding timestamp and user’s action
(click on the ad or not). A click is considered a contagion. The subset of ads
we consider is: {0: 4520, 1: 4280, 2: 1665, 3: 4282}. The resulting dataset is
made of 87,500 sequences of average length of 23.9 steps, for 240,932,401
observed interactions.

Table 9: Experimental results on real-world data — Our model outperforms all
of the baselines in almost every dataset for every evaluation metric. The
standard deviations of the 5 folds cross-validation are negligible. The models
presented in this chapter are underlined.

RSS JS div. BCF1 Twitter IR-RBF $0.001\,467\,6$ $0.000\,058\,2$ $0.983\,202$
IR-EXP $0.001\,135\,9$ $0.000\,048\,8$ $0.986\,201$ ICIR $0.013\,710\,0$
$0.000\,629\,3$ $0.961\,401$ 2-5 Naive $0.016\,086\,6$ $0.000\,725\,2$
$0.937\,949\,9$ CoC $0.001\,676\,5$ $0.000\,067\,2$ $0.957\,223\,0$ IMMSBM
$0.014\,730\,5$ $0.000\,682\,9$ $0.954\,292\,3$ PD IR-RBF $1.126\,79$
$0.007\,583$ $0.978\,903$ IR-EXP $1.552\,556$ $0.008\,668\,6$ $0.966\,056$
ICIR $3.535\,85$ $0.018\,225$ $0.938\,12$ 2-5 Naive $3.6527$ $0.019\,146\,6$
$0.945\,45$ CoC $1.240\,857\,3$ $0.008\,087\,6$ $0.973\,645\,3$ IMMSBM
$20.3773$ $0.087\,010\,4$ $0.767\,153$ Ads IR-RBF $0.004\,338$ $0.000\,043\,2$
$0.981\,43$ IR-EXP $0.003\,016$ $0.000\,029\,7$ $0.985\,24$ ICIR $0.098\,337$
$0.000\,848\,1$ $0.965\,88$ 2-5 Naive $0.145\,311\,1$ $0.001\,262\,9$
$0.912\,58$ CoC $0.004\,514\,6$ $0.000\,045\,1$ $0.974\,115\,3$ IMMSBM
$0.015\,465$ $0.000\,153\,1$ $0.954\,27$

######  Numerical results

The results on real-world datasets are presented in Tab.9. We see that the
IMMSBM baseline performs poorly on the PD dataset: either considering the time
plays a consequent role in the probability of contagion, or interactions are
not symmetric. Indeed, as we saw in the previous chapter, a core hypothesis of
the IMMSBM is that the effect of exposition A on B is the same as B on A,
whichever is the time separation between them. In a prisoner’s dilemma game
setting, for instance, we expect that a player does not react in the same way
to defection followed by cooperation as to cooperation followed by defection,
a situation for which the IMMSBM does not account. When there are few
entities, the CoC baseline performs as good as IR, but it fails when this
number increases; this is mainly due to the non-convexity of the problem that
does not guarantee convergence towards the optimal solution. Overall, the
InterRate models yield the best results on every dataset.

###  .13 Discussion

####  .13.1 Exponential interaction profiles

In Fig.22, we represent the interaction intensity over time for every pair of
information considered in every corpus fitted with the RBF kernel model. The
intensity of the interactions is the inferred probability of contagion minus
the base contagion probability in any context: $P_{ij}(t)-P_{0}(i)$. We recall
that $P_{0}(i)$ is inferred along with the other parameters; it is in line
with the discussion of Section .6.1 on the necessity to infer the virality
along with interaction terms. Therefore, we can determine the characteristic
range of interactions, investigate recurrent patterns in interactions, whether
the interaction effect is positive or negative, etc.

Overall, we understand why the EXP kernel performs as good as the RBF on the
Twitter and Ads datasets: interactions tend to have an exponentially decaying
influence over time. However, this is not the case in the PD dataset: the
effect of a given interaction is very dependent on its position in the history
(pike on influence at $\Delta t=3$, shift from positive to negative influence,
etc.).

Figure 22: Visualization of the interaction profiles — Intensity of the
interactions between every pair of entities according to their time separation
(one line is one pair’s interaction profile, similar to Fig.19 seen "from the
top"). A positive intensity means that the interaction helps the contagion,
while a negative intensity means it blocks it. Note that we represented
discontinuous times for visualization proposes, but the inferred kernel is
continuous, as in Fig. 19. The key linking numbers on the y-axis to names for
each dataset is provided in the main text, Section .12.2.b.1.

####  .13.2 Recovering state of the art conclusions

In the Twitter dataset, the most substantial positive interactions occur
before $\Delta$t=3. This finding agrees with previous works, which stated that
the most informative interactions within the Twitter URL dataset occur within
the 3 time steps before the possible retweet [87]. We also find that the vast
majority of interactions are weak, matching with previous study’s findings see
Chapter - and [87]. However, it seems that tweets still exert influence even a
long time after being seen, but with lesser intensity.

In the Prisoner’s Dilemma dataset, players’ behaviours are heavily influenced
by the previous situations they have been exposed to. For instance, in the
situation where both players cooperated in the previous round (pairs 2-x, 3rd
section in Fig.22-PD). The probability that the player defects is then
significantly increased if both players cooperated or if one betrayed the
other exactly two rounds before but decreased if it has been two rounds that
players both cooperate.

Finally, we find that the interactions play a lesser role in the clicks on
ads. We observe a slightly increased probability of clicking on every ad after
direct exposure to another one. We also observe a globally decreasing
probability of click when two exposures are distant in time, which agrees with
previous work’s findings [22]. Finally, the interaction profile is very
similar for every pair of ads; we interpret this as a similarity in users’ ads
perception.

We showed that for each of the considered corpus, considering the interaction
profile provides an extended comprehension of choice adoption mechanisms and
retrieves several state-of-the-art conclusions. The proposed graphical
visualization also provides an intuitive view of how the interaction occurs
between entities and the associated trends, hence supporting its relevance as
a new tool for researchers in a broad meaning.

###  .14 Conclusions

######  Modelling the temporal aspect of interactions

We showed in previous works that interactions between entities play a
substantial role in individuals’ actions. However, the temporal aspect of
interactions has been little explored in the literature, despite their
importance in modelling real-world processes [87, 145]. In this chapter, we
filled the gap and introduced an efficient convex model that investigates the
temporal aspect of interactions.

Unlike previous models, our method accounts for both the interaction effects
and their influence over time (the interaction profile). We showed InterRate
yields better results on synthetic and real-world datasets. Therefore, taking
the temporal aspect of interactions provides a finer comprehension of the
processes at stake.

However, the method introduced in this section presents a major flaw. As it
has been stated in the previous chapter, most entities do not interact with
each other –interactions are rare. Therefore, in this chapter, we restricted
to the study of small subsets of nodes which we knew were strongly interacting
together. Our method shows that time is of importance in the modelling of
interactions in spreading processes, but it cannot generalize these
interaction patterns to a class of entities; it lacks clustering. It makes the
results less interpretable, as the model does not look for a trend in
individual entities’ behaviour, in which we are interested.

######  Recovering conclusions on temporal interactions

We showed that InterRate manages to recover several state-of-the-art
conclusions that have been made on the same datasets. On Twitter, most
significant interactions happen at small times, and the majority of
interactions are weak. On the advertisement corpus, our model highlighted the
effects of marketing fatigue, which is a short increase in the probability of
the click on an ad immediately after a first exposition, and a decreasing
probability when expositions get more distant. On the prisoner’s dilemma
dataset, we showed our model recovers clues to a deterministic circular
reasoning between players.

######  Interactions are short

The overall conclusion of this section is that interactions in spreading
processes are brief. Typically, the intensity of an interaction on Twitter
decreases exponentially with the time separating the interacting entities.

######  Towards proper interactions modelling

All the models proposed in both Chapter - and Chapter .8.0..5 make strong
assumptions about the type of entities considered. Typically, an entity has
been identified by a link on Twitter, independently from its content. However,
two different links may be about the same topic and thus be regarded
identically by other entities they interact with. It seems reasonable to
assume that entities carrying identical semantic meanings are regarded in the
same way by other entities they interact with.

As such, a more refined model should comprise three crucial elements to
interaction modelling. Entities should be given a more complete definition
that includes their semantic meaning rather than only their identifier. These
entities must be clustered together according to their dynamics to be able to
recover significant interactions. Interactions should be dynamic, meaning
their value should depend on the time separating the interacting entities. In
the next chapter, we derive such refined model that answers all these
challenges.

##  Chapter Dirichlet-Hawkes Processes - Modelling rare and brief interactions

The conclusions drawn from Chapter - and Chapter .8.0..5 are that interactions
are sparse (we need clusters to model them) and that interactions are brief
(we need to consider time). Besides, we discussed the fact that interacting
entities may have a short lifespan, and that their semantic content must be
taken into account. In this last chapter, we present the steps paving our way
to a single model that answers all of these challenges.

Section .15, we frame our approach as an in-depth modification of an existing
Bayesian prior, the Dirichlet-Hawkes Process.

Section .16, we detail how to get this model by merging Dirichlet Processes
with Hawkes processes, and highlight its main limits.

Section .17, as an indirect way to overcome these limits, we explore
alternative forms of Dirichlet Processes and end up with an expression that
alleviates a major feature of the Dirichlet Processes, their “rich-get-richer”
property. We show that such a hypothesis is not always a relevant modelling
choice. We propose the Powered Dirichlet Process as a way to directly control
the “rich-get-richer” assumption.

Section .18, we then incorporate the newly proposed Powered Dirichlet Process
into the standard Dirichlet-Hawkes process to create the Powered Dirichlet-
Hawkes Process. We show our formulation yields significantly better results
than state-of-the-art models when temporal information or textual content is
weakly informative and alleviates the hypothesis that textual content and
temporal dynamics are always perfectly correlated. Our approach eventually
allows us to correctly model self-interactions within a cluster.

Section .19, we extend our approach to the multivariate case, which allows us
to explore not only self-interactions, but interactions between all clusters.
We present the challenges that arise from such an extension and how we
overcome them.

Section .20, finally, we perform large-scale experiments on a real-world
Reddit dataset, made of all of the headlines published on news subreddits
throughout the whole year 2019. We confirm the findings of previous chapters
and conclude that interactions play a minor role in this particular dataset.

Published works:

  * •

[99] in Section .18

  * •

[100] in Section .20

###  .15 Introduction

####  .15.1 How to properly model interactions

In the two previous chapters, we underlined several challenges that arise when
modelling interactions between entities in information spread. In Chapter -,
we showed that interactions between pairs of entities are sparse –few entities
interact significantly. In Chapter .8.0..5, we showed that pair-interaction
between entities are brief –entities have to be close in time for significant
interactions to happen.

These conclusions invite us to rethink the definition of what we consider an
entity. In previous chapters, an entity is identified using a unique ID. It
can be a link, a word, or a song, but it take multiple forms. This is a strong
assumption since several of such entities can carry an identical semantic
meaning. Let link A222https://www.youtube.com/watch?v=dQw4w9WgXcQ and link
B333https://www.youtube.com/watch?v=oHg5SJYRHA0 point to the same content. The
previous definition of entities would consider them as distinct despite their
semantic meaning being strictly identical. This could become problematic in
contexts where such entities can express the same thing in several manners.
Typically on Twitter or Reddit, opinions on recent news are expressed in
various forms but may express the same thing. Another advantage of a more
complete definition is the case of repeated information. Imagine a Twitter
account that provides a daily weather forecast; each of the reports would be
considered as a different entity whereas users are expected to react in
similar ways to this or that forecast (sunny, rainy, cloudy, etc.), but not to
every one of them.

Therefore, we need to account for entities’ content. This is even more crucial
in situations where pieces of information appear and disappear at a high rate,
such as on online media platforms. The half-life of a Tweet is around 18
minutes, around 30 minutes for a Facebook post, 19h on Instagram, 6 days on
Youtube, etc., which are all short in some perspective. Short lifespans
provide less information on which to learn the interacting processes. However,
aggregating these pieces of information together using both their content and
dynamics would provide enough data to study interactions.

To summarize, we need to develop a model that creates clusters of entities
based on both their content and temporal interactions.

####  .15.2 Objective

Our final goal is to model the temporal interaction between entities in real-
world large-scale datasets. Nowadays, online information is generated at an
unprecedented rate. At the time of writing, every minute, 500,000 comments are
posted on Facebook, 400 hours of videos are uploaded on Youtube, and 500,000
tweets are published on Twitter.

As we stated in Chapter -, we need to cluster this mass of entities together
to make sense of this mass of information. As we saw in the introduction,
Section .15.1, this clustering must take entities ’ content into account
–typically their textual content. Many clustering algorithms are based on text
similarity, that is, how similar the words of two published documents are [19,
8, 107].

However, in Chapter .8.0..5, we saw that these clusters must also include a
temporal dimension to reflect the interacting processes at stake. Another
variable to account for is thus the time of publication [18, 28].

####  .15.3 Proposed approach

Many clustering models that claim to model dynamic clusters do not in fact
explicitly account for time. At a given time $t$, they sample a subset of
recent observations according to a temporal sampling function and then learn a
static model for this time-step using the selected data only [18, 2, 141].
However, sampling observations over time implies defining a sampling function
that might not correctly model the temporal dynamics at stake. In general,
time is not used by the model but instead only reduces the data provided to
static models. It has been argued that such modelling is not fit to account
for the arrival of documents in continuous-time settings [30].

In [30], the authors combine techniques of standard textual clustering with
temporal point processes. The idea is to infer the time-sampling functions for
data selection jointly with the clustering model with which they are
associated. They derive the Dirichlet-Hawkes process (DHP) prior for
clustering document streams by using jointly textual and temporal information
in the cluster inference. In this model, clusters are self-stimulated. It
means that each entity they comprise gets associated with a temporal function
representing the probability of a new entity appearing at all times. This can
be interpreted as the diagonal of the interaction matrix in Fig. 12, but in
its temporal version –with each case being associated with an interaction
profile, see Section .13.

This model [30] seems to fit our task very well. However, we cannot use it as
such, as it suffers from some limitations and assumptions. For instance, it
has been argued that this method cannot handle limit cases where text is less
informative –e.g., short texts, overlapping vocabularies [141].

####  .15.4 Workflow

In this chapter, we will first detail the Dirichlet-Hawkes Process (DHP) prior
introduced in [30] and point out its limits (Section .16). In particular, we
will show that their inclusion of the temporal dimension results from an
arbitrary choice, and that tuning the influence of time allows us to recover
better results on challenging datasets.

To answer this problem, we will first question the Dirichlet Process (DP) on
which DHP is built (Section .17). We demonstrate that alternative, more
flexible variations of DP are possible and that they allow for better
modelling performances; we call this alternative process the Powered Dirichlet
Process (PDP).

We then develop the Powered Dirichlet-Hawkes Process (PDHP) as the
reformulation of the DHP model in terms of PDP (Section .18). We show this
novel, more flexible formulation yields significantly better results than the
original DHP and alleviate several hypotheses the original model made.

Finally, in Section .19, we extend the PDHP to the multivariate case (MPDHP):
we allow clusters to have temporal interactions with each other, and not only
with themselves. This is equivalent to a temporal version of every matrix
presented in Page 12-Fig. 12 –and not only their diagonal as for DHP.

The final form of the proposed approach is able to model interactions that
are:

  * •

sparse –by clustering entities together (Chapter -)

  * •

temporal –by associating each cluster an interaction profile, (Chapter
.8.0..5)

  * •

content sensitive –by considering documents instead of entities identifiers

We finally conduct a large-scale study of the multivariate temporal
interactions between clusters of documents (entities with a content) in
Section .20. We use 12 months of Reddit data specific to news subreddits, and
conclude that the impact of interactions is small in this specific dataset.
Reassuringly, we also confirm the findings of previous chapters on the
sparsity and persistence of interactions in social media.

###  .16 State of the art and limits

####  .16.1 A brief overview of temporal clustering of textual documents

The use of temporal dimension in document clustering has been studied on many
occasions; a notable spike of interest happened in 2006. Many authors tackled
the problem of inferring time-dependent clusters from models based on LDA [18,
127, 63]. However, most of these models are parametric, meaning the number of
clusters is fixed at the beginning of the algorithm. Depending on the
considered time range and the dataset, the number of clusters needs to be
fine-tuned with several independent runs, making them hardly usable for many
real-world applications. In all three references cited, the authors mention
that a non-parametric version of the model might be derivable.

In 2008, A. Ahmed et al. proposed the Recurrent Chinese Restaurant Process
(RCRP) as an answer to this problem [2]. Instead of considering a fixed-size
dataset, this model can handle a stream of documents arriving in chronological
order, and the number of clusters is automatically updated. In this model,
time is split into episodes to capture the temporal aspect of cluster
formation; it considers an integer count of publications within a given time
window. A later version of the model from 2010, the Distance-Dependent Chinese
Restaurant Process (DD-CRP), tries to alleviate this approximation by
replacing fixed-time episodes with a continuous-time sampling function [17].
However, the model still considers integer counts with only their distribution
over time changing. The temporal dimension is not explicitly modelled, but
instead used as a filter for the data fed to the model. Such slicing (even
based on continuous sampling functions) can induce strong bias in the temporal
modelling. One of these biases is illustrated in Fig. 23, where observations
in the same time slice can be further in time than observations in different
ones. Thus, the model is not designed to consider every temporal information
in a continuous-time setting.

Figure 23: Slicing the data to consider time can introduce bias – Here, events
in a same time slice can be further in time than observations in different
ones. Explicitly modelling time using temporal processes helps getting rid of
this bias.

In 2015, N. Du et al. answered this problem by combining the Dirichlet process
with the Hawkes process, used to model the appearance of events in a
continuous-time setting. The key idea is to replace the counts of a Dirichlet
process with the intensity function of the Hawkes process. The resulting
Dirichlet-Hawkes process (DHP) is then used as a prior for clustering
documents appearing in a continuous-time stream. The inference is realized
with a Sequential Monte-Carlo (SMC) algorithm. Following DHP, two articles
have been published extending the idea: the Hierarchical Dirichlet-Hawkes
process [83] and the Indian Buffet Hawkes process [118]. Another work proposed
an EM algorithm for the inference [134] (it uses a heuristic method to update
the number of clusters and cannot handle a stream of documents).

####  .16.2 Dirichlet-Hawkes Process

#####  .16.2.a Dirichlet Process

The Dirichlet Process is typically used in clustering models. It naturally
yields a partition over a possibly infinite number of clusters. It is used as
a prior on entities’ membership–such as the ones introduced in Section .6.1.b.

A well-known metaphor for the Dirichlet process is referred to as “Chinese
restaurant”. The corresponding process is named “Chinese Restaurant Process”
(CRP). It can be illustrated as follows: when a $n^{th}$ client arrives in a
Chinese restaurant, she will sit at one of the $K$ already occupied tables
with a probability proportional to the number of persons already sitting at
this table. She can also go to a new table in the restaurant and be the first
client to sit there with a probability inversely proportional to the total
number of clients already sitting at other tables. It can be written formally
as:

$CRP(C_{i}=c|\alpha,C_{1},C_{2},...,C_{i-1})=\begin{cases}\frac{N_{c}}{\alpha+N}\text{
if c = 1, 2, ..., K}\\\ \frac{\alpha}{\alpha+N}\text{ if c = K+1}\end{cases}$
(38)

where $c$ is the cluster chosen by the $n^{th}$ customer, $N_{k}$ is the
population of cluster $k$, $K$ is the number of already occupied tables and
$\alpha\in\mathcal{R}^{+}$ the concentration parameter. When the number of
clients goes to infinity, this process is equivalent to a draw from a
Dirichlet distribution over an infinite number of clusters with a uniform
concentration parameter $\alpha$. The form of Eq.38 is helpful to understand
the underlying dynamics of the process and the contribution of seminal works
we will detail now. It can be shown that the expected number of clusters after
$N$ observations evolves as $\log N$ [6].

The two best-known variations of the regular Dirichlet process that address
the “rich-get-richer” property control are the seminal Pitman-Yor process [93]
and the Uniform process [124]. Each of them can be expressed in a similar form
as Eq.38, and will be detailed in Section .17.

#####  .16.2.b Hawkes Process

Figure 24: A realization of a Hawkes process

Hawkes processes are typically used to model sequences of events in continuous
time, under the assumption that the probability that new events happen is
conditioned by the realization of earlier events. Typically, it can be used to
model retweet cascades, as the more retweets there are, the most likely other
retweets are to appear [24].

A Hawkes process is defined as a self-stimulating temporal point process.
Point processes are fully characterized by an intensity function $\lambda(t)$,
which is related to the probability $P(t_{events}\in[t;t+\Delta t])$ of an
event happening between $t$ and $t+\Delta t$ by ${\lambda(t)=\lim_{\Delta
t\rightarrow 0}\frac{P(t_{events}\in[t;t+\Delta t])}{\Delta t}}$. In the case
of Hawkes processes, $\lambda(t)$ is defined conditionally on all the events
that happened at times lower than $t$, that it the history up to $t$, noted
$\mathcal{H}_{<t}$. We provide a synthetic realization of a Hawkes process as
an illustration in Fig. 24. In our setup, we define one Hawkes process for
each cluster, independent from the others. The intensity of the Hawkes process
associated with cluster $c$ is defined as:

$\lambda_{c}(t|\mathcal{H}_{<t,c})=\sum_{\mathcal{H}_{<t,c}}\vec{\alpha_{c}}^{T}\cdot\vec{\kappa}(t_{i,c})$
(39)

where $t_{i,c}$ is the time of the $i^{th}$ observed in cluster $c$,
$\mathcal{H}_{<t,c}=\\{t_{i,c}|t_{i,c}<t\\}_{i=1,2,...}$ is the history of
events in cluster $c$ up to $t$, $\vec{\alpha_{c}}$ is a vector of
coefficients, $\vec{\kappa}(t)$ is a vector of kernel functions with the same
dimension as $\vec{\alpha}$ and $\cdot$ represents the dot product. The kernel
functions are defined at the beginning of the algorithm and are not modified
afterwards. We will later infer the weights vector $\vec{\alpha}$ to determine
which entries of the kernel vector are the most relevant for a given
situation. This technique has become standard in Hawkes processes modelling
and used in several occasions [28, 143].

The likelihood of a combination of $C$ independent Hawkes processes can be
written:

$\begin{split}\mathcal{L}(&\vec{\lambda}|\mathcal{H}_{<T,c})=\prod_{c\in
C}\mathcal{L}_{c}(\lambda_{c}|\mathcal{H}_{<T,c})\\\
&=\prod_{c}e^{-\int_{0}^{T}\lambda_{c}(t)dt}\prod_{t_{i,c}}\lambda_{c}(t_{i}|\mathcal{H}_{<t_{i,c},c})\\\
&=e^{-\sum_{c}\int_{0}^{T}\lambda_{c}(t|\mathcal{H}_{<t,c})dt}\prod_{t_{i,c^{\prime}},c^{\prime}=c}\lambda_{c}(t_{i,c^{\prime}}|\mathcal{H}_{<t_{i,c^{\prime}},c})\end{split}$
(40)

where $T$ is the upper time of the considered observation window, going from
$0$ to $T$. From now on, the dependence of Hawkes intensity functions on the
history will be implicit for clarity of presentation; we define
$\lambda(t):=\lambda(t|\mathcal{H}_{<t})$.

#####  .16.2.c Dirichlet-Hawkes Process – Expression

As the merging of Dirichlet processes and Hawkes processes, Dirichlet-Hawkes
processes aim to model self-stimulating clusters in continuous time. In their
definition of the DHP, the authors of [30] substitute the counts $N_{k}$ of
the DP (Eq. 38) with the inferred Hawkes intensities (Eq. 39), resulting in
the following form for the Dirichlet-Hawkes prior:

$P(C_{i}=c|t_{i},\lambda_{0},\mathcal{H}_{<t_{i},c})=\begin{cases}\frac{\lambda_{c}(t_{i})}{\lambda_{0}+\sum_{c^{\prime}}\lambda_{c^{\prime}}(t_{i})}\text{
if c$\leq$C}\\\
\frac{\lambda_{0}}{\lambda_{0}+\sum_{c^{\prime}}\lambda_{c^{\prime}}(t_{i})}\text{
if c=C+1}\end{cases}$ (41)

where $t_{i}$ is the arrival time of document $i$. The expression in Eq. 41 is
used as a prior that accounts for the publication dynamics in Bayesian
modelling. It is used as an a priori for a model that explicitly consider the
textual content of a document.

In Eq. 41, the authors consider a time-independent intensity function
$\lambda(t)=\lambda_{0}$. This process is used as the Dirichlet-Hawkes
equivalent of the concentration parameter $\alpha_{0}$ in a Dirichlet process
(see Eq. 63). It corresponds to a background Poisson process which determines
the rate at which new Hawkes processes are launched –that is, at which new
dynamic clusters are created.

Given the existence of the underlying Poisson process of parameter
$\lambda_{0}$, the temporal likelihood of the processes associated with all
clusters can be written:

$\begin{split}\mathcal{L}(&\vec{\lambda}|\mathcal{H}_{<T,c})=\mathcal{L}(\lambda_{0})\prod_{c}\mathcal{L}_{c}(\lambda_{c})\\\
&=e^{-\int_{0}^{T}\lambda_{0}dt}\prod_{c}e^{-\int_{0}^{T}\lambda_{c}(t)dt}\prod_{t_{i,c}}\lambda_{c}(t_{i})\\\
&=e^{-\lambda_{0}T-\sum_{c}\int_{0}^{T}\lambda_{c}(t)dt}\prod_{t_{i,c^{\prime}},c^{\prime}=c}\lambda_{c}(t_{i,c^{\prime}})\end{split}$
(42)

Note that $\mathcal{L}(\lambda_{0})=e^{-\int_{0}^{T}\lambda_{0}dt}$ because no
event is ever assigned to the Poisson process; the product over the history of
events runs over 0 events, which equals 1 by convention. We recall that the
Hawkes intensity dependence on the history of events is implicit;
$\lambda(t):=\lambda(t|\mathcal{H}_{<t})$.

#####  .16.2.d Textual modelling

We choose to model the textual content of documents as the result of a
Dirichlet-Multinomial distribution. This model is purposely simple to ease the
understanding, but can easily be replaced by a more complex one. A more
complete textual modelling is out of the scope of this presentation, which
focuses on the definition of the DHP prior. A document will be associated with
a given cluster according to word count in every cluster and words count in
the document only. The generative process is as follows:

$\theta_{i}\sim Dir(\theta_{0})\ \ \ \ ;\ \ \ \ \omega_{v,i}\sim
Mult(\theta_{i})$ (43)

where $\theta_{i}$ is the cluster of document $i$, and $\omega_{v,i}$ is the
$v^{th}$ word of document $i$. Let
$\mathcal{L}_{txt}(\vec{C}_{<i,c}|N_{<i,c},\theta_{0})$ be the marginal joint
distribution of every document’s cluster allocation up to the $i^{th}$ one.
The likelihood of the $i^{th}$ document belonging to cluster $c$ can then be
expressed as:

$\begin{split}\mathcal{L}(C_{i}=c|N_{<i,c},n_{i},\theta_{0})&=P(n_{i}|C_{i}=c,N_{<i,c},\theta_{0})\\\
&=\frac{\mathcal{L}_{txt}(\vec{C}_{<i,c}|N_{<i,c},\theta_{0})}{\mathcal{L}_{txt}(\vec{C}_{<i-1,c}|N_{<i,c},\theta_{0})}\\\
&=\frac{\frac{\cancel{\Gamma(\theta_{0})}}{\Gamma(N_{c}+n_{i}+\theta_{0})}\prod_{v}\frac{\Gamma(N_{c,v}+n_{i,v}+\theta_{0,v})}{\cancel{\Gamma(\theta_{0,v})}}}{\frac{\cancel{\Gamma(\theta_{0})}}{\Gamma(N_{c}+\theta_{0})}\prod_{v}\frac{\Gamma(N_{c,v}+\theta_{0,v})}{\cancel{\Gamma(\theta_{0,v})}}}\\\
&=\frac{\Gamma(N_{c}+\theta_{0})}{\Gamma(N_{c}+n_{i}+\theta_{0})}\prod_{v}\frac{\Gamma(N_{c,v}+n_{i,v}+\theta_{0,v})}{\Gamma(N_{c,v}+\theta_{0})}\end{split}$
(44)

where $N_{c}$ is the total number of words in cluster $c$ from observations
previous to $i$, $n_{i}$ is the total number of words in document $i$,
$N_{c,v}$ the count of word $v$ in cluster $c$, $n_{i,v}$ the count of word
$v$ in document $i$ and $\theta_{0}=\sum_{v}\theta_{0,v}$.

####  .16.3 Limits

A common feature of all the models we mentioned in this section is that they
use a non-parametric Dirichlet process (DP) prior (describe in Eq. 38) or
variations built on it, such as DHP (Eq. 41) and HDHP. Yet, on several
occasions, it has been pointed out that there are no specific reasons to use
this process in particular and that alternative forms might work better
depending on the dataset. In [128], the author relaxes several conditions
associated with DP and shows that alternative priors are an equally valid
choice in Bayesian modelling. In [124], the authors derive the Uniform process
(UP) and show that it performs better on a document clustering task. In
Section .17, we generalize UP and DP within a more general framework: the
Powered Dirichlet process (PDP). We show it performs better than DP on several
datasets. As we show in Section .17 and Section .18, considering alternative
definitions of the DP significantly leads to significantly different results.

As for DHP itself, it does not work well when the textual information within
documents conveys little information, that is when the text is short [141] or
when vocabularies overlap significantly. To answer this problem, the authors
develop an approach based on Dirichlet process mixtures, which is not designed
for continuous-time document streams – the temporal aspect comes from a
sampling function as in [2, 17]. There are other limiting cases for DHP, for
instance when temporal information conveys little information (few
observations, overlapping temporal intensities) or when documents within
textual clusters do not follow the same temporal dynamics. To overcome those
limitations, we develop the Powered Dirichlet-Hawkes process in the next
section.

In addition, we uncover in Section .18 new limiting cases in which DHP fails,
typically when publication times convey little information (overlapping Hawkes
intensities, few observations). We also show there are cases where different
documents generated from the same textual cluster do not follow the same
temporal dynamics (they are associated with a different $\lambda(t)$), which
the DHP is not designed to handle. For instance, an article published by a
popular newspaper is unlikely to have the same influence on subsequent similar
articles (temporal dynamics) as the same article published by a less popular
newspaper. Textual content is not perfectly correlated to publication times.

###  .17 Powered Dirichlet Process – Alleviate the “rich-get-richer”
assumption

####  .17.1 Introduction

The limits of DHP raised in Section .16.3 can be overcome by redefining the
Dirichlet Process (DP) it is built on. From a broader perspective, most
existing works based on Dirichlet Processes can be revisited by considering
alternative forms for Dirichlet-like priors. These priors have been used
extensively in Bayesian clustering over the last years. A non-exhaustive list
of application includes medicine, [50], natural language processing [19, 140],
genetics [104, 65, 85], recommender systems [4, 42], sociology [49, 27], etc.

The key idea of Bayesian clustering is to simulate a corpus of independent
observations by drawing them from a set of latent variables (clusters). Those
clusters are each associated with a probability distribution on the
observations, whose parameters are drawn from a prior distribution, as we will
formulate mathematically later. Now, an often desirable property of Bayesian
models is to make them nonparametric. In our case, it means that both the
number of clusters and their associated distributions are inferred. A very
popular prior on clusters distributions that allows this is the Dirichlet
process. It includes a chance for a new cluster to be created in the prior
probability of a distribution (often when an observation is not likely to be
explained by existing clusters). Otherwise, the observation is associated a
priori with an existing cluster with a probability proportional to that
cluster’s population. Note that the model the prior is associated with might
use this prior information, but might as well ignore it by design.

However, the Dirichlet process (and the related Pitman-Yor process) prior
comes with a strong hypothesis on the way observations are allocated to
various clusters: the rich-get-richer property [35]. A new observation a
priori belongs to a cluster with a probability proportional to the number of
observations already present in the cluster (see Eq.38); large clusters have a
greater chance to get associated with new observations. This modelling implies
a strong assumption on the way data is generated. It has already been pointed
out [128] that there is a need for more flexible priors.

####  .17.2 Motivation

The need for alternative priors is particularly relevant in the case of
imbalanced data and scale-dependent clustering. A cluster made of fewer
entities might go unnoticed due to a rich-get-richer prior. As an example of
the imbalance problem, consider a case where data is processed sequentially
–which is often the case when it comes to the Dirichlet process prior. The
first observation from a new cluster would then have a much larger a priori
probability to belong to a populated but irrelevant cluster, than to open a
new one (this probability decreases as $\frac{1}{N_{obs}}$). This typically
happens when sampling topics from news streams [124, 137]. In the case of
scale-dependent clustering, a similar problem arises. Consider clustering
people pinpointed on a map. Tiny clusters (at the scale of cities, for
instance) might go unnoticed at larger scales (countries, for instance). To
spot city clusters on a world map, the “rich-get-richer” assumption becomes
irrelevant and a “rich-get-no-richer” prior would be preferred [124]; the
optimal solution might as well be in-between these two priors, as in Fig.28.
We design a method to bridge the variety of possible priors between the
Dirichlet process and the Uniform process in a continuous fashion. By
generalizing existing works, our method shows there exist Dirichlet-based
priors that exhibit a yet unexplored class of behaviours, such as “rich-get-
less-richer”, “rich-get-more-richer” and “poor-get-richer”.

Little effort has been put into exploring alternative forms of priors for
nonparametric Bayesian modelling. In this section, we address this problem by
deriving a more general form of the Dirichlet process that explicitly controls
the importance of the “rich-get-richer” assumption. Explicitly, we derive the
Powered Chinese Restaurant Process (PCRP) that generalizes state-of-the-art
works such as UP [124] and DP. We show that controlling the “rich-get-richer”
prior of simple models yields better results on synthetic and real-world
datasets.

This work is motivated by the need to control the importance of the “rich-get-
richer” assumption in Dirichlet process (DP) priors. Developing such a more
permissive prior would impact many works based on the vanilla Dirichlet
Process. In our case, we are particularly interested in the implication of DHP
for the reasons discussed in Section .16.3.

The “rich-get-richer” property of the DP may not always be the most suitable
prior for modelling a given dataset. The usual motivation for using a DP prior
is that a new observation has a probability of being assigned to any cluster
proportional to its population (or intensity function, as in [30]) in the
absence of external information (such as inter-points distance in case of
spatial clustering, for instance). However, this assumption might be flawed in
several cases.

Typically, most state-of-the-art works rely on tuning a parameter $\alpha$
(see Eq.38) to get the “right” number of clusters (this parameter shifts the
distribution of the number of clusters as $\mathbb{E}(K|N)\propto\alpha\log N$
with $K$ the number of clusters and $N$ the number of observations). However,
we argue this is a bad practice. Imagine sampling topics in a news stream:
there is no specific reason for topics to appear at a rate $\alpha\log N$ as
in the regular DP. The logarithmic dependence on $N$ cannot be tuned using
$\alpha$. Such prior is then unfit to describe the data correctly as the
number of observations grows; worst, it can lead the model it is coupled with
on the wrong track. Moreover, when considering observations streams [124,
137], there is usually no specific a priori reason for a new observation to
belong to a cluster with a probability depending linearly on the cluster’s
size, as in the regular Dirichlet and Pitman-Yor processes (and variants). The
order of appearance then plays too important of a role. Typically, newer
observations might be associated with an existing large cluster, despite being
radically different from the points it comprises, due to a too influential
“rich-get-richer” assumption. To alleviate those assumptions, we develop a
more general form of the DP process allowing a natural control of the “rich-
get-richer” property.

####  .17.3 Background

#####  .17.3.a Previous works

######  Dirichlet process

The standard Dirichlet process has already been detailed in Section .16.2.a.
For completeness, we recall its expression as a Chinese Restaurant Process:

$CRP(C_{i}=c|\alpha,C_{1},C_{2},...,C_{i-1})=\begin{cases}\frac{N_{c}}{\alpha+N}\text{
if c = 1, 2, ..., K}\\\ \frac{\alpha}{\alpha+N}\text{ if c = K+1}\end{cases}$
(45)

######  Uniform process

A first process that breaks the “rich-get-richer” property is the Uniform
process. It has been used on some occasions [104, 65] without focus on the
prior itself. It has later been formalized and studied in comparison with the
regular Dirichlet and Pitman-Yor processes [124]. It can be written as
follows:

$UP(C_{i}=c|\alpha,C_{1},C_{2},...,C_{i-1})=\begin{cases}\frac{1}{\alpha+K}\text{
if c = 1, 2, ..., K}\\\ \frac{\alpha}{\alpha+K}\text{ if c = K+1}\end{cases}$
(46)

This formulation completely gets rid of the “rich-get-richer” property. The
probability of a new client joining an occupied table is a uniform
distribution over the number of occupied tables; it does not depend on the
tables’ population. In [124], it has been shown that the expected number of
tables evolves with $N$ as $\sqrt{N}$. Removing the “rich-get-richer” property
leads to a flat prior. As we show later, our formulation allows to retrieve
such flat priors and thus generalizes the Uniform Process.

######  Pitman-Yor process

Following the Chinese Restaurant process metaphor, the Pitman-Yor process [93,
62] proposed to incorporate a discount parameter when a client opens a new
table. Mathematically, the process can be formulated as:

$PY(C_{i}=c|\alpha,\beta,C_{1},C_{2},...,C_{i-1})=\begin{cases}\frac{N_{c}-\beta}{\alpha+N}\text{
if c = 1, 2, ..., K}\\\ \frac{\alpha+\beta K}{\alpha+N}\text{ if c =
K+1}\end{cases}$ (47)

The introduction of the parameter $\beta>0$ increases the probability of
creating new clusters. A table with a small number of customers has
significantly fewer chances to gain new ones, while the probability of opening
a new table increases significantly. It can be shown that the number of tables
evolves with the number of clients $N$ as $N^{\beta}$ [115, 124, 43]. However,
this process does not control the arguable “rich-get-richer” hypothesis [128],
since the relation to the population of a table remains linear; it only shifts
this dependence of a value $\beta$. It makes so by creating clusters based on
the number of existing clusters and the total number of observations, but not
according to the population of already existing clusters. Those play the same
role in the Pitman-Yor process as in the DP. The Pitman-Yor process thus comes
with two limitations. Firstly, since $\beta>0$, it cannot modify the process
to generate fewer clusters. Secondly, the discount parameter does not modify
the linear dependence on previous observations for cluster allocations — rich
still get richer; the prior is as peaky on large clusters as before. The
present section offers to address those two limitations.

######  Other extensions

Another similar prior, the Power-law Indian Buffet Process, has been proposed
so that a realization would yield a number of clusters obeying a power-law as
the number of observations increases [121]. This formulation can be seen as a
generalization of the Pitman-Yor process; it adds an additional parameter that
sums with $N$ in the denominator of Eq. 47. However, the posterior probability
for a new customer to belong to a cluster depends linearly on each cluster’s
size, and the “rich-get-richer” hypothesis is preserved.

Finally, the Generalized Gamma Process proposed a similar discount idea to
increase the probability of opening new clusters in [78]. The proposed prior
([78]-Eq.4) modifies a cluster’s probability to get chosen by subtracting a
constant term from each cluster’s population. Thus, the “rich-get-richer”
property is not alleviated in their approach either, since the dependence on
the cluster’s population is still linear. As for the PY process, this
formulation only allows to increase the number of clusters and does not
alleviate the “rich-get-richer” hypothesis.

#####  .17.3.b Contributions

In the next section, we derive the Powered Dirichlet Process (PDP) that allows
controlling the “rich-get-richer” property. The process is also referred to as
the Powered Chinese Restaurant Process (PCRP) due to its formulation being
close to the vanilla metaphor. This new process generalises state-of-the-art
works. Such generalization allows to define unexplored classes of a priori
hypotheses: poor-get-richer, rich-get-no-richer (Uniform process), rich-get-
less-richer, rich-get-richer (DP), and rich-get-more-richer. In doing so, we
define the Powered Dirichlet-Multinomial distribution. We detail some key
properties of the Powered Dirichlet Process (convergence, expected number of
clusters). Finally, we show that controlling the “rich-get-richer” prior of
simple models yields better results on synthetic and real-world datasets.

####  .17.4 The model

#####  .17.4.a The Dirichlet-Multinomial distribution

As explained earlier, the Dirichlet distribution yields a collection of
positive variables whose sum equals 1. The Multinomial distribution yields a
sum of counts in each of the $K$ “boxes” (here clusters) following $N$ draws
from an identical distribution over those boxes. We recall the definition of
the Dirichlet distribution and of the Multinomial distribution:

$Dir(\vec{p}|\vec{\alpha})=\frac{\prod_{k}p_{k}^{\alpha_{k}-1}}{B(\vec{\alpha})}\
\ \ \,\ \ \ \
Mult(\vec{N}|N,\vec{p})=\frac{\Gamma(\sum_{k}N_{k}+1)}{\prod_{k}\Gamma(N_{k}+1)}\prod_{k}p_{k}^{N_{k}}$
(48)

with $\vec{N}=(N_{1},N_{2},...,N_{K})$ where $N_{k}$ is the integer number of
draws assigned to cluster $k$, $N=\sum_{k}N_{k}$ the total number of draws,
$\Gamma(x)=(x-1)!$ is the gamma function, and
${B(\vec{x})=\prod_{k}\Gamma(x_{k})/\Gamma(\sum_{k}x_{k})}$ is the beta
function.

The Dirichlet-Multinomial distribution merges both distributions: the
probabilities for each “box” in the Multinomial distribution are sampled once
from a Dirichlet distribution. As we will show later, the Dirichlet process
can be derived from this distribution. The Dirichlet-Multinomial distribution
is defined as follows:

$\begin{split}p(\vec{N}|\vec{\alpha},n)&=\int_{\vec{p}}p(\vec{N}|\vec{p},n)p(\vec{p}|\vec{\alpha})d\vec{p}\\\
&=\frac{(n!)\Gamma(\sum_{k}\alpha_{k})}{\Gamma(n+\sum_{k}\alpha_{k})}\prod_{k=1}^{K}\frac{\Gamma(N_{k}+\alpha_{k})}{(N_{k}!)\Gamma(\alpha_{k})}\\\
\text{where }&\vec{p}\sim Dir(\vec{p}|\vec{\alpha})\,;\,\vec{N}\sim
Mult(\vec{N}|n,\vec{p})\end{split}$ (49)

In Eq.49, we sample $n$ values over a space of $K$ distinct clusters each with
probability $\vec{p}=(p_{1},p_{2},...,p_{K})$, using a Dirichlet prior with
parameter $\vec{\alpha}=(\alpha_{1},\alpha_{2},...,\alpha_{K})$. To derive the
Dirichlet process equation, we must compute a new observation’s conditional
distribution to belong to any cluster given the allocation of all the previous
random variables when $K\rightarrow\infty$.

#####  .17.4.b Powered conditional Dirichlet prior