On Modeling and Predicting Individual Paper Citation Count over Time

1 Introduction and related work

Integral to the success of scientific research is the impact of works. Paper citation and its derivatives e.g. g-index [Egghe, 2006], H-index [Hirsch, 2005; Acuna et al., 2012] have become popular measures to gauge the journals, scholars, labs, departments, and institutes [Fuyuno and Cyranoski, 2006], despite their well-known lack of predictive power to future impact [Wang et al., 2013]: current citations and the derived metrics can only capture past accomplishments.

A candidate's potential future impact e.g. his/her citation count often plays a more important role for policy/decision making concerning with recruitment, promotion and funding, because the ultimate question is: Who will be the most successful in this position, with this fellowship? When an earlycareer candidate is selected for a tenure-track position, it is an investment. In those institutions with low tenure rate, this can amount to an outright bet on one scientist who acquire a start-up package up to millions of dollars [Stephan, 2012].

The scientific literature is turning into an unbounded collection such that it becomes intimidating to have a thorough comprehension on relevant papers even in one area. To find frontier research materials, there is also a need for identifying the pertinent and influential work in a setting where a plenty of papers emerge each day, before they become widely recognized. Also, researchers may re-think if their research is on an exciting path or a dead end that will end the careers prematurely. A (reliable) prediction model can serve as a selfevaluation tool to streamline their research agenda.

⇤Corresponding author is Junchi Yan. The work is partially supported by STCSM (15JC1401700, 14XD1402100), NSFC (61527804, 61221001, 61129001), the 111 Program (B07022).

We agree that the main way of predicting a researcher's future impact is peer-assessment, but also think that algorithmic approaches could be valuable complementary ways, especially for junior scientists representing a group closer to the typical case in which algorithmic approaches will be applied in real academic hiring decisions, under an appropriate mechanism. It is felt that 'pipeline' leaks in the later career decision points, especially confounded with the subjective gender bias in academic career [Ginther and Kahn, 2004].

As a widely recognized metric [Wang et al., 2013] to scientific impact, however, predicting an individual paper's citation count over time is (arguably) very difficult. For instance, a seminal work may start-up by a small number of follow-up papers that builds up to a pioneering work within a field, and it takes a long time before they generate greater impact. Or a researcher may work on a hot topic, and publish a novel method related to this topic, which immediately draw the community's wide attention. Or simply different papers by the same author can have significant citation variation due to various reasons such as the topic, timing, fields, etc. Such heterogeneous citation curves call for advanced models.

We give an overview about the general problem - scientific impact analysis and prediction, and then focus on the literature on (long-term) individual-level paper citation prediction, which is (arguably) more challenging and has become an emerging applied research topic [Wang et al., 2013]. Scientific impact analysis and prediction Since predicting individual paper's citation count looks very challenging, and the skewed distribution of citations often obeys a powerlaw [Dong et al., 2015] or log-normal form [Radicchi et al., 2008], many researchers resort to other more accessible impact analysis and prediction problems. For instance, [Petersen et al., 2014] perform a longitudinal analysis to measure the effect of the central author's reputation on the paper citation rate. [Pan and Fortunato, 2014] give a formal definition concerning the author-wise impact metric Author Impact Factor (AIF) and perform empirical study to verify its capability to capture the trends and variations of the impact of the scientific output of scholars in time. However, they are prescriptive methods with no capability of prediction.

For scientific impact prediction, [Dong et al., 2015] formulates and addresses the classification problem of whether a paper will influence an author's h-index within a pre-defined timeframe, whereby various factors/features are fed into the classification model including publication topic, venue, author's affiliation etc. Similarly, [Acuna et al., 2012] present an approximate formula obtained from linear regression with elastic net regularization, to estimate the future h-index of life scientists of which the factors relate to current h-index, number of written papers of authors, years since first published article etc. [Penner et al., 2013] identify the flaws of the above linear regression model and empirically show this model tend to overestimate the scientist's future impact. Specifically, they suggest that h-index, or other cumulative metrics are inappropriate targets of regression based predictive models in that they contain trivial correlation by construction (see 'Discussion' in [Penner et al., 2013]). In contrast, in this paper, our point process based model does not suffer from this limitation, and directly estimates the transient citation counts at any future time point or period. [Stern, 2014] empirically study the social science top-ranked journals and discover that 'half of the papers in the top 20% in 2012 were already in the top 20% in the year of publication (2006)'.

