# Systematic errors in estimates of $R_{t}$ from symptomatic cases in the
presence of observation bias

Guido Sanguinetti1,2 1SISSA, Trieste, Italy
2School of Informatics, University of Edinburgh, UK
<EMAIL_ADDRESS>
(23 November 2020)

###### Abstract

We consider the problem of estimating the reproduction number $R_{t}$ of an
epidemic for populations where the probability of detection of cases depends
on a known covariate. We argue that in such cases the normal empirical
estimator can fail when the prevalence of cases among groups changes with
time. We propose a Bayesian strategy to resolve the problem, as well as a
simple solution in the case of large number of cases. We illustrate the issue
and its solution on a simple yet realistic simulation study, and discuss the
general relevance of the issue to the current covid19 pandemic.

##  1 Introduction

The instantaneous reproduction number $R_{t}$ is an effective empirical
measurement of the velocity with which an epidemic is spreading through a
population. A value of $R_{t}$ greater than 1 denotes an epidemic which is
growing exponentially at that time, while values smaller than 1 are witnessed
in declining phases of the epidemic. Because of its immediate interpretation
and its independence of detailed modelling assumptions, estimate of the value
of $R_{t}$ are of immense importance to public health experts and policy
makers, and are normally an essential consideration in determining measures to
fight the spread of an epidemic.

$R_{t}$ is defined as the expected number of secondary infections at time $t$
from each infected case, or equivalently the ratio of the number of new
infected cases at time $t$ to the average infectiousness of individuals at
that time (see Section 2 for a more precise definition). As such, it returns a
local (in time), model-free exponential approximation to the epidemic
dynamics, providing precious information to population health scientists and
policy makers. While alternative approximations which do not assume
exponential dynamics exist [2], in practice most epidemiological studies
utilise the exponential approximation approach, as witnessed in the manifold
studies of the dynamics of the recent covid19 pandemic, e.g. [6, 7]. Because
the ascertainment of dates of infection is in general difficult if not
impossible, in practice $R_{t}$ is estimated on the basis of new cases whose
symptom onset happens at time $t$, as suggested in several publications, see
e.g. [3]. Therefore, estimates of $R_{t}$ are describing the dynamics of the
epidemic in the subset of the infected population which develops sufficiently
severe symptoms as to be detected as infected. Such a restriction is often
considered as a positive feature, as it avoids the unavoidable uncertainty
linked to detection of asymptomatic and paucisymptomatic cases, which is
heavily confounded by testing strategies which may be highly time-varying.

Mathematically, computing $R_{t}$ from symptomatic cases only is fully
justified under the assumption that new symptomatic cases constitute a
constant fraction $\alpha\leq 1$ of the total (unknown) number of new cases in
a day. Such an assumption guarantees that the overall dynamics of the epidemic
among the symptomatic cases are identical to the true epidemic, just rescaled
by a factor $\alpha$. However, empirical evidence from the covid19 pandemic in
western Europe suggests that this assumption is not always justified.

It is well known that infection by the new coronavirus results in a large
fraction of asymptomatic cases, however such fraction is heavily dependent on
a number of additional covariates. In particular, it is well known that
younger people are far less likely to develop symptoms than older people [4],
therefore if the progression of the disease is different between different age
groups, it is plausible that the ratio $\alpha$ of symptomatic cases to total
cases will not remain constant, but it will follow the spread of the epidemic
across the age groups. Such a situation has been widely observed during the
autumnal second wave of the infection, when the virus spread rapidly among
younger people in holiday resorts in the summer, and gradually started
spreading to older people during the autumn, with a gradual increase of the
average age of the cases detected witnessed.

In this brief paper, I show that estimating $R_{t}$ from symptomatic cases
only can lead to mistaken conclusions when the probability of developing
symptoms depends on a covariate which is not constant in time. I also show how
the problem can be fixed by introducing latent variables, and how an extremely
simple correction can provide very accurate maximum likelihood estimates of
the true $R_{t}$ when the number of cases is large and the dependence of the
probability of developing symptoms on the covariate is known. Numerical
examples on a simple yet informative model confirm the potential relevance of
these considerations.

##  2 Mathematical formulation of the problem and proposed solutions

###  2.1 Definitions and problem statement

The instantaneous reproduction number $R_{t}$ is defined mathematically by the
equation

$I_{t}=R_{t}\sum_{\tau=1}^{t}p(\tau)I_{t-\tau}$ (1)

where $I_{t}$ is the number of infected individuals at time $t$, and $p(\tau)$
is the distribution over generation times, effectively accounting for
incubation times. In practice, the distribution over generation times is
assumed known from epidemiological studies; for example, the Italian Istituto
Superiore di Sanitá (ISS) uses estimates of generation times distribution
obtained from the early Covid19 outbreak in Lombardy [1] to compute $R_{t}$
estimates in its weekly monitoring. Equation (1) defines an instantaneous rate
of change of the infected population independently of the stochastic
generative process underpinning the epidemic. Practically one simply combines
equation (1) with a Poisson or Gaussian noise distribution on the output
$I_{t}$ to provide straightforward Bayesian or maximum likelihood estimates of
the parameter of interest $R_{t}$.

Equation (1) is linear in the number of new daily infections $I_{t}$,
therefore it still holds when the vector of new infections is rescaled by a
constant. Therefore, under the assumption that symptomatic cases are a
constant fraction of new cases, estimates of $R_{t}$ obtained using only
symptomatic cases will be identical to estimates obtained using all new
infections. Such an assumption will however fail if the fraction of
symptomatic cases depends on time.

The Covid19 pandemic in Western Europe, during its second wave in the autumn
of 2020, has shown a constant gradual increase in the median age of cases. As
an example, the Italian ISS reported a median age of Covid19 positives of 29
in the period 17-23 August 2020, which increased to 47 in the period from mid
October to mid November 2020111Weekly data on median age of cases appear to be
unavailable after September. Data from ISS reports
https://www.epicentro.iss.it/coronavirus/sars-cov-2-dashboard.. The fraction
of Covid19 cases presenting clinically relevant symptoms is estimated to vary
from less than 20% in young individuals to over 70% in elderly people [4].
Therefore, the symptomatics rate must have changed in the period August-
November, invalidating the assumptions underlying the procedure to estimate
$R_{t}$. In particular, part of the increase in symptomatic cases will be the
result of an increased rate of symptomatic cases in older individuals, leading
therefore to an inflation in the estimates of $R_{t}$.

###  2.2 Latent variable model

A natural solution to the problem is to retain the original definition of
$R_{t}$ in terms of infected numbers, and to introduce observable variables
$S_{t}(l)$ denoting number of symptomatic cases at time $t$ with covariate
value equal to $l$222This treatment is equally suited to handling discrete or
continuous covariates, however for practical convenience it might be easier to
group continuous covariates such as age into discrete classes.. The
observables are obtained from the latent number of infected cases with
covariate $l$ $I_{t}(l)$ through a Binomial observation model

$S_{t}(l)\sim\mathrm{Binom}\left(I_{t}(l),\pi_{l}\right)$ (2)

where $\pi_{l}$ is the symptomatic rate for infected cases in group $l$.
Assuming this rate to be known (for example from epidemiological studies such
as [4]), then the variables $I_{t}(l)$ are all independent a posteriori333This
is because the renewal equations (1) and (3) are deterministic conditioned on
$I_{t}$ and $I_{t}(l)$ respectively. and, given a suitable prior, they can be
estimated independently via e.g. the Metropolis-Hastings (M-H) algorithm or
any other Bayesian sampler444I am not aware of a conjugate prior distribution
over the number of trials to a binomial likelihood where the success rate is
known.. To obtain a posterior distribution over the parameter of interest
$R_{t}$, one simply needs to rewrite equation (1) by summing over the
covariate values

$\sum_{l}I_{t}(l)=R_{t}\sum_{\tau,l}\pi_{l}p(\tau)I_{t-tau}(l)$ (3)

and plugging samples of $I_{t}(l)$ in equation (3)555We assume that generation
times are independt of the covariate $l$..

If the symptomatic rate $\pi_{l}$ is not known, then it can also be estimated
within a Bayesian framework. An efficient solution could be to assign each
$\pi_{l}$ independent Beta priors and run a block Metropolis-within-Gibbs
scheme: conditioned on observables and $\pi_{l}$ values, the $I_{t}(l)$ can be
sampled in parallel via independent M-H samplers. Given sampled values of
$I_{t}(l)$, the posterior on each $\pi_{l}$ value can be obtained analytically
and sample values of $\pi_{l}$ could be drawn in a straightforward way. It is
however probable that informative priors on $\pi_{l}$ would be needed for the
estimates to have an acceptable level of uncertainty.

###  2.3 Gaussian approximation and analytical corrections

Assuming that the numbers of infected individuals for each group $I_{t}(l)$
are sufficiently large, and that the symptomatic rates $\pi_{l}$ are also not
too close to zero or one, it is then possible to approximate the binomial
distribution with a Gaussian conditional on $I_{t}(l)$, so that

$p\left(S_{t}(l)|I_{t}(l),\pi_{l}\right)\simeq\mathcal{N}\left(\pi_{l}I_{t}(l),\pi_{l}(1-\pi_{l})I_{t}(l)\right).$
(4)

Introducing $\delta_{t}(l)=I_{t}(l)-\frac{S_{t}(l)}{\pi_{l}}$, we can rewrite
the exponent in (4) as

$-\frac{\left(S_{t}(l)-\pi_{l}I_{t}(l)\right)^{2}}{\pi_{l}(1-\pi_{l})I_{t}(l)}=-\frac{\pi_{l}}{1-\pi_{l}}\frac{\delta_{t}(l)^{2}}{\frac{S_{t}(l)}{\pi_{l}}\left(1+\frac{\pi_{l}\delta}{s_{t}(l)}\right)}$

from which it is clear that the maximum likelihood value for $I_{t}(l)$ is
obtained for $\delta=0$, i.e. when $I_{t}(l)=\frac{S_{t}(l)}{\pi_{l}}$.
Plugging this estimator in the equation for $R_{t}$ (3) yields the following
estimator

$R_{t}^{MLE}=\frac{\sum_{l}S_{t}(l)\pi_{l}^{-1}}{\sum_{\tau=1}^{t}\sum_{l}S_{t-\tau}(l)\pi_{l}^{-1}p(\tau)}$
(5)

##  3 Numerical illustration

To illustrate the impact of violating the assumption of constant probability
of symptoms onset, we analyse a simple epidemic model which simulates the
spread of a disease among two equal sized populations with weak interactions.
The two populations have equal size of 200.000 individuals, and differ only in
the probability that an infected individual will develop symptoms: infected
individuals in one population, Y, will develop symptoms with probability 0.3,
while individuals in the other, O, will develop symptoms with probability 0.8.
These values were chosen to mimic the values reported for COVID19 in young and
old individuals [4].

Infection dynamics are almost identical in the two sub-populations, according
to a stochastic SI model with $R_{0}=1.4$. To simplify further, we assumed
that infected individuals remain infectious for only 1 time step, i.e. there
is no incubation period and can only infect individuals the day after they
were infected. Individuals in the Y population can also infect with low
probability individuals in the O population. The resulting epidemic dynamics
over 100 simulation steps are shown in Figure 1 top left, resulting in two
nearly identical waves of infections with a delay of approximately 20 time
units. The observed dynamics of symptomatic cases is however very different,
with a first low peak of symptomatic Y cases and a later, much larger peak of
O cases (Fig. 1 top right). The probability of showing symptoms transitions
sharply at around time 30 from the Y probability (0.3) to the O probability
(0.8), reflecting the change in prevalence of the infection among the two
populations (Fig 1 bottom left). Correspondingly, around time 30 we see a
significant deviation of the $R_{t}$ estimate computed on symptomatics only
(Fig 1 bottom right blue line) from the ground truth estimate computed on all
cases (red line). The Gaussian MLE correction proposed in equation (5) is in
almost perfect agreement with the true $R_{t}$ (cyan line).

Figure 1: Example results in a simple SI epidemic model with two populations
with different probabilities of symptoms onset (young, $p(s_{y}|I_{y})=0.3$,
and old $p(s_{o}|I_{o})=0.8$): top left, new daily infections for young (blue)
and old (red) populations. An early wave of infection in the young population
slowly trickles to the old population, which develops a nearly identical wave
approximately 20 days later. Top right: observed new daily symptomatic cases,
showing an apparently much larger infection among the older. Bottom left,
fraction of new symptomatic cases over total new cases, showing a clear rapid
increment when the infection among the old population starts to dominate.
Bottom right, estimated $R_{t}$ using symptomatic cases (blue), total cases
(red) and corrected maximum likelihood estimates (cyan), in the time interval
when the fraction of symptomatic cases changes. Standard estimates based on
symptomatic cases only overestimate $R_{t}$ by approximately 10%, while the
corrected maximum likelihood estimate is nearly identical to the true $R_{t}$
value consistently throughout the period.

While the analysis of this simple model reveals a marked deviation of the
estimated $R_{t}$ from the true $R_{t}$, the effect appears relatively modest,
approximately around 10% relative error. One possible explanation for this is
that the simple model does not involve incubation times, so that the rate of
change of the fraction of symptomatics is slow relative to the fast contagion
time.

I therefore explored a different scenario where the generation time
distribution is zero everywhere except for 0.5 at days -3 and -4 (i.e., there
is an incubation time of 3 days and infectiousness lasts two days). I set the
$R_{0}$ to 2.4 to yield similar epidemic dynamics to the simple case in Figure
1. The deviation of the estimates from symptomatic cases from the true values
appear more marked, presumably because the nonzero incubation time leads to a
sharper apparent increment in symptomatic cases due to the shift in
population. Once again, the simple correction from the Gaussian approximation
seems to rescue the problem to a large extent.

Quantifying the errors in estimates across 100 independent simulations, we see
that for both the simple model with no incubation and the more complex model
the Gaussian approximation yields a very significant improvement in estimates
of the true $R_{t}$, with error reduction ranging from a factor 7 in the
simple model to over an order of magnitude for the more complex model.

Figure 2: Errors in estimates of $R_{t}$ over 100 independent simulations:
left boxplot of differences in estimates from symptomatics only (left box) and
corrected by Gaussian approximation (right box) in model without incubation
time; right similar as left panel but on model with incubation time

##  4 Discussion

In this short paper, I have identified and discussed a possible shortcoming in
the current method for estimating the reproduction number $R_{t}$ of an
epidemic. While estimating $R_{t}$ from symptomatic cases only is certainly a
practical solution to avoid the difficulties in estimating total infection
numbers in epidemics with many asymptomatic individuals, I argue both
theoretically and numerically that this strategy can yield erroneous results
when the proportion of symptomatic cases is not constant in time. I also
showed that the problem can be rescued to a large extent when the change of
symptomatic cases is the result of a shift of the epidemic across populations
with different (known) characteristics.

A major limitation of this work is the assumption that the time dependence of
the fraction of symptomatics results only from the population structure and
how the epidemic spreads across the population. In reality, myriad effects
might lead to even larger errors, from local difficulties in recording and
reporting symptomatic cases due to health services being overwhelmed, to
changes in testing policies (for example leading to the identification of more
or less individuals with light symptoms). Nevertheless, I would argue that,
despite its limitations, the proposed analysis has practical merit, since
changes in median age of infected individuals have been recorded widely in the
current COVID19 pandemic, leading to errors that could easily be corrected. It
is worth pointing out that even fractional errors in estimates of $R_{t}$ can
have huge practical consequences, since the reproduction number is the main
quantitative indicator adopted by governments in deciding containment measures
which have frequently enormous economic impact. This issue is even more
pressing because frequently estimates of $R_{t}$ at local level are adopted,
where factors such as infections spreading in care homes (which have by
definition a high median age and hence a large fraction of symptomatics) can
completely dominate the recorded numbers of symptomatic cases.

Methodologically, this work points to a more substantial usage of latent
variable models in epidemiological studies. Latent variable models are a
popular and well established area of research in statistics and machine
learning; to my knowledge, the proposal of using latent variables to estimate
$R_{t}$ in the current covid19 pandemic has been explored in [5], where
however the focus was reconstructing total numbers of infected individuals
from total numbers of positive tests. It is widely believed that machine
learning methods might greatly help monitoring and predicting the spread of
pandemics, however their application is critically limited by data
availability. Regrettably, relevant real data to test models is scarce: in the
case of this work, while median age of cases and total recorded cases are
widely available, I have not been able to obtain real data on symptomatic
cases per day and their age categories, therefore preventing a more meaningful
evaluation of the proposed methodology.