Paper citation prediction As shown in Fig.1, individual level paper citation prediction is challenging. A line of methods formulate the paper citation prediction problem into a regression task, and examine which covariates are effective input features. [Yan et al., 2011; Chakraborty et al., 2014] extract author-centric attributes (e.g. productivity, co-author network, influence), paper-specific features (e.g. team-size, reference count), and venue-centric features (e.g. venue rank) to build a supervised regression model. [Yu et al., 2012] study the problem of predicting the linkage i.e. citation between a pair of papers. However, their approach is based on link prediction and cannot predict the dynamic citation count at any time point. In the seminal Science paper [Wang et al., 2013], the authors propose a point process based behavioral model which tries to capture the dynamics of the individual paper citations. Its intensity function is a multiplication of three factors: i) fitness term, which is paper's intrinsic value being different from paper to paper; ii) the aging effect over time; iii) the reinforcement term indicating the well-documented fact that highly cited papers are more visible and are more likely to be cited again than less-cited contributions. They employ maximum likelihood estimation for individual paper to infer its set of parameters, which is at the risk of over-fitting as also observed in [Wang et al., 2014b]. To mitigate this issue, [Shen et al., 2014] adopt a Bayesian treatment by using a conjugate prior for the fitness parameter θd. This prior is not paper-specific nor flexible to capture paper's arbitrary profile.