## References

  * [1] D Cereda, M Tirani, F Rovida, V Demicheli, M Ajelli, P Poletti, F Trentini, G Guzzetta, V Marziano, A Barone, M Magoni, S Deandrea, G Diurno, M Lombardo, M Faccini, A Pan, R Bruno, E Pariani, G Grasselli, A Piatti, M Gramegna, F Baldanti, A Melegaro, and S Merler.  The early phase of the COVID-19 outbreak in Lombardy, Italy, 2020.
  * [2] Gerardo Chowell, Cécile Viboud, Lone Simonsen, and Seyed M Moghadas.  Characterizing the reproduction number of epidemics with early subexponential growth dynamics.  Journal of The Royal Society Interface, 13(123):20160659, 2016.
  * [3] Anne Cori, Neil M Ferguson, Christophe Fraser, and Simon Cauchemez.  A new framework and software to estimate time-varying reproduction numbers during epidemics.  American journal of epidemiology, 178(9):1505–1512, 2013.
  * [4] Nicholas G Davies, Petra Klepac, Yang Liu, Kiesha Prem, Mark Jit, Rosalind M Eggo, CMMID COVID-19 working group, et al.  Age-dependent effects in the transmission and control of covid-19 epidemics.  Nature Medicine, 26:1205–1211, 2020.
  * [5] Andrea de Simone and Marco Piangerelli.  A bayesian approach for monitoring epidemics in presence of undetected cases.  Chaos Solitons & Fractals, 08 2020.
  * [6] Seth Flaxman, Swapnil Mishra, Axel Gandy, H Unwin, Helen Coupland, T Mellan, Harisson Zhu, Tresnia Berah, J Eaton, P Perez Guzman, et al.  Report 13: Estimating the number of infections and the impact of non-pharmaceutical interventions on covid-19 in 11 european countries.  2020\.
  * [7] Kin On Kwok, Florence Lai, Wan In Wei, Samuel Yeung Shan Wong, and Julian WT Tang.  Herd immunity–estimating the level required to halt the covid-19 epidemics in affected countries.  Journal of Infection, 80(6):e32–e33, 2020.


# A Fast Algorithm for Ranking Users by their Influence in Online Social
Platforms