Despite the recent advances in scientific impact prediction and more specifically, paper citation prediction, it is still unclear and even controversial on the reliability and bound of prediction accuracy of a long-term citation prediction model - see the comments [Wang et al., 2014b] and response [Wang et al., 2014a] published in the Science communication papers (http://www.science.com/) after the pioneering work [Wang et al., 2013] and improvement [Shen et al., 2014].

Contribution of the paper This paper is aimed to provide in-depth findings on a recently released real-world dataset.

Specifically, we propose a novel point process model for long-term paper citation prediction, which is also quite general in applicability. Our approach captures the Matthew effect [Merton, 1968] (or accumulated advantage, richer get richer and the poor get poorer) and the recency effect of past citations. In particular, it can help better address the common but unresolved 'second-acts' scenarios in [Wang et al., 2013] 1 (a.k.a. 'Sleep Beauty' [Ke et al., 2015]). The covariates w.r.t the author and paper are also incorporated in the intensity function to improve interpretability and mitigates overfitting.

We also provide an empirical analysis of the predictive power and interpretability of the learned point process model on the public Microsoft Academic Graph [Sinha et al., 2015]. Our method consistently outperforms the state-of-the-arts approaches [Wang et al., 2013; Shen et al., 2014].

We think more importantly, this work provides a new investigation on how effective algorithmic citation prediction can be devised, regarding the recent arguments appear on Science communication [Wang et al., 2014b; 2014a].

2 Model and algorithm

2.1 Model formulation

The received citation count of an individual paper d during time period [0↪ T] is characterized by a time-stamped sequence ¶tdi ♦ni=0 when a citation occurs - which we dub it as event in the setting of point process in this paper: 0 = td0  td1  ▷▷▷  tdi  ▷▷▷  tdn  T . The goal is to model and predict the future citation count over an arbitrary time window given the historical citations and other available covariates.

It is clear that papers having been cited frequently tend to accumulate more citations, especially for recent citations. It is also clear that, with time, even the most novel paper loses its popularity. Some papers, however, seem to have an inherent 'quality' that can be interpreted as a community's recognition to the work. Building on a foundation of the above observations, we derive our prediction model in three regards. Intrinsic popularity The quality of the paper is its intrinsic factor contributing to its popularity. To some extent, the quality can be measured by its paper/author-specific covariates, such as the H-Index [Hirsch, 2005] of the author, the field where the paper is published etc. In line with [Yan et al., 2011; 2012; Chakraborty et al., 2014], we extract a set of covariates for each paper as listed in Table 1. We also plot three of their scatters regarding with the citation count in Fig.2 and some of them exhibit strong correlations. Therefore, we can use these covariates to regress the associated coefficients of the intrinsic popularity. Here the Lasso (`1 norm) can be used to induce sparse coefficients to mitigate overfitting.

1 The 'second-acts' e.g. the citation burst for superconductivity papers after the discovery of high-temperature superconductivity in the 1980s, or delayed impact, like the citation explosion to Erd˝os and R´enyi's work 40 years after their publication [Barab´asi and Albert, 1999], following the emergence of network science [Redner, 2005].

Impact decaying over time A general and common trend is that the paper's attractiveness fades away over time. This can be explained by the fact that a topic goes through its lifecycle and ends up with an out-of-dated status, or still being a hot topic, but the novelty is incorporated in subsequent work that dilutes its impact and relevance to other work.

Recency-sensitive citation triggering Previous methods [Wang et al., 2013; Shen et al., 2014] ignore the time-stamp, and aggregate all past citations to model the intensity, which might be less effective to capture the citation dynamics - see our experiments. We propose to use a self-triggering process a.k.a. Hawkes processes [Hawkes, 1971; Hawkes and Oakes, 1974] which favor more on recent citations and the effect time-decaying window can be controlled by, for computational effectiveness, a Laplace kernel [Yan et al., 2013; 2015]. The 'recency-sensitive' model can naturally address those papers with spiking citation curve being still remain unresolved in [Wang et al., 2013; Shen et al., 2014].

Hence, we define the citation count intensity of paper d:

where sd(t) = (1↪ sd1↪ sd2↪ ▷▷▷↪ sdK) is a row vector encoding the K paper-specific covariates for paper d and β is a column vector for the coefficients. e-w1dt is the aging function accounting for attractiveness decrease since its publication. ↵d is the triggering strength of each citation before current time point t, with the decaying effect e-w2d(t-ti) . For point process, the estimated citation count can in general be computed by integrating Eq.1 over a specified future time period.

2.2 Discussion on peer methods

As an emerging problem, the most relevant work to ours is the Reinforced Poisson Process (RPP) model as presented in [Wang et al., 2013] and [Shen et al., 2014], whereby the latter adds a conjugate prior on the fitness of an individual to the former work [Wang et al., 2013] published in Science studying the problem of long-term individual citation dynamics. We also mention the form of Hawkes processes [Hawkes, 1971; Hawkes and Oakes, 1974] as our model is partly originated from this type of point process. We also solve the learning problem via a tailored ADMM based algorithm.

Reinforced Poisson Process - RPP The seminal work published in Science [Wang et al., 2013] begins to study the fundamental problem for the predictability of long-term citation.

Their point process model involves three individual paperspecific parameters: the relative fitness λi capturing a paper's importance relative to to other papers; immediacy θi governing the time for a paper to reach its citation peak; and a longevity parameter σi accounting for the decay rate of its popularity. By solving a master equation associated with the intensity function, they directly give the equation for the number of citations cti (λi↪ θi↪ σi) at time t (see Eq.2 and Eq.3 in that paper). The parameters are then estimated by leastsquare-fitting, given historical citation data for paper i.

To overcome the problem that maximum likelihood parameter estimation suffers from overfitting, especially for relatively small sample size as we need to train one model for each paper by its citations, [Shen et al., 2014] adopt a Bayesian treatment by adding a conjugate prior on the relative fitness parameter θd via a gamma distribution, for each paper and showcase superior results compared with [Wang et al., 2013] for prediction accuracy. While it incurs some doubts frommore recent study [Wang et al., 2014b]. Hawkes Process The intensity of self-exciting Hawkes pro-

cess [Hawkes, 1971; Hawkes and Oakes, 1974] is given by:

Here a paper-specific parameter θd for paper d is used instead of the parameterized linear regression term βsd by our approach, and other parameters have the similar meaning compared to Eq.1. There are some recent studies/applications on Hawkes processes, [Zhao et al., 2015] use this model to predict the popularity of Twitter, where popularity intensity is determined by a stochastic infectiousness process (Cox process). [Zhou et al., 2013] introduce the low-rank sparsity on the infect matrix formed by the mutually-exciting Hawkes model. [Luo et al., 2015] propose a multi-task learning variation for the mutually-exciting Hawkes model and [Yan et al., 2015] adopt it for sales pipeline modeling.

Compared with the above methods, especially [Wang et al., 2013; Shen et al., 2014], the bullets of our approach are:

i) Introducing paper-specific covariates The paperspecific covariates are involved via parameterizing the quality term (see Eq.1). In RPP and Hawkes process, for each paper, their quality term θd is modeled by one their own parameter respectively, which we think is the main reason for overfitting. One shall note that in their methods, each of these parameters is learned from its past citations of each paper independently. Meanwhile, the covariates are not used which otherwise can play a bridge to cross-distributing the information over training samples i.e. papers. For instance, for those papers with very short observation window and with few citations, it is difficult to interpret and predict the behavior of such papers without exploring the covariates used in Table 1. In this sense, involving the informative covariates will mitigate the over-fitting problem and meanwhile help improve the interpretability of our model. Note in the improved RPP [Shen et al., 2014], they only impose a prior on the global distribution of θd, without using the covariates to parameterize the prior and such valuable information is ignored.

ii) Modeling citation recency by self-exciting kernel We account for citation recency by modeling triggering effect in continuous time space. This feature is inspired from the Hawkes process. To our best knowledge, this is the first time for adapting this component in paper citation prediction. More importantly, we find this recency-weighted triggering model is more appropriate for the citation dynamics, especially for those 'second-acts' and 'delayed-impact' phenomena that once appear in citation history [Redner, 2005].

iii) Additive intensity model We model the relation of the first two components with the third one by an additive composition rather than multiplication used in [Wang et al., 2013; Shen et al., 2014]. In general these two forms have their respective strengths (see more details in Chapter 4 in [Aalen et al., 2008]), in the analogous context of the multiplicative Cox and the additive Aalen functions [Aalen et al., 2008].

Specifically, the additive model decouples the temporal aging from the triggering effect while the multiplicative couples each other. We simplify the temporal aging term by a general decaying kernel for the difficulty to capture various and unknown citation life-cycle patterns. Our additive mechanism can isolate the adverse effect by this coarse design. Moreover, it is mathematically easier and more efficient to learn the additive model than a multiplicative one [Vu et al., 2011] 2 .

2 One technical issue is our additive model does not automatically guarantee the non-negativeness of the first term in Eq.1. Thus we normalize the covariates sd to [0,1] and make them almost always greater than zero, and β are ensured to be positive according to Eq.8.

2.3 Model learning and prediction

The length of time interval between two consecutive citations follows an inhomogeneous Poisson process. Therefore, given that the (i - 1)th citation arrives at ti-1, the probability that the ith citation arrives at ti follows

Then the log-likelihood of time-stamped sequence is:

By plugging Eq.1 into the above function and adding sparsity regularization ♣♣β♣♣1, for Gd(t) = R t 0 gd(t)dt we reach:

where gd(t) = e-w(t-tj) is the triggering kernel in Eq.1. Adding `1 norm renders Eq.2 non-differentiable. We apply the idea of Alternating Direction Method of Multipliers (ADMM) [Boyd et al., 2011] to convert the optimization problem to several sub-problems that are easier to solve. The optimization problem in Eq.2 can be rewritten as the following equivalent form by introducing an auxiliary variable z:

The corresponding augmented Lagrangian of the problem is:

where u is the scaled dual variables corresponding to the constraint β = z, and ⇢ is the penalty parameter, which is usually used as the step size in updating the dual variable. Solving the above augmented Lagrangian problem using ADMM algorithm involves the following sub-problem:

Sλ◁⇢ is soft thresholding [Donoho and Johnstone, 1995].

To update β and ↵ in Eq.5 efficiently, we adopt EM framework to solve this convex problem. The EM step is as follows. Let pki denotes probability that feature k triggers event ti and the pij denotes the probability that event ti triggers event tj.

We empirically iterate the expectation step (Eq.6, 7) maximization step (Eq.8, 9) until convergence:

where B = PNd=1 sdkGd(T) + ⇢(uk - zk). In fact, in our tests, we always find our method converges to a stationary point though its convergence property is not theoretically proved in the paper.

We update w1↪w2 for paper d (below d is omitted for notational simplicity) by gradient descent:

After learning the parameters, we simulate the Hawkes process by Ogata's thinning algorithm [Ogata, 1981] and estimate the predicted citations before time t, denoted by cd (t).

3 Experiments and discussion

3.1 Experimental settings

Dataset and compared methods We perform citation count prediction on the real-world dataset: Microsoft Academic Graph [Sinha et al., 2015] of which the papers are well collected, complete and authorized. We select publications in Computer Science, which consists of 3,539,403 papers authored by 1,598,575 researchers. Two networks are constructed: the co-author collaboration network with 1,598,575 vertices and the other is citation network with time-stamped directed link, indicating when the citation is received.

We further use papers published during 1969-1989 from the so-called 'main' Computer Science venues (refer to http://libra.msra.cn/), including 1,240 journals and 547 conference series, resulting in a total of 47,293 papers. Similar to the protocol in [Wang et al., 2013; Shen et al., 2014], we use papers with more than 5 citations during the first 5 years after publication as training data and predict their citations in the next 10 years. Other fields, Physics, Biology and Chemistry are also evaluated. The improved RPP [Shen et al., 2014] based on [Wang et al., 2013], the Hawkes model [Hawkes, 1971] are implemented and tuned to their best performance.

Two metrics used in [Shen et al., 2014] are also used:

Mean Absolute Percentage Error (MAPE) It measures the average deviation between predicted and true popularity over N papers. Denoting with cd (t) the predicted number of citations for a paper d up to time t and with rd (t) its real number of citations, MAPE is given by 1N PNd=1 ∣∣∣ cd (t)-rd (t) rd (t) ∣∣∣.

Accuracy It measures the fraction of papers correctly predicted for a given error tolerance ✏. Hence the accuracy of popularity prediction on N papers is 1N PNd=1 ∣∣∣d : ∣∣ cd (t)-rd (t) rd (t) ∣∣  ✏∣∣∣. [Shen et al., 2014] set ✏ = 0▷1 on their dataset. We find in our test, our methods always outperforms regardless ✏ and we set ✏ = 0▷3.

3.2 Results and further discussion

MAPE and accuracy They are given in Fig.3 where each column for the first two rows shows the results for Computer Science papers published in journal, conference, and IJCAI respectively. Our method (denoted by OURS) consistently outperforms across different observation time window (7, 10 years). The third row on Biology, Physics and Chemistry reveal that our method performs robust across fields.

Time complexity The time cost for RPP, Hawkes and our method is 0.355, 0.805, 0.860 seconds per iteration. However, our efficient EM framework renders OURS converges with less iterations as illustrated in Fig.4. The total consumption time for RPP, Hawkes and our method is 42.6, 24.1 and

17.2 minutes in average for training the Computer Science models. This suggests our method is efficient and accurate. Visual comparison We also compare the learned intensities of RPP and our method. Fig.5 shows the learned intensity distribution and citation events. Event intensity of RPP relies strongly on the time decaying function while ours can flexibly capture the triggering effect of recent citations.

Study on 'second-acts' papers We conduct experiments on the so-called 'second-acts' phenomenon that papers receive plenty of citations in their late-stage of life cycle. This type of papers is also called in Sleep Beauty by [Ke et al., 2015]. There is a continuous spectrum of delayed recognition where both the hibernation period and the awakening intensity are taken into account. Fig.6 compares both the real and cumulative predicted citations for sleep beauties type of papers, where the observation window is 10 years. 249 papers are chosen by setting i) less than 20 citations in the first 10 years since publication, and ii) larger than 70 citations in the next 10 years. Note the MAPE is worse than the results in Fig.3 while the Hawkes model and our method performs relatively better. We visualize real and predicted citations of two concrete examples: solid line indicates paper [Kahn and Roth, 1971], and dash line [Eklundh, 1986]. RPP increases linearly regardless of the booming citations in the late stage (around the end of the 10 year observation time window) while our method follows the trend more timely and closely.

Interpretability of covariates By using the sparsity regularization (set λ = 2 in Eq.1), we can select the most important and interpretable features. Table.1 ranks the covariates by the amplitude of coefficients. The most important factors are author's authority, such as H-index, author rank, and venue's rank, which relate to the novelty of the scientific works.

4 Conclusion

We present an individual paper citation prediction model. Empirical results suggest that its utility for prediction and interpretability. It also gives an independent study on the argument for how effective algorithmic citation prediction approaches can be devised among Science communication [Wang et al., 2014b; 2014a]. We empirically find robust methods is achievable for individual paper citation prediction by appropriate modeling in line with [Wang et al., 2014a].