Nouamane Arhachoui, Esteban Bautista, Maximilien Danisch, Anastasios
Giovanidis  Sorbonne Université, CNRS – LIP6, F-75005 Paris, France
{nouamane.arhachoui, esteban.bautista-ruiz<EMAIL_ADDRESS>
###### Abstract

Measuring the influence of users in social networks is key for numerous
applications. A recently proposed influence metric, coined as $\psi$-score,
allows to go beyond traditional centrality metrics, which only assess
structural graph importance, by further incorporating the rich information
provided by the posting and re-posting activity of users. The $\psi$-score is
shown in fact to generalize PageRank for non-homogeneous node activity.
Despite its significance, it scales poorly to large datasets; for a network of
$N$ users, it requires to solve $N$ linear systems of equations of size $N$.
To address this problem, this work introduces a novel scalable algorithm for
the fast approximation of $\psi$-score, named Power-$\psi$. The proposed
algorithm is based on a novel equation indicating that it suffices to solve
one system of equations of size $N$ to compute the $\psi$-score. Then, our
algorithm exploits the fact that such a system can be recursively and
distributedly approximated to any desired error. This permits the
$\psi$-score, summarizing both structural and behavioral information for the
nodes, to run as fast as PageRank. We validate the effectiveness of the
proposed algorithm, which we release as an open source Python library, on
several real-world datasets.

###### Index Terms:

Social Networks, Influence, Algorithms, Scalability, PageRank.

††publicationid: pubid: IEEE/ACM ASONAM 2022, November 10-13, 2022
978-1-6654-5661-6/22/$31.00 © 2022 IEEE

##  I Introduction

###  I-A Context

Online Social Platforms (OSPs) have become ubiquitous in society. Their
proliferation, as well as their massive number of users, has made the
quantification of the influence of users in OSPs a task of utmost importance.
For example, spotting the so-called influencers of a network is crucial for
social scientists to better understand the dynamics that determine new trends
[1] or the roots of political polarization [2]. It is also key for companies
in order to develop a better marketing of their products [3]. Machine learning
algorithms are another example, where identifying a reduced set of users whose
features are present over the entire network is essential to reduce the number
of parameters and to avoid the curse of dimensionality [4, 5].

Social platforms usually provide a user with a wall and a news-feed. In the
wall, users can post messages. In the news-feed, users see messages that other
users, called their leaders (or followees), have posted on their wall. A user
can then re-post a message from its news-feed into its wall, allowing posts to
diffuse through the network [6]. A fundamental problem is then to quantify the
influence of users. Namely, the degree to which their posts spread in the
network and are visible to others.

In order to address the aforementioned problem, classical centrality metrics
[7] permit to assign a ranking to the nodes of the network according to their
structural importance. However, such metrics are unsatisfactory to assess
influence in OSPs, as they do not consider the posting and sharing activity of
users in their scores. The relevance of this information is highlighted in
[8], where it is shown that users with the highest number of followers (in-
degree) may not necessarily be the ones that generate the highest number of
retweets or impressions. This means that the in-degree, which is correlated
with the popularity of a user, is not the only factor that determines
influence. Similarly, the authors in [9] address the question “Are social
links valid indicators of real user interaction?” and analyze a large user
trace from Facebook to showcase the difference. As a result, it is necessary
to go beyond pure structural metrics.

To overcome this limitation, a new metric, coined as $\psi$-score, has been
recently introduced in [10].

The $\psi$-score results from an Online Social Platform model, also presented
in [10], in which the diffusion of information is taken into account with
activity features for each user. In the general case where users post and re-
post content at a different frequency from each other, the OSP model
accurately incorporates this user activity to the structural information
related to the graph, thus providing a more meaningful ranking metric of user
influence in OSPs. This metric was shown in [10, Theorem 5] to equal PageRank
in the special case of homogeneous activity, where all users create the same
amount of posts in a given time interval and share with the same rate within
this time interval.

###  I-B Goal and contributions

The problem we address in this paper is that the $\psi$-score calculation from
[10] does not scale to large graphs since it needs the resolution of $N$
linear systems of $N$ equations each, $N$ being the number of users in the
network, whereas the PageRank would only require to solve one system with $N$
unknowns.

We thus propose a new algorithm based on an equation that revamps these $N$
system into a single one, still with $N$ equations. Our algorithm, which
accepts a recursive and distributed implementation, allows to approximate the
$\psi$-score of the OSP with a similar speed as the traditional approximation
for PageRank via the power method. Hence, the $\psi$-score can be used as a
more expressive alternative to PageRank in the task of ranking users by their
influence.

Finally, we provide an open source software
library111https://github.com/NouamaneA/psi-score in order to allow anyone to
implement the $\psi$-score in their own projects.

###  I-C Related Works

Centrality measures have been systematically used in various fields of
application of network science such as social networks [11], transportation
networks [12], communication networks [13], biological networks [14] or even
political networks [15]. The most used metrics are degree, betweenness,
closeness, eigenvector centrality, and PageRank, among others [7, 16].

These metrics, which are practically used to rank nodes of a network, only
consider the network topology. For example, PageRank [16] is derived through a
random walk with teleportation on the studied graph, taking at each step
uniform decisions about which node to consider next.

There are several existing algorithms such as Push [17] and the use of
Chebyshev polynomials [18] that exploit the random walk interpretation of
PageRank in order to accelerate its calculation. The Push method [17] loops on
the nodes by allowing each of them to update the approximation of its own
PageRank value and updates its neighbors’ values. The method proposed in [18]
uses the Chebyshev polynomials to approximate the PageRank vector in a more
efficient way compared to the known power-method introduced with the metric
[16]. Although the literature to accelerate graph eigenvalue-based measures is
rich, in this paper we will need to focus in what follows on the specific
structure of the $\psi$-score with its current algebraic limitations in order
to propose an original way to make it scale as fast as PageRank does. Having
established in this paper a power-iteration inspired algorithm that scales,
this opens future perspectives to study other methods inspired by the Push
method [17] or Chebyshev polynomials [18] to achieve even better performance.

###  I-D Outline of the paper

The paper is structured as follows. Section II presents the recently proposed
$\psi$-score and discusses its relation to classical centrality metrics.
Section III contains our main contribution: a new algorithm for the fast
approximation of the $\psi$-score that scales to massive datasets. We show
that the algorithm runtime is comparable to PageRank, thus providing a more-
diverse alternative. Section IV briefly presents our software contribution: an
open source Python library for the $\psi$-score computation. Section V
evaluates the performance of our new algorithm on several real-world networks.
Conclusions and future work are discussed in Section VI.

TABLE I: Notations used in the paper Notation | Description
---|---
$\psi_{i}$ |  $\psi$-score of user $i$
$\bm{\psi}$ |  $\psi$-score column-vector
$p_{i}^{(n)}$ | impressions of $i$ on the Newsfeed of $n$
$q_{i}^{(n)}$ | (influence) impressions of $i$ on the Wall of $n$
$\mathbf{p}_{i}$ | vector with every $p_{i}^{(n)}$
$\mathbf{q}_{i}$ | vector with every $q_{i}^{(n)}$
$\mathbf{P}$ | matrix with $\mathbf{p}_{i}$ as column $i$
$\mathbf{Q}$ | matrix with $\mathbf{q}_{i}$ as column $i$
$\mathbf{b}_{i}$ | input vector of normalised posting activity
$\mathbf{d}_{i}$ | vector with the posting ratio of a user $i$
$\mathbf{B}$ | matrix with $\mathbf{b}_{i}$ as column $i$
$\mathbf{D}$ | diagonal matrix with $\mathbf{d}_{i}$ as column $i$
$\mathbf{1}$ | column-vector of $\mathbb{R}^{N}$ full of ones, i.e. $(1\;1\;\cdots\;1)^{T}$

##  II The $\psi$-score

Figure 1: A small example of the Online Social Platform from the point of view
of user $i$.

Recently, a new score of users’ influence over their social network (platform)
has been introduced in [10]. The calculation of this score involves
information not only about the users’ positions inside the graph structure,
but additionally their posting and sharing activity. These are two extra
important features that enrich previous graph-based metrics: when a user posts
more frequently he/she should have a higher influence because his/her posts
appear more often on the newsfeed of others; furthermore, when a user shares
content rarely, very few posts (other than of his/her own creation) pass
through to reach his/her followers.

More formally, let $\mathcal{G}=(\mathcal{N},\mathcal{E})$ be an unweighted
and directed graph where $\mathcal{N}$ denotes the set of vertices of
cardinality $|\mathcal{N}|=N$ and $\mathcal{E}$ refers to the set of edges of
cardinality $|\mathcal{E}|=M$. Vertices model the users of social platforms
and edges the follower-leader relations: $(i,j)\in\mathcal{E}$ encodes that
user $i$ follows user $j$. Let $\lambda$ be the vector of posting activity
over the set of vertices. That is, $\lambda^{(n)}$ is the number of posts per
unit of time (posting frequency) that user $n$ creates on his/her wall.
Similarly, let $\mu$ be the vector of re-posting/sharing activity over the set
of vertices. Then, $\mu^{(n)}$ refers to the frequency with which user $n$
visits his/her news-feed and chooses one of the current entries uniformly at
random to re-post it on his/her wall.

From now on in this section we will focus on posts of origin $i$. The authors
of [10] show that the aforementioned graph and activity features can be used
to compute the following two vectors associated with the influence of user
$i\in\mathcal{N}$:

  * •

$\mathbf{p}_{i}=(p_{i}^{(1)}\;p_{i}^{(2)}\;\cdots\;p_{i}^{(N)})^{T}$ where,
for $n\in\mathcal{N}$, $p_{i}^{(n)}$ is the expected percentage of posts
originating from user $i$ on the news-feed of user $n$.

  * •

$\mathbf{q}_{i}=(q_{i}^{(1)}\;q_{i}^{(2)}\;\cdots\;q_{i}^{(N)})^{T}$ where,
for all $n\in\mathcal{N}$, $q_{i}^{(n)}$ is the expected percentage of posts
originating from user $i$ on the wall of user $n$.

A central point in [10] is that $q_{i}^{(n)}$ can be interpreted as the
influence of user $i$ on user $n$, meaning that the more posts of origin $i$
appear on the wall of $n$, the larger the influence that $i$ has on $n$. This
is used to define the influence of user $i$ over the entire network as the
average influence of $i$ on all users. This notion is formalized by the so-
called $\psi$-score of user $i$ which is defined as

$\psi_{i}=\frac{1}{N}\sum\limits_{n\in\mathcal{N}}{q_{i}^{(n)}}.$ (1)

The expression in (1) involves all entries of the vector $\mathbf{q}_{i}$ for
the calculation of the influence score $\psi_{i}$ of user $i$. However,
$\mathbf{q}_{i}$ is not a given information and to find $q_{i}^{(n)}\;\forall
n\in\mathcal{N}$ we need to resolve a system of equations. Specifically, each
$q_{i}^{(n)}$ can be derived from the $p_{i}^{(n)}$ thanks to the following
equations:

$\displaystyle\text{for }n=j\neq i,\;(\lambda^{(j)}+\mu^{(j)})q_{i}^{(j)}$
$\displaystyle=$ $\displaystyle\mu^{(j)}p_{i}^{(j)}$ (2)
$\displaystyle\text{for }n=i,\;(\lambda^{(i)}+\mu^{(i)})q_{i}^{(i)}$
$\displaystyle=$ $\displaystyle\lambda^{(i)}+\mu^{(i)}p_{i}^{(i)}$ (3)

Moreover, the vector $\mathbf{p}_{i}$ is the solution of a fixed-point problem
(see [10], section III). For a post of origin $i$ in the Newsfeed of $j$ we
have the following balance equation:

$\sum\limits_{k\in\mathcal{L}^{(j)}}(\lambda^{(k)}+\mu^{(k)})p_{i}^{(j)}=\lambda^{(i)}\mathds{1}_{\\{i\in\mathcal{L}^{(j)}\\}}+\sum\limits_{k\in\mathcal{L}^{(j)}}\mu^{(k)}p_{i}^{(k)}$
(4)

Since there are no self-loops and a user cannot be its own leader or follower,
we have the following balance equation for user $i$’s Newsfeed:

$\sum\limits_{k\in\mathcal{L}^{(i)}}(\lambda^{(k)}+\mu^{(k)})p_{i}^{(i)}=\sum\limits_{k\in\mathcal{L}^{(i)}}\mu^{(k)}p_{i}^{(k)}$
(5)

The system of equations (2)-(3) and (4)-(5) can be written in matrix form with
the following linear system:

$\displaystyle\mathbf{p}_{i}$ $\displaystyle=$
$\displaystyle\mathbf{A}.\mathbf{p}_{i}+\mathbf{b}_{i}$ (6)
$\displaystyle\mathbf{q}_{i}$ $\displaystyle=$
$\displaystyle\mathbf{C}.\mathbf{p}_{i}+\mathbf{d}_{i}$ (7)

where,
$\mathbf{A}\in\mathbb{R}_{+}^{N\times
N}:\;a_{ji}=\frac{\mu^{(i)}}{\sum\limits_{\ell\in\mathcal{L}^{(j)}}(\lambda^{(\ell)}+\mu^{(\ell)})}\mathds{1}_{\\{i\in\mathcal{L}^{(j)}\\}},\\\
\mathbf{b}_{i}\in\mathbb{R}_{+}^{N}:\;b_{ji}=\frac{\lambda^{(i)}}{\sum\limits_{\ell\in\mathcal{L}^{(j)}}(\lambda^{(\ell)}+\mu^{(\ell)})}\mathds{1}_{\\{i\in\mathcal{L}^{(j)}\\}},\\\
\mathbf{C}\in\mathbb{R}_{+}^{N\times
N}:\;c_{ji}=\frac{\mu^{(j)}}{\lambda^{(j)}+\mu^{(j)}}\mathds{1}_{\\{j=i\\}},\\\
\mathbf{d}_{i}\in\mathbb{R}_{+}^{N}:\;d_{ji}=\frac{\lambda^{(i)}}{\lambda^{(i)}+\mu^{(i)}}\mathds{1}_{\\{j=i\\}}.\\\
$

To solve the fixed-point problem of (6), we can use the following iteration
(see [10], Theorem 4):

$\mathbf{p}_{i}(t)=\mathbf{A}\mathbf{p}_{i}(t-1)+\mathbf{b}_{i}$ (8)

where $t$ is the current iteration. With any initialization of
$\mathbf{p}_{i}(0)$, the $\mathbf{p}_{i}(t)$ converges towards the solution
$\mathbf{p}_{i}$ of (6) as $t\rightarrow\infty$ because $\mathbf{A}$ is a sub-
stochastic matrix. This method is less costly than using matrix inversion to
solve the system.

input : Origin User $i$, $N$ number of users, $N\times N$ matrix $\mathbf{A}$,
vector $\mathbf{b}_{i}$, $\mathbf{p}$-tolerance $\varepsilon$

output : vector $\mathbf{p}_{i}$

$\mathbf{p}_{i}\leftarrow\mathbf{b}_{i}$;

$t\leftarrow 0$;

$gap\leftarrow 1$;

while  _ $(gap>\varepsilon)$ _ do

$\mathbf{p}_{i}^{old}\leftarrow\mathbf{p}_{i}$;

$\mathbf{p}_{i}\leftarrow\mathbf{A}\mathbf{p}_{i}^{old}+\mathbf{b}_{i}$;

$gap\leftarrow\left\lVert\mathbf{p}_{i}-\mathbf{p}_{i}^{old}\right\rVert$;

$t\leftarrow t+1$;

end while

return $\mathbf{p}$;

Algorithm 1 Power-NF from [10]: Power method for the Newsfeed probabilities
$\mathbf{p}_{i}$ related to user $i$.

To rank all users, we need to calculate $\psi_{i}$ for all $i\in\mathcal{N}$,
i.e. we need the vector
$\bm{\psi}=(\psi_{1},\;\psi_{2},\;\cdots,\;\psi_{N})^{T}\in\mathbb{R}_{+}^{N}$.
Thus, the iterative algorithm in (8) has to be applied $N$ times to solve the
linear system (6) for each origin $i\in\mathcal{N}$. After that, we can map
each vector $\mathbf{p}_{i}$ to $\mathbf{q}_{i}$ through (7) and finally using
(1) we can derive each $\psi_{i}$.

As shown in [10, Theorem 5], in the homogeneous case when all users in the
network have the same posting and re-posting activity, i.e. $\forall
n\in\mathcal{N},\;\lambda^{(n)}=\lambda$ and $\mu^{(n)}=\mu$, it holds
$\bm{\psi}=\bm{\pi}$ where $\bm{\pi}$ is the PageRank vector with a damping
factor of $\alpha=\frac{\mu}{\lambda+\mu}$.

Problem Statement: The current computation of the $\psi$-score is very slow
(compared e.g. to PageRank), as it requires to solve $N$ systems of $N$
equations each, where $N$ is the number of users. This is especially
problematic when we want to calculate the $\psi$-score in real networks with
very large number of users.

Hence, given a directed social graph $\mathcal{G}=(\mathcal{N},\mathcal{E})$
where each node has information about its activity of posting and sharing, we
aim for an algorithm that computes the $\psi$-score for all nodes in the graph
as fast as PageRank, which is the classically used alternative.

##  III Proposed Method

In this section, we present our main contribution: a distributed algorithm for
the fast approximation of the $\psi$-score, to any given requirement of error-
tolerance. To attain this goal, we first derive a new expression for the
$\psi$-score that reduces its computational complexity from solving $N$ linear
systems of size $N$, to just one linear system of size $N$ summarizing all
users. Then, we show that power-iteration rules can be used to approximate the
solution of such linear system, something that can furthermore profit by
distributed computations. This way, we obtain an algorithm for the derivation
of the $\psi$-score that matches the empirical complexity of PageRank, while
we benefit by the expressiveness of the score due to the incorporation of
richer information.

###  III-A Expressing with a single linear system

For the sake of notation clarity, we will express the $\psi$-score of the
network in terms of the matrices $\mathbf{P,Q,B,D}$, which are defined in
Table I. This matrix notation allows us to express the vector of $\psi$-scores
as

$\displaystyle\bm{\psi}^{T}=\frac{1}{N}\mathbf{1}^{T}\mathbf{Q}$ (9)

Thus, the problem of computing the $\psi$-score amounts to computing a matrix
$\mathbf{Q}$, whose columns require knowledge of the vectors $\mathbf{p}_{i}$.
Finding each $\mathbf{p}_{i}$ amounts to solving a linear system as shown in
Eq. (6). Note in (6)-(7) that irrespective of the user origin $i$, all
equations use the same matrices $\mathbf{A}$ and $\mathbf{C}$. The solution of
each system can be expressed as (see [10, Lemma 2]):

$\mathbf{p}_{i}=(\mathbf{I}-\mathbf{A})^{-1}\mathbf{b}_{i}=\sum\limits_{t=0}^{\infty}{\mathbf{A}^{t}}\mathbf{b}_{i}$
(10)

The linearity of (10) implies that the solution for all users can be written
in a single matrix form as

$\mathbf{P}=(\mathbf{I}-\mathbf{A})^{-1}\mathbf{B}=\sum\limits_{t=0}^{\infty}{\mathbf{A}^{t}}\mathbf{B}$
(11)

Then, by leveraging (7) in matrix form and the series development in (11), we
can further develop Eq. (9) as follows:

$\displaystyle\bm{\psi}^{T}$
$\displaystyle=\frac{1}{N}\mathbf{1}^{T}\mathbf{Q}$
$\displaystyle=\frac{1}{N}\mathbf{1}^{T}(\mathbf{CP}+\mathbf{D})$
$\displaystyle=\frac{1}{N}\mathbf{1}^{T}\left[\mathbf{C}(\mathbf{I}-\mathbf{A})^{-1}\mathbf{B}+\mathbf{D}\right]$
$\displaystyle=\frac{1}{N}\mathbf{1}^{T}\left[\mathbf{C}\left(\sum\limits_{t=0}^{\infty}\mathbf{A}^{t}\right)\mathbf{B}+\mathbf{D}\right]$
$\displaystyle=\frac{1}{N}\left[\mathbf{1}^{T}\mathbf{C}\left(\sum\limits_{t=0}^{\infty}\mathbf{A}^{t}\right)\mathbf{B}+\mathbf{1}^{T}\mathbf{D}\right]$
$\displaystyle=\frac{1}{N}\left[\left(\sum\limits_{t=0}^{\infty}\mathbf{1}^{T}\mathbf{C}\mathbf{A}^{t}\right)\mathbf{B}+\mathbf{1}^{T}\mathbf{D}\right]$

Now, if we let $\mathbf{c}^{T}:=\mathbf{1}^{T}\mathbf{C}$ and
$\mathbf{d}^{T}:=\mathbf{1}^{T}\mathbf{D}$, we have that $\mathbf{c}$ and
$\mathbf{d}$ are vectors in $\mathbb{R}_{+}^{N}$ (the diagonals of the
matrices $\mathbf{C}$ and $\mathbf{D}$). These derivations allow us to state
our main theoretical result, which is the following expression for the
$\psi$-score:

$\displaystyle\bm{\psi}^{T}=\frac{1}{N}\left[\left(\sum\limits_{t=0}^{\infty}\mathbf{c}^{T}\mathbf{A}^{t}\right)\mathbf{B}+\mathbf{d}^{T}\right]$
(12)

Eq. (12) indicates that it is not necessary to solve $N$ linear systems of
size $N$ in order to compute the $\psi$-score. Instead, it suffices to solve
one linear system of size $N$ which is expressed as an infinite sum of vector
matrix products, and can even be calculated distributedly. The drawback is
that we lose with this new expression the calculation of intermediate detailed
influence quantities $\mathbf{p}_{i}$ and $\mathbf{q}_{i}$. These are not
necessary to be explicitly known for the calculation of the $\psi$-score as
(12) shows, but they do contain valuable information that could be useful for
certain practical applications. Accelerating the detailed calculation of all
by-products of the $\psi$-score is an important challenge for future work.

###  III-B Convergence of the sum

In this subsection, we demonstrate that the infinite sum in Eq. (12) can be
evaluated by power-iteration rules. The sum always converges (for the specific
expression of $\mathbf{A}$) and can thus be truncated to a finite number of
terms to find an approximation of the $\psi$-score. To show this, let us
define:

$\displaystyle\mathbf{s}_{t}^{T}=\sum\limits_{\tau=0}^{t}\mathbf{c}^{T}\mathbf{A}^{\tau},$
$\displaystyle\mathbf{s}=\lim_{t\rightarrow\infty}\mathbf{s}_{t},$ (13)

where $T$ denotes (matrix or vector) transpose, $t$ is the iteration index and
$\tau$ the running summation index. We stress that Eq. (13) converges if the
spectral radius bound of $\mathbf{A}$ is in the open unit interval. This
convergence is guaranteed by Lemmas 1 and 2 of [10] which use the fact that
$\mathbf{A}$ is a sub-stochastic matrix.

###  III-C Truncation of the sum

We now proceed to study the approximation error after truncating the sum to
the first $t$ terms. For this, we express the truncation of the sum in the
following recursive way:

$\mathbf{s}_{t}^{T}=\mathbf{s}_{t-1}^{T}\mathbf{A}+\mathbf{c}^{T}$ (14)

where $\mathbf{s}_{0}=\mathbf{c}$. This allows us to define the following gap
parameter

$\varepsilon_{t}=\left\lVert\mathbf{s}_{t}^{T}-\mathbf{s}_{t-1}^{T}\right\rVert$
(15)

which is useful for the termination rule of the algorithm. We do not specify
which norm to use in the gap formula because all the work done in this part is
valid with any norm. Nevertheless, in our implementation, we choose to use the
$L^{1}$-norm as it is the common way to deal with series approximation with
power-iterations [16]. We stop the recursion if the trajectory towards the
fixed point between two successive iterations does not change, i.e. when the
following condition is satisfied:

$\varepsilon_{t}\leq\varepsilon$ (16)

An important question to address is therefore how the tolerance set to
terminate the series summation $\mathbf{s}_{t}$ in (13) affects the
approximation of the $\psi$-score in (12). To address this point, let us
define $\bm{\psi}_{t}$ as the approximation of $\bm{\psi}$ by replacing (14)
into (12) using $t$ terms. Then, we can define the gap between two
approximations of the $\psi$-score with different truncations as:

$\delta_{t}=\left\lVert\bm{\psi}_{t}^{T}-\bm{\psi}_{t-1}^{T}\right\rVert$ (17)

Our goal is thus to express $\delta_{t}$ as a function of $\varepsilon_{t}$ or
to find an upper-bound relation between them to not only set a convergence
condition for the vector $\mathbf{s}_{t}$ but more importantly for the vector
$\bm{\psi}_{t}$. We derive this relationship by expressing the difference
$\bm{\psi}_{t}-\bm{\psi}_{t-1}$ in terms of $\mathbf{s}_{t-1}$ and
$\mathbf{s}_{t}$ as follows:

$\displaystyle\bm{\psi}_{t}^{T}-\bm{\psi}_{t-1}^{T}$
$\displaystyle=\frac{1}{N}(\mathbf{s}_{t}^{T}\mathbf{B}+\mathbf{d}^{T})-\frac{1}{N}(\mathbf{s}_{t-1}^{T}\mathbf{B}+\mathbf{d}^{T})$
$\displaystyle=\frac{1}{N}(\mathbf{s}_{t}^{T}-\mathbf{s}_{t-1}^{T})\mathbf{B}$
(18)

Eq. (III-C) then allows to see that:

$\displaystyle\delta_{t}$
$\displaystyle=\frac{1}{N}\left\lVert(\mathbf{s}_{t}^{T}-\mathbf{s}_{t-1}^{T})\mathbf{B}\right\rVert$
$\displaystyle\leq\frac{1}{N}\left\lVert\mathbf{s}_{t}^{T}-\mathbf{s}_{t-1}^{T}\right\rVert\left\lVert\mathbf{B}\right\rVert$
$\displaystyle\delta_{t}$
$\displaystyle\leq\frac{\varepsilon_{t}\left\lVert\mathbf{B}\right\rVert}{N},$
(19)

With this upper bound, we can thus select a termination tolerance
$\varepsilon$ that guarantees that the $\psi$-score trajectory did not change
more than $\delta$. Namely, if we terminate the algorithm by satisfying the
condition $\varepsilon_{t}\left\lVert\mathbf{B}\right\rVert\leq\varepsilon$,
then we are certain that $\delta_{t}$ is lower than $\frac{\varepsilon}{N}$.
We coin our algorithm for the fast approximation of the $\psi$-score as
Power-$\psi$. It is summarized in Algorithm 2.

input : $N$ number of users, $N\times N$ matrices $\mathbf{A}$ and
$\mathbf{B}$, two vectors $\mathbf{c}$ and $\mathbf{d}$,
$\mathbf{s}$-tolerance $\varepsilon$

output : vector $\bm{\psi}$ with the $\psi$-score of all users

$\mathbf{s}\leftarrow\mathbf{c}$;

$B\\_norm\leftarrow\left\lVert\mathbf{B}\right\rVert$;

$t\leftarrow 0$;

$gap\leftarrow 1$;

while  _ $(gap>\varepsilon)$ _ do

$\mathbf{s}_{old}\leftarrow\mathbf{s}$;

$\mathbf{s}^{T}\leftarrow\mathbf{s}_{old}^{T}\mathbf{A}+\mathbf{c}$;

$gap\leftarrow B\\_norm\left\lVert\mathbf{s}_{old}-\mathbf{s}\right\rVert$;

$t\leftarrow t+1$;

end while

$\bm{\psi}^{T}\leftarrow\frac{1}{N}\left(\mathbf{s}^{T}\mathbf{B}+\mathbf{d}^{T}\right)$;

return $\bm{\psi}$;

Algorithm 2 Power-$\psi$: Power iteration based algorithm for the $\psi$-score
vector.

###  III-D Relationship between $\psi$-score and PageRank

For the sake of completeness, in this subsection we discuss the relation
between equation (14) in the homogeneous activity case and PageRank [16].
Then, we highlight the difference between our proposed method and the standard
power-iteration algorithm for PageRank.

Let us define $\bm{\pi}$ as the PageRank vector and
$\mathbf{W}=\mathbf{D}_{out}^{-1}\mathbf{L}$ the random walk transition matrix
where $\mathbf{L}$ is the adjacency matrix and $\mathbf{D}_{out}$ is the
diagonal matrix in which each diagonal entry $(i,i)$ is the out-degree of node
$i$. It was shown in [10, Theorem 5] that in the homogeneous activity case
(where $\lambda^{(i)}=\lambda$ and $\mu^{(i)}=\mu$, $\forall i$) the
$\psi$-score coincides with PageRank. This is, $\bm{\psi}=\bm{\pi}$ with
$\alpha=\frac{\mu}{\lambda+\mu}$ being PageRank’s damping factor.

In the case of homogeneous activity, the vectors and matrices involved in the
$\psi$-score expression of Eq. (12) take the following form:

$\displaystyle\mathbf{A}$ $\displaystyle=\alpha\mathbf{W}$
$\displaystyle\mathbf{B}$ $\displaystyle=(1-\alpha)\mathbf{W}$
$\displaystyle\mathbf{c}$ $\displaystyle=\alpha\mathbf{1}$
$\displaystyle\mathbf{d}$ $\displaystyle=(1-\alpha)\mathbf{1}.$

Then, after substitution in Eq. (14), it becomes:

$\mathbf{s}_{t}^{T}=\alpha\mathbf{s}_{t-1}^{T}\mathbf{W}+\alpha\mathbf{1}^{T}$
(20)

By substituting the above expression in (12) using the homogeneous values for
$\mathbf{B}$ and $\mathbf{d}$ we get

$\displaystyle\bm{\psi}_{t}^{T}$ $\displaystyle=$
$\displaystyle\frac{1}{N}\left(\mathbf{s}_{t}^{T}(1-\alpha)\mathbf{W}+(1-\alpha)\mathbf{1}^{T}\right)$
(21) $\displaystyle=$
$\displaystyle\frac{1}{N}\left((\alpha\mathbf{s}_{t-1}^{T}\mathbf{W}+\alpha\mathbf{1}^{T})(1-\alpha)\mathbf{W}+(1-\alpha)\mathbf{1}^{T}\right)$
$\displaystyle=$
$\displaystyle\alpha\left(\frac{1}{N}(\mathbf{s}_{t-1}^{T}(1-\alpha)\mathbf{W}+(1-\alpha)\mathbf{1}^{T})\right)\mathbf{W}+$
$\displaystyle+$ $\displaystyle\frac{(1-\alpha)}{N}\mathbf{1}^{T}$
$\displaystyle=$
$\displaystyle\alpha\bm{\psi}_{t-1}^{T}\mathbf{W}+\frac{(1-\alpha)}{N}\mathbf{1}^{T}$

For comparison, the above iteration actually evolves as the power-method for
PageRank

$\bm{\pi}_{t}^{T}=\alpha\bm{\pi}_{t-1}^{T}\mathbf{W}+\frac{1-\alpha}{N}\mathbf{1}^{T}.$
(22)

Notice that the proposed Power-$\psi$ algorithm in the general non-homogeneous
case uses power-iteration to first approximate the series $\mathbf{s}$ up to
some truncation step $t$ based on the tolerance required; after convergence it
multiplies by $\mathbf{B}$ and adds $\mathbf{d}$ at the very end, to avoid
unnecessary matrix vector multiplications. This is the main difference in
implementation compared to the power-iteration with PageRank. As shown above,
both methods essentially converge to the same value in the homogeneous case,
but the critical detail is the tolerance, which in the case of $\psi$-score is
on the $\mathbf{s}$-vector, whereas in the case of PageRank is on the
PageRank-vector itself.

##  IV The psi-score Python package

The psi-score project222https://github.com/NouamaneA/psi-score is a software
library that we have developed in Python. It is open source and licensed under
the terms of the MIT license, allowing free usability.

This package is a practical tool for the community to easily use the
$\psi$-score metric in Network Science projects without the need to develop
each algorithm individually. Both methods mentioned in this work are
accessible.

Inspired by scikit-network [19] and scikit-learn [20], the project has a
similar Application Programming Interface (API) in order to be user-friendly.

##  V Numerical Evaluation

In this section, we aim to evaluate the performance of our Power-$\psi$
algorithm (Algorithm 2). We employ the following two metrics to assess
complexity: $(i)$ the number of matrix-vector multiplications to reach a
targeted tolerance; and $(ii)$ the actual runtime of the algorithms in a
machine with the following characteristics: 256GB of RAM, 4 CPUs (with 6
cores) at 3GHz. The former allows to assess and compare the complexity of
distributed algorithms irrespective of the technology where the algorithms are
implemented. The latter allows to get insights on the scalability of the
algorithms in a real-world use case.

We run three series of experiments. The first one aims to validate the ability
of Power-$\psi$ to approximate the true $\psi$-score by comparing the
approximation error Power-$\psi$ does with the current state-of-the-art
alternative. The second one aims to compare the performance of Power-$\psi$
with Power-NF and PageRank in terms of required number of matrix-vector
multiplications. The third experiment shows that the proposed Power-$\psi$
scales to large graphs. It aims to compare the computation time for the three
mentioned methods.

Each experiment has been run twice: $(i)$ one time to compute the $\psi$-score
with heterogeneous activity for the nodes, i.e.
$\lambda^{(i)}\neq\lambda^{(j)}$ and $\mu^{(i)}\neq\mu^{(j)}$ for at least one
pair of nodes $i,j$ (in practice, the activity frequencies have been chosen
uniformly at random within the interval $(0,1)$); $(ii)$ the second time we
run the experiments to compute the $\psi$-score algorithm with homogeneous
activity $\lambda^{(n)}=\lambda=0.15$ and $\mu^{(n)}=\mu=0.85$ for all $n$
(case that reduces to PageRank as discussed in Section III-D) and compare it
to the standard power-iteration for PageRank with
$\alpha=\frac{\mu}{\lambda+\mu}=0.85$ for consistency in our comparisons.

In the experiments we vary the quantity $\varepsilon$ related to tolerance.
Note here that each one of the algorithms (Power-$\psi$, Power-NF and
PageRank’s power-method) have a different notion of tolerance, depending on
the criterion of convergence each algorithm applies. Specifically, for
PageRank we refer to $\pi$-tolerance
$\left\lVert\bm{\pi}_{t}^{T}-\bm{\pi}_{t-1}^{T}\right\rVert\leq\varepsilon$
because the convergence is based on the change in PageRank values; for Power-
NF we refer to $p$-tolerance
$\left\lVert\mathbf{p}_{i}(t)-\mathbf{p}_{i}(t-1)\right\rVert\leq\varepsilon$
because the convergence is based on the change in the news-feed $p$-values;
for Power-$\psi$ we refer to $s$-tolerance
$\left\lVert\mathbf{s}_{t}^{T}-\mathbf{s}_{t-1}^{T}\right\rVert$ because as
shown above the convergence is based on the change in the series $s$-values.
Note in the case of Power-NF and Power-$\psi$ further processing needs to be
done after convergence in order to eventually derive the resulting
$\psi$-scores. Hence, when we set some value $\varepsilon$ for the x-tolerance
in a specific method, we are interested to evaluate the resulting relative
error in the result $\bm{\psi}_{\varepsilon}$ of this method compared to the
true value $\bm{\psi}_{true}$, defined as follows

$error=\frac{\left\lVert\bm{\psi}_{true}-\bm{\psi}_{\varepsilon}\right\rVert_{2}}{\left\lVert\bm{\psi}_{true}\right\rVert_{2}}$
(23)

Table II lists the datasets used for our experiment, which are available in
the Konect repository [21]. Twitter and Facebook graphs, where users have
posting and re-posting activity, are natural application scenarios for our
method. The two citation datasets (DBLP and HepPh arXiv) are also relevant
application cases for our algorithm. This is because paper publications per
author can be seen as equivalent to posting activity, and paper citations by
author can be seen as equivalent to reposting activity. In all datasets, only
the graph topology is available, thus user activity has been generated either
at random for $(i)$ or chosen for $(ii)$ in order to have
$\frac{\mu}{\lambda+\mu}=0.85$.

In order to be able to assess the error with respect to a true $\psi$-score,
for experiments 1 and 2 we focus entirely on the DBLP citation network [22]
that has a moderate size and it is thus realistic to exactly solve the linear
system. Experiment 3 addresses the remainder of datasets from Table II.

Our implementation is available on GitHub333https://github.com/NouamaneA/fast-
algorithm-psiscore as well as the code for numerical evaluation.

TABLE II: Datasets used for evaluation Dataset name | Type | #Nodes | #Edges | References
---|---|---|---|---
DBLP | Citation Network | 12 591 | 49 743 | [21, 22]
Twitter | Social Network | 465 017 | 834 797 | [21, 23]
Facebook | Social Network | 63 731 | 817 035 | [21, 24]
HepPh arXiv | Citation Network | 34 546 | 421 578 | [21, 25]

###  V-A Experiment 1

In the first experiment, we assess the approximation error made with
Power-$\psi$. We compare it with the approximation errors made with the other
methods Power-NF and PageRank’s power-method. This experiment seeks to show
that our proposed algorithm can approximate the $\psi$-score vector (also the
PageRank vector in the homogeneous case) at least as well as existing methods.

We compute the true $\psi$-score vector: once for the heterogeneous case $(i)$
and a second time for the homogeneous case. Then, we run each algorithm with
several targeted tolerance from $10^{-9}$ to $10^{-1}$. Finally we calculate
the approximation error with equation (23) for each realization of each
method.

The results are shown for heterogeneous activity $(i)$ in Fig. 2 and for
homogeneous activity $(ii)$ in Fig. 3. We fix the same level of $\mathbf{p}$\-
and $\mathbf{s}$-tolerance for both algorithms Power-NF and Power-$\psi$ as
well as PageRank power iteration and derive the error from (23). We see in
Fig. 2 and Fig. 3, for the same targeted tolerance, that the error made with
Power-$\psi$ method is lower than the approximation error made with Power-NF
and PageRank, which validates the convergence of our method to the true value.

###  V-B Experiment 2

The second experiment evaluates the number of matrix-vector multiplications
that are required to compute the $\psi$-score and compare the methods with the
same precision.

We run the proposed algorithm and the state-of-the-art algorithm presented in
[10] with different tolerance parameters from $10^{-9}$ from $10^{-1}$. We
measure the number of matrix-vector multiplications done for each method as
well as the relative error calculated with (23) in order to compare the
algorithms within a given precision. Figure 4 shows the result in the case
where we compute the $\psi$-score with $(i)$ heterogeneous activity and Figure
5 shows the result when computing $\psi$-score in the specific case of $(ii)$
homogeneous activity that gives the PageRank vector.

We see clearly that the proposed algorithm Power-$\psi$ outperforms the old
one Power-NF by several orders of magnitude. We can also observe that the
difference in the number of matrix-vector multiplications compared to PageRank
is very small, which renders the two methods computationally equivalent, at
least for the power-iteration method.

###  V-C Experiment 3

In the third experiment, we seek to demonstrate that our proposed algorithm
scales to real-world networks with several hundreds of thousands of edges, at
least as well as PageRank does. For this purpose, we chose one targeted
tolerance (for space reasons) and we measured the computation time after
running each method. Table III shows the results we obtained in the case where
we compute the $\psi$-score with heterogeneous activity $(i)$ and Table IV
shows the results in the case of homogeneity $(ii)$ when the $\psi$-score
vector equals the PageRank vector.

Here again, the proposed Power-$\psi$ outperforms by far the state-of-the-art
Power-NF. We can also conclude that the proposed algorithm is close to
PageRank’s power-method in terms of computation time. In some cases,
Power-$\psi$ takes longer but that is because the algorithm needs to perform
other operations after the convergence of the approximated series.

Figure 2: Experiment 1 $(i)$ with DBLP: Precision assessment of the proposed algorithm compared to the state-of-the-art method. Heterogeneous case where users do not necessarily have the same activity. Figure 3: Experiment 1 $(ii)$ with DBLP: Precision assessment of the proposed algorithm compared to the state-of-the-art method. Homogeneous case where all users have the same activity and $\psi$-score=PageRank. Figure 4: Experiment 2 $(i)$ with DBLP: Comparison of the proposed algorithm with the state-of-the-art method in terms of number of matrix-vector multiplications. Heterogeneous case where users do not necessarily have the same activity. Figure 5: Experiment 2 $(ii)$ with DBLP: Comparison of the proposed algorithm with the state-of-the-art method and PageRank’s power-method in terms of number of matrix-vector multiplications. Homogeneous case where all users have the same activity and $\psi$-score=PageRank. TABLE III: Experiment 3 $(i)$: Computation time evaluation with $\varepsilon$ set to $10^{-9}$. Heterogeneous case where users do not necessarily have the same activity. Dataset | Power-NF [10] | Power-$\psi$
---|---|---
DBLP | 17.805 sec | 0.029 sec
Facebook | 1764.226 sec | 0.307 sec
Twitter | 14526.039 sec | 0.634 sec
HepPh | 272.358 sec | 0.622 sec
TABLE IV: Experiment 3 $(ii)$: Computation time evaluation with $\varepsilon$ set to $10^{-9}$. Homogeneous case where all users have the same activity and $\psi$-score=PageRank. Dataset | PageRank [16] | Power-NF [10] | Power-$\psi$
---|---|---|---
DBLP | 0.023 sec | 20.775 sec | 0.034 sec
Facebook | 0.308 sec | 2253.302 sec | 0.454 sec
Twitter | 0.584 sec | 17411.146 sec | 0.806 sec
HepPh | 0.361 sec | 360.769 sec | 0.908 sec

##  VI Conclusion and Future Work

We have proposed a new method to compute the $\psi$-score, which is influence
of users in a social network (or in any network with nodes having a posting a
sharing activity). We have shown that the proposed algorithm converges faster
than the state-of-the-art method. The algorithm enables scalability for real-
world datasets for the $\psi$-score.

We have also provided an open source software library that implements the
existing algorithm for the $\psi$-score calculation as well as the novel
approach.

The proposed method gives directly the influence of every user over the entire
network. But in some cases, one might need to use the influence of a user on a
specific user, which is giver by the OSP model but skipped by the proposed
Power-$\psi$ to go faster to the point. In the future, we are interested in
further study on the $\psi$-score in order to explore the possibility to get
these other information potentially valuable for some applications.

## Acknowledgements

The work of N.A. is funded directly by Sorbonne University, CNRS-LIP6
laboratory by an internal project. The research of A.G. is supported by the
French National Agency of Research (ANR) through the FairEngine project under
Grant ANR-19-CE25-0011. E.B. is funded by the French National Agency of
Research (ANR) under the FiT LabCom grant.

## References

  * [1] D. Kempe, J. Kleinberg, and E. Tardos, “Maximizing the spread of influence through a social network,” in _Proceedings of the Ninth ACM SIGKDD International Conference on Knowledge Discovery and Data Mining_ , ser. KDD ’03. New York, NY, USA: Association for Computing Machinery, 2003, p. 137–146. [Online]. Available: https://doi.org/10.1145/956750.956769
  * [2] M. S. Alvim, B. Amorim, S. Knight, S. Quintero, and F. Valencia, “A multi-agent model for polarization under confirmation bias in social networks,” in _International Conference on Formal Techniques for Distributed Objects, Components, and Systems_. Springer, 2021, pp. 22–41.
  * [3] K. Bok, Y. Noh, J. Lim, and J. Yoo, “Hot topic prediction considering influence and expertise in social media,” _Electronic Commerce Research_ , vol. 21, no. 3, pp. 671–687, 2021.
  * [4] C. Musto, P. Lops, M. de Gemmis, and G. Semeraro, “Context-aware graph-based recommendations exploiting personalized pagerank,” _Knowl. Based Syst._ , vol. 216, p. 106806, 2021. [Online]. Available: https://doi.org/10.1016/j.knosys.2021.106806
  * [5] V. Cabannes, L. Pillaud-Vivien, F. Bach, and A. Rudi, “Overcoming the curse of dimensionality with laplacian regularization in semi-supervised learning,” _Advances in Neural Information Processing Systems_ , vol. 34, 2021.
  * [6] S. Goel, D. J. Watts, and D. G. Goldstein, “The structure of online diffusion networks,” in _Proceedings of the 13th ACM Conference on Electronic Commerce_ , ser. EC ’12. New York, NY, USA: Association for Computing Machinery, 2012, p. 623–638. [Online]. Available: https://doi.org/10.1145/2229012.2229058
  * [7] Z. Wan, Y. Mahajan, B. W. Kang, T. J. Moore, and J.-H. Cho, “A survey on centrality metrics and their network resilience analysis,” _IEEE Access_ , vol. 9, pp. 104 773–104 819, 2021.
  * [8] M. Cha, H. Haddadi, F. Benevenuto, and K. Gummadi, “Measuring user influence in twitter: The million follower fallacy,” _Proceedings of the International AAAI Conference on Web and Social Media_ , vol. 4, no. 1, pp. 10–17, May 2010. [Online]. Available: https://ojs.aaai.org/index.php/ICWSM/article/view/14033
  * [9] C. Wilson, A. Sala, K. P. N. Puttaswamy, and B. Y. Zhao, “Beyond social graphs: User interactions in online social networks and their implications,” _ACM Trans. Web_ , vol. 6, no. 4, nov 2012. [Online]. Available: https://doi.org/10.1145/2382616.2382620
  * [10] A. Giovanidis, B. Baynat, C. Magnien, and A. Vendeville, “Ranking online social users by their influence,” _IEEE/ACM Trans. Netw._ , vol. 29, no. 5, pp. 2198–2214, 2021. [Online]. Available: https://doi.org/10.1109/TNET.2021.3085201
  * [11] C. Correa, T. Crnovrsanin, and K.-L. Ma, “Visual reasoning about social networks using centrality sensitivity,” _IEEE Transactions on Visualization and Computer Graphics_ , vol. 18, no. 1, pp. 106–120, 2012.
  * [12] R. Puzis, Y. Altshuler, Y. Elovici, S. Bekhor, Y. Shiftan, and A. Pentland, “Augmented betweenness centrality for environmentally aware traffic monitoring in transportation networks,” _J. Intell. Transp. Syst._ , vol. 17, no. 1, pp. 91–105, 2013. [Online]. Available: https://doi.org/10.1080/15472450.2012.716663
  * [13] A. Tizghadam and A. Leon-Garcia, “Betweenness centrality and resistance distance in communication networks,” _IEEE Network_ , vol. 24, no. 6, pp. 10–16, 2010.
  * [14] M. Ashtiani, A. Salehzadeh-Yazdi, Z. Razaghi-Moghadam, H. Hennig, O. Wolkenhauer, M. Mirzaie, and M. Jafari, “A systematic survey of centrality measures for protein-protein interaction networks,” _bioRxiv_ , 2017. [Online]. Available: https://www.biorxiv.org/content/early/2017/10/09/149492
  * [15] P. R. Miller, P. S. Bobkowski, D. Maliniak, and R. B. Rapoport, “Talking politics on facebook: Network centrality and political discussion practices in social media,” _Political Research Quarterly_ , vol. 68, no. 2, pp. 377–391, 2015. [Online]. Available: https://doi.org/10.1177/1065912915580135
  * [16] L. Page, S. Brin, R. Motwani, and T. Winograd, “The pagerank citation ranking: Bringing order to the web,” 1998.
  * [17] J. Whang, A. Lenharth, I. Dhillon, and K. Pingali, “Scalable data-driven pagerank: Algorithms, system issues, and lessons learned,” 08 2015, pp. 438–450.
  * [18] E. Bautista and M. Latapy, “A local updating algorithm for personalized pagerank via chebyshev polynomials,” _Soc. Netw. Anal. Min._ , vol. 12, no. 1, p. 31, 2022. [Online]. Available: https://doi.org/10.1007/s13278-022-00860-5
  * [19] T. Bonald, N. de Lara, Q. Lutz, and B. Charpentier, “Scikit-network: Graph analysis in python,” _Journal of Machine Learning Research_ , vol. 21, no. 185, pp. 1–6, 2020. [Online]. Available: http://jmlr.org/papers/v21/20-412.html
  * [20] L. Buitinck, G. Louppe, M. Blondel, F. Pedregosa, A. Mueller, O. Grisel, V. Niculae, P. Prettenhofer, A. Gramfort, J. Grobler, R. Layton, J. VanderPlas, A. Joly, B. Holt, and G. Varoquaux, “API design for machine learning software: experiences from the scikit-learn project,” in _ECML PKDD Workshop: Languages for Data Mining and Machine Learning_ , 2013, pp. 108–122.
  * [21] J. Kunegis, “KONECT – The Koblenz Network Collection,” in _Proc. Int. Conf. on World Wide Web Companion_ , 2013, pp. 1343–1350. [Online]. Available: http://dl.acm.org/citation.cfm?id=2488173
  * [22] M. Ley, “The DBLP computer science bibliography: Evolution, research issues, perspectives,” in _String Processing and Information Retrieval, 9th International Symposium, SPIRE 2002, Lisbon, Portugal, September 11-13, 2002, Proceedings_ , ser. Lecture Notes in Computer Science, A. H. F. Laender and A. L. Oliveira, Eds., vol. 2476. Springer, 2002, pp. 1–10. [Online]. Available: https://doi.org/10.1007/3-540-45735-6_1
  * [23] M. De Choudhury, Y. Lin, H. Sundaram, K. Candan, L. Xie, and A. Kelliher, “How does the data sampling strategy impact the discovery of information diffusion in social media?” in _ICWSM 2010 - Proceedings of the 4th International AAAI Conference on Weblogs and Social Media_ , ser. ICWSM 2010 \- Proceedings of the 4th International AAAI Conference on Weblogs and Social Media, Dec. 2010, pp. 34–41, 4th International AAAI Conference on Weblogs and Social Media, ICWSM 2010 ; Conference date: 23-05-2010 Through 26-05-2010.
  * [24] B. Viswanath, A. Mislove, M. Cha, and K. P. Gummadi, “On the evolution of user interaction in facebook,” in _Proceedings of the 2nd ACM Workshop on Online Social Networks_ , ser. WOSN ’09. New York, NY, USA: Association for Computing Machinery, 2009, p. 37–42. [Online]. Available: https://doi.org/10.1145/1592665.1592675
  * [25] J. Leskovec, J. Kleinberg, and C. Faloutsos, “Graph evolution: Densification and shrinking diameters,” _ACM Trans. Knowl. Discov. Data_ , vol. 1, no. 1, p. 2–es, mar 2007. [Online]. Available: https://doi.org/10.1145/1217299.1217301