Sequential Modeling of Topic Dynamics with Multiple Timescales
TOMOHARU IWATA, NTT Communication Science Laboratories
TAKESHI YAMADA, NTT Science and Core Technology Laboratory Group
YASUSHI SAKURAI and NAONORI UEDA, NTT Communication Science Laboratories

19
We propose an online topic model for sequentially analyzing the time evolution of topics in document collections. Topics naturally evolve with multiple timescales. For example, some words may be used consistently
over one hundred years, while other words emerge and disappear over periods of a few days. Thus, in the
proposed model, current topic-specific distributions over words are assumed to be generated based on themultiscale word distributions of the previous epoch. Considering both the long- and short-timescale dependency yields a more robust model. We derive efficient online inference procedures based on a stochastic
EM algorithm, in which the model is sequentially updated using newly obtained data; this means that pastdata are not required to make the inference. We demonstrate the effectiveness of the proposed method interms of predictive performance and computational efficiency by examining collections of real documentswith timestamps.
Categories and Subject Descriptors: H.2.8 [Database Management]: Database Applications—Datamining; I.2.6 [Artificial Intelligence]: Learning; I.5.1 [Pattern Recognition]: Models—Statistical
General Terms: Algorithms
Additional Key Words and Phrases: Topic model, time-series analysis, online learning, multiscale
ACM Reference Format:
Iwata, T., Yamada, T., Sakurai, Y., and Ueda, N. 2011. Sequential modeling of topic dynamics with multiple
timescales. ACM Trans. Knowl. Discov. Data 5, 4, Article 19 (February 2012), 27 pages.
DOI = 10.1145/2086737.2086739 http://doi.acm.org/10.1145/2086737.2086739

1. INTRODUCTION

Great interest is being shown in developing topic models that can analyze and summarize the dynamics of document collections, such as scientific papers, news articles, andblogs [Ahmed and Xing 2010; AlSumait et al. 2008; Blei and Lafferty 2006; Canini et al.
2009; Iwata et al. 2009; Nallapati et al. 2007; Wang et al. 2008; Wang and McCallum2006; Watanabe et al. 2011; Wei et al. 2007; Zhang et al. 2010]. A topic model is a hierarchical probabilistic model, in which a document is modeled as a mixture of topics,
and a topic is modeled as a probability distribution over words. Topic models are successfully used in a wide variety of applications including information retrieval [Bleiet al. 2003], collaborative filtering [Hofmann 2003], and visualization [Iwata et al.
2008] as well as the analysis of dynamics.
In this article, we propose a topic model that permits the sequential analysis of the
dynamics of topics with multiple timescales, which we call the Multiscale Dynamic
Topic Model (MDTM), and its efficient online inference procedures. Topics naturally
Authors’ addresses: T. Iwata, Y. Sakurai, and N. Ueda, NTT Communication Science Laboratories; email:iwata.tomoharu@lab.ntt.co,jp; T. Yamada, NTT Science and Core Technology Laboratory Group.
Permission to make digital or hard copies of part or all of this work for personal or classroom use is granted
without fee provided that copies are not made or distributed for profit or commercial advantage and that
copies show this notice on the first page or initial screen of a display along with the full citation. Copyrights
for components of this work owned by others than ACM must be honored. Abstracting with credit is permitted. To copy otherwise, to republish, to post on servers, to redistribute to lists, or to use any componentof this work in other works requires prior specific permission and/or a fee. Permission may be requested
from Publications Dept., ACM, Inc., 2 Penn Plaza, Suite 701, New York, NY 10121-0701, USA, fax +1 (212)869-0481, or permissions@acm.org.
c 2012 ACM 1556-4681/2012/02-ART19 $10.00

DOI 10.1145/2086737.2086739 http://doi.acm.org/10.1145/2086737.2086739

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:2

T. Iwata et al.

evolve with multiple timescales. Let us consider the topic “politics” in a news article
collection as an example. There are some words that appear frequently over many
years, such as “constitution,” “congress,” and “president.” On the other hand, somewords, such as the names of congress members, may appear frequently over periods oftens of years, and other words, such as the names of bills under discussion, may appear
for only a few days. Thus, in MDTM, current topic-specific distributions over words are
assumed to be generated based on estimates of multiple timescale word distributions of
the previous epoch. Using these multiscale priors improves the predictive performanceof the model because the information loss is reduced by considering both the long- andshort-timescale dependency. Furthermore, by using multiple estimates, the varianceof the inference is reduced when compared with models that have a single timescale.
The online inference and parameter estimation processes can be achieved efficiently
based on a stochastic expectation-maximization (EM) algorithm, in which the model issequentially updated using newly obtained data; past data does not need to be stored
and processed to make new inferences. Some topics may exhibit strong long-timescaledependence, and others may exhibit strong short-timescale dependence. Furthermore,
the dependence may change over time. Therefore, we infer these dependencies foreach timescale, for each topic, and for each epoch. By inferring the dependencies fromthe observed data, MDTM can flexibly adapt to topic dynamics. A disadvantage of
online inference is that it can be more unstable than batch inference. With MDTM,
the stability can be improved by smoothing using multiple estimates with different
timescales.
The remainder of this article is organized as follows. In Section 2, we briefly review
related work. In Section 3, we formulate a topic model for multiscale dynamics, anddescribe its online inference procedures. In Section 4, we demonstrate the effectiveness of the proposed method by analyzing the dynamics of real document collections.
Finally, we present concluding remarks and a discussion of future work in Section 5.
2. RELATED WORK
2.1 Topic Modeling

A number of methods for analyzing the evolution of topics in document collections havebeen proposed, such as the dynamic topic model [Blei and Lafferty 2006], topic overtime [Wang and McCallum 2006], online latent Dirichlet allocation [AlSumait et al.
2008], and topic tracking model [Iwata et al. 2009]. However, none of these methods
take account of multiscale dynamics. For example, the dynamic topic model (DTM)[Blei and Lafferty 2006] depends only on the previous epoch distribution. On the other
hand, MDTM depends on multiple distributions with different timescales. Therefore,
with MDTM, we can model multiple timescale dependence, and so infer the current
model more robustly. Moreover, while DTM uses a Gaussian distribution to accountfor the dynamics, the proposed model uses conjugate priors. Therefore, inference in
MDTM is relatively simple compared with that in DTM.
There are two approaches for incorporating dynamics in topic models. The first
approach models dynamics by defining evolution on the hidden variables [Blei and
Lafferty 2006; Nallapati et al. 2007]. The second approach uses the topics in the pastand current epochs to define a prior for future epochs [Ahmed and Xing 2008, 2010;
Blei and Frazier 2010; Gerrish and Blei 2010]. The proposed model uses the second
approach; however it uses the topics for modeling word distributions with multiple
scales.
Recently, online inference algorithms for topic models have been proposed [Canini
et al. 2009; Hoffman et al. 2010; Sato et al. 2010], which can process one document
at a time. These algorithms are standard topic models but adapted to handle massive
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:3

document collections, including those arriving in a stream, and they do not model
temporal information or dynamics, which is our goal. Another difference is that the
proposed inference algorithm is online in the sense that it learns one epoch at a time,
instead of one document at a time. However, it is straightforward to handle the lattercase by incorporating the algorithm described in Hoffman et al. [2010] and Sato et al.
[2010].
2.2 Multiscale Dynamics

The multiscale topic tomography model (MTTM) [Nallapati et al. 2007] can analyze
the evolution of topics at various resolutions of timescales by assuming nonhomogeneous Poisson processes. In contrast, MDTM models the topic evolution within the
Dirichlet-multinomial framework in the same way as most topic models including latent Dirichlet allocation [Blei et al. 2003]. Another advantage of MDTM over MTTM isthat it can make inferences in an online fashion. Therefore, MDTM can greatly reduce
both the computational cost and the memory requirement because past data need not
be stored. Online inference is essential for modeling the dynamics of document collections in which large numbers of documents continue to accumulate at any given
moment, such as news articles and blogs. This is because it is necessary to adapt to
the new data immediately for topic tracking, and it is impractical to prepare sufficient
memory capacity to store all past data.
The multiscale analysis of time-series data such as wavelets is remotely related.
AWSOM [Papadimitriou et al. 2003] is one of the first streaming methods for forecasting and is designed to discover arbitrary periodicities in single time sequences.
Sakurai et al. [2005] proposed BRAID, which builds a multilevel window structure andefficiently detects lag correlations between data streams. Singular value decomposition (SVD) is used for analyzing multiscale patterns in streaming data [Papadimitriou
et al. 2005] as well as topic models. However, since SVD assumes Gaussian noise, it isinappropriate for discrete data such as document collections [Hofmann 1999].
3. PROPOSED METHOD
3.1 Preliminaries

In the proposed model, documents are assumed to be generated sequentially at eachepoch. Suppose we have a set of D t documents at the current epoch, t, and each docuNt,dment is represented by wt,d = {wt,d,n}n=1, i.e. the set of words in the document, where
1 ≤ d ≤ D t. Our notation is summarized in Table I. We assume that epoch t is a discrete variable, and we can set the time period for an epoch arbitrarily at, for example,
one day or one year.
Before introducing the proposed model, we review latent Dirichlet allocation (LDA)[Blei et al. 2003; Griffiths and Steyvers 2004], which forms the basis of the proposedmodel. LDA assumes the following generative process of words in a document. Eachdocument has topic proportions θt,d. For each of the Nt,d words in the document, topic
z t,d,n is chosen from the topic proportions, and then word wt,d,n is generated from a
topic-specific multinomial distribution over words φzt,d,n . Topic proportions θt,d andword distributions φz , are assumed to be generated according to symmetric Dirichlet distributions. Figure 1 (a) shows a graphical model representation of LDA, where
the shaded and unshaded nodes indicate observed and latent variables, respectively.
3.2 Model

We consider a set of multiple timescale distributions over words for each topic to incorporate multiple timescale properties. To account for the influence of the past at
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:4

T. Iwata et al.
Table I. Notation
SymbolDescription

Dt
Nt,d
W
wt,d,nnumber of documents at epoch t
number of words in the dth document at epoch t
number of unique words
nth word in the dth document at epoch t,
wt,d,n ∈ {1, · · · , W}
number of topics
topic of the nth word in the dth document at epoch t, z t,d,n ∈ {1, · · · , Z }number of scales
number of past-periods to consider in scale s
multinomial distribution over topics for the dth document at epoch t,Zθt,d = {θt,d,z }z=1
, θt,d,z ≥ 0, z θt,d,z = 1
multinomial distribution over words for the zth topic at epoch t,W, φt,z,w ≥ 0, w φt,z,w = 1
φt,z = {φt,z,w }w=1
multinomial distribution over words for the zth topic with scale s at epoch t,
 (s)
(s)
(s)
(s)
W
= {ξt,z,w}w=1, ξt,z,w≥ 0, w ξt,z,w
=1
ξt,z

Z
z t,d,n
S
Ls
θt,d
φt,z
(s)
ξt,z

Fig. 1. Graphical models of, (a) latent Dirichlet allocation, and (b) the multiscale dynamic topic model.

different timescales on the current epoch, we assume that current topic-specific word
distributions φt,z are generated according to the multiscale word distributions at the
(s)
(s)
(s)
S
W
previous epoch {ξt−1,z}s=1. The multiscale word distribution ξt−1,z= {ξt−1,z,w}w=1
represents a distribution over words of topic z with scale s at epoch t − 1, and is defined asfollows.
tt =t−L +1 Nt ,z,w
(s)
,
(1)
ξt,z,w =  t swt =t−L s +1 Nt ,z,w
where L s is the length of scale s and Nt,z,w is the number of times word w was assigned to topic z at epoch t. In particular, we use the following asymmetric Dirichlet
distribution for the prior of the current word distribution φt,z , in which the Dirichlet
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:5

(s)
Fig. 2. Multiscale word distributions at epoch t with S = 4. Each histogram shows ξt−1,z, which is a multinomial distribution over words with timescale s.

parameter is defined such that the mean of φt,z becomes proportional to the weighted
sum of multiscale word distributions at the previous epoch,
 S


(s)
φt,z ∼ Dirichlet
λt,z,sξt−1,z ,(2)
s=0

where λt,z,s is a weight for scale s in topic z at epoch t, and λt,z,s > 0. By estimatingSfor each epoch, for each topic, and for each timescale using the curweights {λt,z,s}s=0rent data as described in Section 3.3, MDTM can flexibly respond to the influence ofthe previous short- and long-timescale distributions on the current distribution. In this(s)Sway, the estimated multiscale word distributions {ξt−1,z}s=1
at the previous epoch constitute the hyperparameters in the current epoch. Their estimation will be describedin Section 3.4.
There are many different ways of setting the scales, but for generality and sim(s)indicates the word distribution from t − 2s−1 + 1 to t,
plicity, we set them so that ξt,z
(s=1)
is equivalent to the estimate
where a larger s represents longer timescale, and ξt,z(s=0)
= W −1of a unit-time word distribution φt,z . We use a uniform word distribution ξt,z,w
for scale s = 0. This uniform distribution is used to avoid the zero probability problem. Figure 2 shows multiscale word distributions in this setting. Multiscale word(s)
distributions ξt,z
are likely to become smoother as the timescale increases and exhibit
more peaks as the timescale decreases. By using the information presented in thesevarious timescales as the prior for the current distribution with weights, we can infer
the current distribution more robustly. In Figure 2, we use timescales of a geometricprogression: 1, 2, 4, . . . , 2 S−1 . Instead of using 2s−1 epochs for scale s, we could use any
number of epochs. For example, if we know that the given data exhibit a periodicity
e.g., of one week and one month, we can use the scale of one week for s = 1 and one
month for s = 2. In such a case, we can still estimate the parameters in a similar way
to that used with the algorithm described in Section 3.4. Typically, we do not know
the periodicity of given data in advance, we therefore consider a simple scale settingin this article.
In LDA, topic proportions θt,d are sampled from a Dirichlet distribution. To capture the dynamics of topic proportions with MDTM, we assume that the Dirichlet paZdepend on the previous parameters. In particular, we use therameters αt = {αt,z }z=1
following Gamma prior for a Dirichlet parameter of topic z at epoch t,
αt,z ∼ Gamma(γ αt−1,z , γ ),(3)

where the mean is αt−1,z , and the variance is αt−1,z /γ . By using this prior, the meanis the same as that at the previous epoch unless otherwise indicated by the new data.
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:6

T. Iwata et al.

Parameter γ controls the temporal consistency of the topic proportion prior αt,z . Whenγ is high, αt,z is likely to be close to αt−1,z .
Assuming that we have already calculated the multiscale parameters at epoch
(s)
S
Z
Z
}s=0
}z=1
and αt−1 = {αt−1,z }z=1
, and given parameters, γ andt − 1, t−1 = {{ξt−1,zSZt = {{λt,z,s}s=0 }z=1 , MDTM is characterized by the following process for a set of docDt
uments Wt = {wt,d}d=1
at epoch t.
(1) For each topic z = 1, · · · , Z :(a) Draw topic proportion prior
αt,z ∼ Gamma(γ αt−1,z , γ ),(b) Draw word distribution(s)),
φt,z ∼ Dirichlet( s λt,z,sξt−1,z(2) For each document d = 1, · · · , D t:(a) Draw topic proportions
θt,d ∼ Dirichlet(αt),(b) For each word n = 1, · · · , Nt,d:i.Draw topic
z t,d,n ∼ Multinomial(θt,d),ii.Draw word
wt,d,n ∼ Multinomial(φt,zt,d,n ).
Figure 1(b) shows a graphical model representation of MDTM.
3.3 Online Inference

We present an online inference algorithm for MDTM, that sequentially updates themodel at each epoch using the newly obtained document set and the multiscale modelof the previous epoch. The information in the data up to, and including, the previous
epoch is aggregated into the previous multiscale model. The online inference and parameter estimation can be efficiently achieved by a stochastic EM algorithm [Andrieu
et al. 2003], in which the collapsed Gibbs sampling of latent topics, Zt, and the maximum a posteriori (MAP) estimation of hyperparameters, αt and t, are alternately
performed.
We assume a set of documents Wt at current epoch t, and estimates of parametersfrom the previous epoch αt−1 , t−1 , γ and t are given. The joint distribution on the setof documents, the set of topics, and the topic proportion priors given the parameters
are defined as follows.
P(Wt, Zt, αt |αt−1 , γ , t−1 , t) = P(αt|αt−1 , γ )P(Zt|αt)P(Wt|Zt, t−1 , t),Nt,d D t}d=1 represents a set of topics.
where Zt = {{z t,d,n}n=1(4)

The first term on the right handside of (4) is as follows using (3).

P(αt|αt−1 , γ ) =γ αt−1,z −1
 γ γ αt−1,z αt,zexp(−γ αt,z )(γ αt−1,z )

z

,

(5)

where (·) is the gamma function. We can integrate out the multinomial distribuDt
Z
tion parameters in MDTM, {θt,d}d=1
and {φt,z }z=1
, by taking advantage of Dirichletmultinomial conjugacy. The second term is calculated by
P(Zt|αt) =Nt,dDt  

d=1

P(z t,d,n|θt,d)P(θt,d|αt)dθt,d,n=1

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

(6)Sequential Modeling of Topic Dynamics with Multiple Timescales
Dtand we have the following equation by integrating out {θt,d}d=1
.
 
Dt

 z  Nt,d,z + αt,z
z αt,z
P(Zt|αt) =, Nt,d + z αt,zz  αt,z
d

19:7

(7)

where Nt,d,z is the number of words in the dth documentZ assigned to topic z at epocht, and Nt,d = z Nt,d,z . Similarly, by integrating out {φt,z }z=1
, the third term is given as
follows.
(s)
w  Nt,z,w +s λt,z,sξt−1,z,w

s λt,z,sP(Wt|Zt, t−1 ,  t) =
,
(8)


(s)
 Nt,z + s λt,z,s

λ
ξ
z
t,z,s
w
s
t−1,z,w
where
Nt,z,w is the number of times word w was assigned to topic z at epoch t, andNt,z = w Nt,z,w .
The inference of the latent topics Zt can be efficiently computed by using collapsed
Gibbs sampling. Let j = (t, d, n) for notational convenience, and z j be the assignmentof a latent topic to the nth word in the dth document at epoch t. Then, given thecurrent state of all but one variable z j, a new value for z j is sampled from the following
probability.
(s)
Nt,d,k\ j + αt,k Nt,k,w j\ j + s λt,k,sξt−1,k,w jP(z j = k|Wt, Zt\ j, αt, t−1 , t) ∝,(9)Nt,d\ j + z αt,zNt,k\ j + s λt,k,s
where \ j represents the count yielded by excluding the nth word in the dth document.
See Appendix A for the derivation.
The parameters αt and t are estimated by maximizing the joint distribution (4).
The fixed-point iteration method described in Minka [2000] can be used for maximizing
the joint distribution as follows.
 old old
old
Nt,d,z + αt,z−αt,zγ αt−1,z − 1 + αt,zd
αt,z ←
 old 
,
(10)

 old
Nt,d + z αt,z −γ+ dz  αt,z (·) is a digamma function defined by (x) = ∂ log∂ x(x) , and old (s ) (s) old (s )
N
−
ξ
+
λ
ξ
λ
ξ




t,z,w
w t−1,z,w
s t,z,s t−1,z,w
s t,z,s t−1,z,w
 old.

λt,z,s ← λold
 old
t,z,s
Nt,z + s λt,z,s −s λt,z,s

where

(11)

See Appendix B and Appendix C for the derivation. By iterating Gibbs sampling with
(9) and maximum likelihood estimation with (10) and (11), we can infer latent topics
while optimizing the parameters. Since MDTM uses the past distributions as thecurrent prior, the label switching problem [Stephens 2000] is not likely to occur whenthe estimated λt,z,s is high, which implies that current topics strongly depend on theprevious distributions. Label switching can occur when the estimated λt,z,s is low. By
allowing low λt,z,s, which is estimated from the given data at each epoch and each topic,
MDTM can adapt flexibly to changes even if existing topics disappear and new topics
appear in midstream.
The time complexity of one iteration of our Gibbs sampling is O(D t Nt Z ), where Nt
is the average number of words in documents at epoch t, and it does not depend on thenumber of scales. The time complexity of one iteration of hyperparameter estimation
increases linearly with the number of scales.
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:8

T. Iwata et al.

3.4 Efﬁcient Estimation of Multiscale Word Distributions

By using the topic assignments obtained after iterating the stochastic EM algorithm,(s)represents the probability
we can estimate multiscale word distributions. Since ξt,z,w
s−1
of word w in topic z from t − 2 + 1 to t, the estimation is as follows.
t(s)N̂t,z,wt =t−2s−1 +1 N̂t ,z,w(s)
ξt,z,w =  (s) =  
,
(12)
t
w
t =t−2s−1 +1 N̂t ,z,ww N̂t,z,w(s)where N̂t,z,w is the expected number of times word w was assigned to topic z from
epochs t − 2s + 1 to t, and N̂t,z,w is the expected number of times word w is assigned
to topic z at epoch t. The expected number is calculated by N̂t,z,w = Nt,z φ̂t,z,w , where
φ̂t,z,w is a point estimate of the probability of word w in topic z at epoch t. Although we
integrate out φt,z,w , we can recover its point estimate as follows.
(s)Nt,z,w + s λt,z,sξt−1,z,w
.
(13)
φ̂t,z,w =Nt,z + s λt,z,s
While we can use the actual number of times Nt,z,w in (12), we use the expected number(s=1)
to be the estimate of φt,z,w as follows.
of times N̂t,z,w to constrain the estimate of ξt,z,wN̂t,z,w(s=1)
ξt,z,w== φ̂t,z,w .
w N̂t,z,w

(14)

(s)

Note that the value N̂t,z,w can be updated sequentially from the previous value(s)N̂t−1,z,w as follows.
(s)

(s)

N̂t,z,w ← N̂t−1,z,w + N̂t,z,w − N̂t−2s−1 ,z,w .

(15)

(s)

Therefore, N̂t,z,w can be efficiently updated through just two additions rather than 2s−1
additions. We refer to it as the naive updating method.
L EMMA 3.1. The naive updating method requires O(2 SZ W) memory for updating
multiscale word distributions.
(s)

P ROOF. The naive updating method subtracts N̂t−2S−1 ,z,w for updating N̂t,z,w . Therefore, it has to maintain O(2 S) values of N̂t−2S−1 ,z,w from t − 2 S−1 to t − 1. We need theO(2 S) values for each of Z topics and each of W words. Thus, it would take O(2 SZ W)
space.
Since the memory requirement increases exponentially with the number of scales
in the naive updating method, this requirement prevents us from modeling longtimescale dynamics. Thus, we consider approximating the update by decreasing the
update frequency for long-timescale distributions as in Algorithm 1, which is linearwith respect to the number of scales.
L EMMA 3.2. The proposed updating method requires O(SZ W) memory for updating
multiscale word distributions.
(s)

P ROOF. The proposed updating method can update N̂t,z,w using the previous value
(s−1)

(s−1)

at scale s − 1, N̂t−1,z,w , and the current value at scale s − 1, N̂t,z,w . When we update
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:9

(s)

ALGORITHM 1: Algorithm for the approximate update of N̂t,z,w .
(1)

N̂t,z,w ← N̂t,z,w
for s = 2, · · · , S do
if t mod 2s−1 = 0 then
(s)

(s−1)

(s)

(s)

(s−1)

N̂t,z,w ← N̂t,z,w + N̂t−1,z,w
else
N̂t,z,w ← N̂t−1,z,w
end
end

Fig. 3. Graphical models of the online multiscale dynamic topic model.
(s−1)

from s = 1 to s = S, we must store only O(S) values of N̂t−1,z,w at the previous epoch.
We need the O(S) values for each of Z topics and each of W words. Thus, it would takeO(SZ W) space.
(s)

Figure 4 shows approximate updating of N̂t,z,w with S = 3 from t = 4 to t = 8. Eachbox represents N̂t ,z,w , where the number represents epoch t . The boxes in the bottom
row indicate that the N̂t,z,w value is newly calculated at epoch t. Each row at each epoch(s)represents N̂t,z,w as sum of the corresponding boxes, and the shaded boxes indicate that(s)the values are updated from the previous values. N̂t,z,w is updated every 2s−1 epoch.
(s)

For example, value N̂t,z,w at scale s = 1 is updated every epoch, the value at scale s = 2
is updated every two epochs, and the value at scale s = 3 is updated every four epochs.
The value at scale s = 1 can be updated by replacing it with the newly calculated value.
The value at scale s = 2 can be updated by adding the newly calculated value to theprevious value at scale s = 1; see epochs t = 6 or t = 8 in Figure 4. The value at scales = 3 can be updated by summing the newly calculated value, the previous value atscale s = 1 and the previous value at scale s = 2; see epochs t = 8 in Figure 4. As seen in
the preceding, we can calculate the value by using values at the previous epoch withthis approximate updating; we do not need to store the past values before the previous
epoch.
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:10T. Iwata et al.

(s)Fig. 4. Approximate updating of N̂t,z,w from t = 4 to t = 8 with S = 3.
(s)buffer

ALGORITHM 2: Algorithm for the approximate updating of N̂t,z,w using a buffer. N̂t,z,w represents the buffer value.
bufferbufferN̂t,z,w ← N̂t−1,z,w + N̂t,z,w
if (t − 1) mod L s = 0 then(s)bufferN̂t,z,w ← N̂t,z,w
bufferN̂t,z,w ← N̂t,z,w
end

Since the dynamics of a word distribution for a long-timescale are considered tobe slower than that for a short-timescale, this approximation—decreasing the update
frequency for long-timescale distributions—is reasonable. Figure 3 shows a graphicalmodel representation of online inference in MDTM.
For the Dirichlet prior parameter of the word distribution, we use the weighted sumof the multiscale word distributions, as in (2). The parameter can be rewritten as theweighted sum of the word distributions for each epoch as follows.
S(s)λt,z,sξt−1,z,w
=

t−1


λt,z,t φ̂t ,z,w ,(16)t =t−2S−1

s=1

where
N̂t,z,w,
φ̂t,z,w = 
w  N̂t,z,w 
is the expected probability that word w appears in topic z at epoch t, andS
λt,z,s w N̂t ,z,w,
λt,z,t = t−1
 =t−2s−1 N̂t ,z,wwt
s= log (t−t )+1

(17)

(18)

2

is its weight. See Appendix D for the derivation. Therefore, the multiscale dynamic
topic model can be seen as an approximation of a model that depends on the word
distributions for each of the previous epochs. By considering multiscale word distributions, the number of weight parameters t can be reduced from O(2 S Z ) to O(SZ ), andthis leads to more robust inference. Furthermore, the use of multiscaling also reduces
the memory requirement from O(2 S Z W) to O(SZ W) as previously described.
3.5 Enhanced Methods for Approximate Updating
Approximate updating method for arbitrary scale. In the preceding discussion, we set scales(s)
indicates the word distribution from t − 2s−1 + 1 to t. However, we canso that ξt,z
set scales arbitrarily. If we know that the given data exhibit periodicity, e.g., of one
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:11

(s)

Fig. 5. Approximate updating of N̂t,z,w with scale size L s = 3.
(s)b

ALGORITHM 3: Algorithm for the approximate updating of N̂t,z,w using B buffers. N̂t,z,w represents the bth buffer value, tb represents the timing with which the bth buffer is updated, and0 ≤ tb < L s.
for b = 1, · · · , B dobbN̂t,z,w ← N̂t−1,z,w + N̂t,z,w
if (t − tb ) mod L s = 0 then(s)bN̂t,z,w ← N̂t,z,wbN̂t,z,w ← N̂t,z,w

end
end

week and one month, we can use the scale of one week for s = 1 and one month fors = 2. In these cases, we can also efficiently update values by using buffers for storing
intermediate values. Algorithm 2 shows the updating algorithm with scale size L susing a buffer.
L EMMA 3.3. The proposed updating method (arbitrary scale) requires O(SZ W)
memory for updating multiscale word distributions.
(s)

P ROOF. The proposed method for an arbitrary scale can update N̂t,z,w using a single
buffer. We need a buffer for each of S scales, each of Z topics, and each of W words.
Thus, it would take O(SZ W) space.
While this approximation needs a buffer, the memory requirement of this method
is still O(SZ W), which is also linear with respect to the number of scales. Figure 5(s)shows the approximate updating of N̂t,z,w with scale size L s = 3. At epoch t = 3, thenewly calculated value N̂t=3,z,w is stored in the buffer. At epoch t = 4, the buffer value(s)is updated by adding the newly calculated value N̂t=4,z,w . At epoch t = 5, value N̂t,z,w
is updated by using the newly calculated value and the buffer value at the previous
epoch. The buffer is cleared and the newly calculated value N̂t=5,z,w is stored in the(s)buffer. In this way, N̂t,z,w is updated every L s − 1 epochs by using the value at theprevious epoch.
Improving approximation by using multiple buffers. We can improve the approximation byincreasing the update frequency. By using B buffers, the value is updated B times
every L s − 1 epochs, and the rate of unapproximated values is L sB−1 . Algorithm 3 shows
the updating algorithm with scale size L s using B buffers.

L EMMA 3.4. The proposed updating method (arbitrary frequency) requires
O(BSZW) memory for updating multiscale word distributions.
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:12

T. Iwata et al.

(s)Fig. 6. Approximate updating of N̂t,z,w with scale size L s = 4 using two buffers.
Table II.
Average perplexities over epochs. The value in the parenthesis represents the standard deviation
over datasets.
NIPS
PNAS
Digg
Addresses

MDTM

DTM

LDAallLDAoneLDAonline

1754.9 (41.3)
2964.3 (122.0)
3388.9 (37.7)
1968.8 (56.5)

1771.6 (37.2)
3105.7 (146.8)
3594.2 (46.4)
2105.2 (49.7)

1802.4 (36.4)
3262.9 (159.7)
3652.6 (27.1)
2217.2 (75.3)

1822.0 (44.0)
5221.5 (268.7)
5162.9 (43.4)
3033.5 (70.9)

1769.8 (41.5)
3401.7 (149.1)
3500.0 (43.6)
2251.6 (62.0)(s)

P ROOF. The proposed method for an arbitrary frequency can update N̂t,z,w using thevalues in B buffers. We need B buffers for each of S scales, each of Z topics, and eachof W words. Thus, it would take O(BSZ W) space.
By using more buffers, we can approximate the values more closely. This method
allows users to choose the update frequency. Figure 6 shows approximate updatingwith two buffers and scale size L s = 4. The value at each buffer is updated by adding
the newly calculated value. When the size of the buffer becomes the scale size L s − 1,(s)the value is used for calculating N̂t,z,w , and the buffer is updated with the value N̂t,z,w .
For example, at epoch t = 5, the first buffer is updated by adding N̂t=5,z,w . At epoch(s)t = 6, it is used for calculating N̂t=6,z,w by adding it to N̂t=6,z,w , and the buffer is updatedwith N̂t=6,z,w .
4. EXPERIMENTS
4.1 Setting

We evaluated the multiscale dynamic topic model with online inference (MDTM) using
four real document collections with timestamps: NIPS, PNAS, Digg, and Addresses.
The NIPS data consist of papers from Neural Information Processing Systems(NIPS) conferences from 1987 to 1999. There were 1740 documents, and the vocabulary size was 14,036. The unit epoch was set at one year, so there were 13 epochs.
The PNAS data consist of the titles of papers that have appeared in the Proceedings of
the National Academy of Sciences from 1915 to 2005. There were 79,477 documents,
and the vocabulary size was 20,534. The unit epoch was set at one year, so there were
91 epochs. The Digg data consist of blog posts that have appeared on the social news
Web site Digg1 from January 29th to February 20th 2009. There were 108,356 documents, and the vocabulary size was 23,494. The unit epoch was set at one day, so there
were 23 epochs. The Addresses data consist of the State of the Union addresses from1

http://digg.com

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:13

Fig. 7. Perplexities for each epoch.

1790 to 2002. We increased the number of documents by splitting each transcript into
3-paragraph documents as done in Wang and McCallum [2006]. We omitted words
that occurred in fewer than 10 documents. There were 6413 documents, and the vocabulary size was 6759. The unit epoch was set at one year, and excluding the years
for which data were missing, there were 205 epochs. We omitted stop-words from alldata sets.
We compared MDTM with DTM, LDAall, LDAone, and LDAonline. DTM is a dynamic topic model with online inference, which does not take multiscale distributions
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:14

T. Iwata et al.

Fig. 8. Average perplexities with different numbers of topics and their standard deviations.

into consideration; it corresponds to MDTM with S = 1. Note that the DTM used here
models dynamics with Dirichlet priors while the original DTM used Gaussian priors.
LDAall, LDAone, and LDAonline are based on LDA, and therefore do not model the
dynamics. LDAall is an LDA model, which uses all past data for inference. LDAone is
an LDA model, which uses just the current data for inference. LDAonline is an online
learning extension of LDA, in which the parameters are estimated using those of theprevious epoch and the new data [Banerjee and Basu 2007]. The time complexity of
Gibbs sampling for LDAone and LDAonline are O(D t Nt Z ), which is thesame as thatfor the proposed model. That for LDAall is O(D 1,t Nt Z ), where D 1,t = tτ =1 D τ . For afair comparison, the hyperparameters in these LDAs were optimized using stochastic
EM. We set the number of latent topics at Z = 50 for all models, and we iterated an
E- and M-step 500 times for each epoch for all models. In MDTM, we used γ = 1, andwe estimated the Dirichlet prior for topic proportions subject to αt,z ≥ 10−2 to avoid
overfitting. We set the scale size at L s = 2s−1 , and the number of scales so that one of
the multiscale distributions covered the entire period, or S = log2 T + 1 , where T is
the number of epochs. We used Algorithm 1 for updating the parameters.
We evaluated the predictive performance of each model using the perplexity of heldout word.
⎛   test⎞
 test
Nt,d
log
P
w
|t,
d,
D
t
d
n=1
t,d,n
⎠,
Perplexity = exp ⎝−(19)
 test
N
d t,d
test
test
is the number of held-out words in the dth document at epoch t, wt,d,n
is thewhere Nt,dnth held-out words in the document, and Dt represents training samples until epoch

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:15

Fig. 9. Average perplexities of MDTM with different numbers of scales and their standard deviations.

t. A lower perplexity represents higher predictive performance. The word probability
can be calculated as follows,P(w|t, d, Dt ) =θ̂t,d,z φ̂t,z,w ,
(20)
z

where
θ̂t,d,z =Nt,d,z + αt,z.
Nt,d + z αt,z(21)

We used half of the words in 10% of the documents as held-out words for each epoch,
and used the other words as training samples. We created ten sets of training and test
data by random sampling, and evaluated the average perplexity over the ten data sets.
4.2 Results

The average perplexities over the epochs are shown in Table II, and the perplexities
for each epoch are shown in Figure 7. For all data sets, MDTM achieved the lowest
perplexity, which implies that MDTM can appropriately model the dynamics of various
types of data sets through its use of multiscale properties. DTM had a higher perplexity than MDTM because it could not model the long-timescale dependencies. LDAalland LDAonline have high perplexities because they do not consider the dynamics. The
perplexity achieved by LDAone is high because it uses only current data and ignores
past information. The average perplexities over epochs with different numbers of topics are shown in Figure 8. Under the same number of topics, MDTM achieved thelowest perplexities in every case except when Z = 150 and 200 in the NIPS data,
and achieved the second lowest perplexities when Z = 150 and 200 in the NIPS data.
LDAone achieved the lowest perplexity with Z = 150 and 200 in the NIPS data because the number of words in each document in the NIPS data is large compared withother data, and the model can be learned adequately by using data at each epoch. The
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:16

T. Iwata et al.

Fig. 10. Average perplexity of MDTM over iterations on the inference for each epoch in NIPS data.

average perplexities over epochs with different numbers of scales in MDTM are shownin Figure 9. Note that s = 0 uses only the uniform distribution, while s = 1 uses theuniform distribution and the previous epoch’s distribution. The perplexities decreased
as the number of scales increased. This result indicates the importance of considering
multiscale distributions. Although the number of parameters increases as the numberof scales increases, overfitting did not occur. This result suggests that we can use as
many scales as we wish as long as we can bear the computational and memory costs.
Figure 10 shows the perplexities of the proposed model over iterations on inferencein NIPS data, which are averaged over ten different initializations. As the number of
iterations increases, the perplexity decreases, and eventually converges to a certain
point.
Figure 11 shows the average computational time per epoch when using a computerwith a Xeon5355 2.66GHz CPU and a 16GB memory. The computational time for
MDTM is roughly linear against the number of scales. Even though MDTM considers
multiple timescale distributions, its computational time is much smaller than thatof LDAall, which considers a single timescale distribution. This is because MDTM
only uses current samples for inference. In contrast, LDAall uses all the samples for
inference.
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:17

Fig. 11. Average computational time (sec) of MDTM per epoch with different numbers of scales, LDAall,
LDAone, and LDAonline and their standard deviations.

When the data in an epoch do not depend on data in other epochs, LDAone is anappropriate model. When the data for all epochs are generated from a distribution,or the data have no dynamics, LDAall and LDAonline are both appropriate models.
MDTM can flexibly model data that have dynamics, where data in consecutive epochs
exhibit some dependencies. In particular, MDTM with λt,z,s = 0 for s > 0 corresponds
to LDAone, and MDTM with λt,z,s = 0 for all s except when the timescale is infinite,
corresponds to LDAonline.
Figure 12 shows the average estimated λt,z,s with different numbers of scales s in
MDTM. The sum of the values for each epoch and for each topic are normalized to one.
The parameters decrease as the timescale lengthens. This result implies that recent
distributions are more informative as regards estimating current distributions, whichis intuitively reasonable.
In Section 3.5, we described an algorithm for improving the approximation of an
estimation by using multiple buffers. We experimentally investigated how the numberof buffers affects the predictive performance. Figure 13 shows the average perplexities
achieved by MDTM when the number of topics Z = 50, the number of scales S = 1,
and the scale size L s = 10. As the number of buffers increased, the perplexity decreased, although the memory requirement increased. This result suggests that we
can improve predictive performance by using as many buffers as possible.
4.3 Application of MDTM

MDTM can extract the topic evolution in multiple timescales. For experiments inthis section, we used MDTM, where the number of topics Z = 50 and the numberof scales S = 4. Figure 14 shows two topic examples of the multiscale topic evolution in NIPS data analyzed by MDTM. Note that we omit words that appeared in thelonger timescales from the table. In the longest timescale, basic words for the research
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:18T. Iwata et al.

Fig. 12. Average normalized weights λ with different scales estimated in MDTM and their standard
deviations.

Fig. 13. Average perplexities with different numbers of buffers in MDTM with multiple buffers and theirstandard deviations.

field are appropriately extracted, such as “speech,” “recognition,” and “speaker” in the
speech recognition topic, “control,” “action,” “policy,” and “reinforcement” in the reinforcement learning topic. In the shorter timescale, we can see the evolution of trends
in the research. For example, in speech recognition research, phoneme classification
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:19

Fig. 14. Two topic examples of the multiscale topic evolution in NIPS data analyzed by MDTM: (a) speech
recognition, and (b) reinforcement learning topics. The ten most probable words for each epoch, timescale,and topic are shown.

was a common task until 1990, and since 1991 probabilistic approaches such as hidden
Markov models (HMM) have been frequently used.
MDTM can also be used for tracking the popularity of each topic. Figure 15 showsthe estimated αt,z of four topics from the NIPS data analyzed by MDTM, which represents the popularity dynamics for each topic. The most probable words on the longtimescale for topics “neural network” and “probabilistic model” are shown in Table III.
ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:20T. Iwata et al.
Table III. Two Topic Examples in NIPS Data Analyzed by MDTM
neural network

layer, net, neural, nets, system, figure,activation, signal, error, noise, backpropagation

probabilistic model

data, model, training, error, probability,
parameters, set, distribution, models, Gaussian, variance

Fig. 15. Popularity dynamics of four topics in NIPS data analyzed by MDTM.

We can see that the popularity of the topics “neural network” has decreased over time,on the other hand, that of the “probabilistic model” has increased.
The estimated parameters λt,z,t−1 in (18) represent the dependency of the word distribution of topic z between t and t− 1. Therefore, we can analyze the magnitude of thechange in the topic by using λt,z,t−1 . Usually, the change in a topic can be calculated byusing the KL divergence between the current word distribution and the previous word
distribution KL(φt−1,z ||φt,z ). The correlation coefficients between 1/ log(λt,z,t−1 ) and the
KL divergences were 0.625, 0.216, 0.729, and 0.152 for NIPS, PNAS, Digg, and Addresses data sets respectively, which were calculated using all time intervals and all
topics. This result shows that we can analyze the magnitude of changes by using λt,z,t−1
without calculating the KL divergence.
Figure 16 shows the estimated 1/ log(λt,z,t−1 ) of four topics in the NIPS data analyzed
by MDTM, where this value represents how the estimated topics change comparedwith those in the previous epoch. In general, the values in late epochs are low, since
the inference of topics becomes stable when we use more data for the inference. The
value for a “neural network” in 1995 was high because it was at the end of the neuralnetwork boom.
5. CONCLUSION

In this article, we have proposed a topic model with multiscale dynamics andefficient online inference procedures. We have experimentally confirmed that the
proposed method can appropriately model the dynamics in document data by considering multiscale properties, and that it is computationally efficient.
We assumed that the number of topics was known and fixed over time. We can
automatically infer the number of topics by extending the model to a nonparametric
Bayesian model such as the Dirichlet process mixture model [Ren et al. 2008; Tehet al. 2006]. In future work, we plan to determine the unit time interval, the length ofscale, and the number of scales, automatically from the given data. The problem of theunit time interval can be solved by using continuous time dynamic topic models [Wang

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:21

Fig. 16. Topic change dynamics of four topics in NIPS data analyzed by MDTM.

et al. 2008]. The length of scale and the number of scales can be determined using
nonparametric Bayesian methods. The length of scale can be continuous as in [Wangand McCallum 2006], which treats time dependence as a continuous beta random variable. The proposed model can be extended to find influential documents by modelingthe influence for each document as described in Gerrish and Blei [2010]. Since theproposed method is applicable to various kinds of discrete data with timestamps, suchas Web access logs, blogs, and e-mail, we will evaluate the proposed method further byapplying it to other data sets.
APPENDIX
A. DERIVATION OF (9)In this appendix, we give the derivation of (9).
P(z j = k|Wt, Zt\ j, αt, t−1 , t)∝ P(z j = k, w j|Wt\ j, Zt\ j, αt, t−1 , t)= P(z j = k|Zt\ j, αt )P(w j|Wt\ j, z j = k, Zt\ j, t−1 , t).

(22)

The first factor of (22) becomes,P(z j = k|Zt\ j, αt)=P(z j = k|θt,d)P(θt,d|Zt\ j, αt )dθt,d Nt,d,z\ j+αz −1z Nt,d,z\ j + αt,z=
θt,d,kθt,d,zdθt,dz  Nt,d,z\ j + αt,zz  Nt,d,z\ j+αz −1 Nt,d,k\ j+αkz Nt,d,z\ j + αt,z=
θt,d,zθt,d,kdθt,dz  Nt,d,z\ j + αt,z
z =k




z =k  Nt,d,z\ j + αt,z  Nt,d,k\ j + αt,k + 1
z [Nt,d,z\ j + αt,z ]



=

z  Nt,d,z\ j + αt,zz Nt,d,z\ j + αt,z + 1
Nt,d,k\ j + αt,k,=
Nt,d\ j + z αt,z(23)

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:22

T. Iwata et al.

αz −1z (αz )where we useddθ = (in the fourth equation, which is the normalizingz θzz αz )
constant of the Dirichlet distribution, and (x + 1) = x(x) in the fifth equation. In asimilar way, the second factor of (22) becomes,
P(w j|Wt\ j, z j = k, Zt\ j, t−1 , t)=P(w j|z j = k, φt,k )P(φt,k |Wt\ j, Zt\ j, t−1 , t)dφt,k



(s)
(s)
( w [Nt,k,w\ j + s λt,k,sξt−1,k,w])  Nt,k,w\ j+s λt,k,sξt−1,k,w−1=
φt,k,w j
φ
dφt,k

t,k,w
(s)
w (Nt,k,w\ j +s λt,k,s ξt−1,k,w ) w


 
(s)

w Nt,k,w\ j +s λt,k,s ξt−1,k,w=
(s)w  Nt,k,w\ j +s λt,k,s ξt−1,k,w

 

(s)
(s)
Nt,k,w\ j+ s λt,k,sξt−1,k,w−1 Nt,k,w j\ j+ s λt,k,sξt−1,k,w j
φt,k,w
φt,k,w j
dφt,k×w =w j

=



 
w

Nt,k,w\ j +w

Nt,k,w\ j +




(s)
s λt,k,s ξt−1,k,w

(s)s λt,k,s ξt−1,k,w



(s)
(s)
 Nt,k,w j\ j + s λt,k,sξt−1,k,w
Nt,k,w\ j + s λt,k,sξt−1,k,w
+1
j


×


(s)

w Nt,k,w\ j +s λt,k,s ξt−1,k,w + 1(s)
Nt,k,w j\ j + s λt,k,sξt−1,k,w
j
=.
Nt,k\ j + s λt,k,sw =w j (24)

B. DERIVATION OF (10)In this appendix, we give the derivation of (10). Consider the following function.
f (x) = log

(x),(n + x)(25)

where x > 0 and n > 0. The tangent of f (x) at x̂ is as follows.
g(x) = log

(x̂)+ (x̂) −(n + x̂)(n + x̂) (x − x̂).

(26)

Since f (x) is convex downward, the tangent is the lower bound, f (x) ≥ g(x), and we get
the following inequality.
log


(x)
(x̂)
≥ log
+ (x̂) −(n + x)(n + x̂)(n + x̂) (x − x̂).

(27)

Additionally, consider the following function.
(n + x)(n + x̂)x
− log
− x̂[ (n + x̂) − (x̂)] log
(x)
(x̂)
x̂

n−1 

x
x̂
log
,
log(m + x) − log(m + x̂) −=m + x̂
x̂

h(x) = log

m=0

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

(28)

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:23

n−1
where x > 0, x̂ > 0, n > 0, and we used log (n+x)= m=0log(m+ x) and (n+ x) − (x) =(x)n−1 1m=0 m+x in the second equality. The value in each summation is nonnegative.
k(x) = log(m + x) − log(m + x̂) −x̂xlog ≥ 0,m + x̂
x̂

(29)

becausem(x − x̂)
∂k(x)
=
,
∂x
(m + x)x(m + x̂)(30)

and therefore k(x) monotonically decreases (increases) when x < x̂ (x > x̂), and k(x̂) = 0.
By using (28) and (29), we can obtain the following inequality.
log

(n + x̂)
(n + x)≥ log
+ x̂[ (n + x̂) −
(x)
(x̂)

x
(x̂)] log .
x̂

(31)

The joint log likelihood can be written as follows, using (4), (5), (6), and (7).
log P(Wt, Zt, αt|αt−1 , γ , t−1 , t)= log P(αt|αt−1 , γ ) + log P(Zt|αt) + C= (γ αt−1,z − 1) log αt,z − γ αt,z








log 
+
αt,z − log  Nt,d +
αt,z
d

+





z

z



log  Nt,d,z + αt,z − log  αt,z + C ,

(32)

d

where C and C are the constants that do not depend on αt,z . The second term of (32)satisfies the following inequality.






αt,z ) − log  Nt,d +αt,zlog (
d

≥



z



log 

d

+















z

old
αt,z
− log  Nt,d +zNt,d +




old
αt,z

−




z

d




old
αt,z

z



old
αt,z

z

(33)



old
αt,z
− αt,z

,

z


 old
where we used (27) by letting x = z αt,z , n = Nt,z and x̂ = z αt,z. The third term of(32) satisfies the following inequality.
[log (Nt,d,z + αt,z ) − log (αt,z )]
d

≥


old
old
− log  αt,zlog  Nt,d,z + αt,z
d



+

old
Nt,d,z + αt,z

−

old
αt,z



old
αt,z

old
log αt,z − log αt,z


,

(34)

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:24

T. Iwata et al.

oldwhere we used (31) by letting x = αt,z , n = Nt,d,z , and x̂ = αt,z
. By using (33) and (34),
we can get the following lower bound of the joint log likelihood.

log P(Wt, Zt, αt|αt−1 , γ , t−1 , t)≥ (γ αt−1,z − 1) log αt,z − γ αt,z − αt,z
+






Nt,d +

old
Nt,d,z + αt,z

old
αt,z

−




−

z

d




old
αt,z




old
αt,z

z

old
αt,z
log αt,z + C ≡ F(αt,z ),

(35)

d

where C is a constant that does not depend on αt,z . The derivative of F(αt,z ) is asfollows.
 






1
∂ F(αt,z )oldold= (γ αt−1,z − 1)−γ −αt,z −
αt,zNt,d +
∂αt,z
αt,z
z
z
d


old
old
old 1
+
.
(36)
Nt,d,z + αt,z
−
αt,z
αt,z
αt,z
d

The function F(αt,z ) is concave, and the lower bound of the joint likelihood can be maxt,z )imized by finding the value where ∂ F(α∂αt,z = 0. Thus, we can obtain (10).
C. DERIVATION OF (11)In this appendix, we give the derivation of the update rule for λt,z,s, (11), which canalso be derived in the same way as the update rule for αt,z .
The joint log likelihood can be written as follows, using (4) and (8).
log P(Wt, Zt, αt|αt−1 , γ , t−1 , t)= log P(Wt|Zt, t−1 , t) + C








=
λt,z,s − log  Nt,z +
λt,z,s
log 
z

+


z

s



s



log  Nt,z,w +

w





(s)
λt,z,sξt−1,z,w

− log 

s




(s)
λt,z,sξt−1,z,w

+ C,

s

(37)
where C is a constant that does not depend on λt,z,s. By using (27) and (31), we canshow that (37) satisfies the following inequality.
log P(Wt, Zt, αt|αt−1 , γ , t−1 , t)
 






old
old
≥
−λt,z,sλt,z,s −λt,z,sNt,z +
z

+


w





s


s

(s)
λold
t,z,s ξt−1,z,w

+C ≡ F λt,z,s ,




−

s

Nt,z,w +




(s)
λold
t,z,s ξt−1,z,w

(s)
λold
t,z,s ξt−1,z,w log λt,z,ss

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

(38)

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:25

where we used (27) and (31), and C is a constant that does not depend on λt,z,s. The
derivative of F(λt,z,s) is as follows.
∂ F(λt,z,s)∂λt,z,s

 

old
= −
λt,z,s −
+


w




Nt,z +

s




s


(s)
λold
t,z,s ξt−1,z,w



−


λold
t,z,s


Nt,z,w +

s




(s)
λold
t,z,sξt−1,z,w

s

(s)
λold
t,z,s ξt−1,z,w

1.
λt,z,s(39)

The function F(λt,z,s) is concave, and the lower bound of the joint likelihood can be∂ F(λt,z,s)maximized by finding the value where ∂λt,z,s= 0. Thus, we can get (11).

D. DERIVATION OF (16) t−1t−1In this appendix, we give the derivation of (16). Let N̂t−2s−1 ,z =wt =t−2s−1 N̂t ,z,w ,and N̂t,z =.
The
Dirichletpriorparameter
of
the
word
distribution
can
be
N̂
w t,z,w
rewritten as the weighted sum of the word distributions for each epoch using (12), as
follows.
S(s)
λt,z,sξt−1,z,w

s=1

=

S


t−1
λt,z,s

s=1

=

=

=

t =t−2s−1 N̂t ,z,wt−1
N̂t−2s−1 ,z

t−1
S



λt,z,sN̂t ,z,wt−1
s=1 t =t−2s−1 N̂t−2s−1 ,z
t−1


S


t−1


S


λt,z,sN̂t ,z,w
t−1
N̂


S−1
s−1
s= log2 (t−t )+1
t =t−2t−2 ,zt =t−2S−1 s= log2 (t−t )+1

=



t−1


λt,z,t φ̂t ,z,w .

λt,z,sN̂t ,z N̂t ,z,wt−1
N̂t−2s−1 ,z N̂t ,z(40)t =t−2S−1

REFERENCES
A HMED, A. AND X ING, E. P. 2008. Dynamic non-parametric mixture models and the recurrent Chinese
restaurant process: With applications to evolutionary clustering. In Proceedings of the Siam International Conference on Data Mining (SDM).

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

19:26

T. Iwata et al.

A HMED, A. AND X ING, E. P. 2010. Timeline: A dynamic hierarchical Dirichlet process model for recoveringbirth/death and evolution of topics in text stream. In Proceedings of the 26th International Conference
on Uncertainty in Artificial Intelligence (UAI).
A L S UMAIT, L., B ARBARA , D., AND D OMENICONI , C. 2008. On-line LDA: Adaptive topic models for mining
text streams with applications to topic detection and tracking. In Proceedings of the IEEE International
Conference on Data Mining (ICDM). 3–12.
A NDRIEU, C., DE F REITAS, N., D OUCET, A., AND J ORDAN, M. I. 2003. An introduction to MCMC formachine learning. Mach. Learn. 50, 1, 5–43.
B ANERJEE , A. AND B ASU, S. 2007. Topic models over text streams: A study of batch and online unsupervised
learning. In Proceedings of the Siam International Conference on Data Mining (SDM).
B LEI , D. AND F RAZIER , P. 2010. Distance dependent Chinese restaurant processes. In Proceedings of theInternational Conference on Machine Learning (ICML).
B LEI , D. M. AND L AFFERTY, J. D. 2006. Dynamic topic models. In Proceedings of the International Conference on Machine Learning (ICML). 113–120.
B LEI , D. M., N G, A. Y., AND J ORDAN, M. I. 2003. Latent Dirichlet allocation. J. Mach. Learn. Res. 3,
993–1022.
C ANINI , K. R., S HI , L., AND G RIFFITHS, T. L. 2009. Online inference of topics with latent Dirichlet allocation. In Proceedings of the International Conference on Artificial Intelligence and Statistics (AISTATS).
Vol. 5, 65–72.
G ERRISH , S. AND B LEI , D. 2010. A language-based approach to measuring scholarly impact. In Proceedingsof the International Conference on Machine Learning (ICML).
G RIFFITHS, T. L. AND S TEYVERS, M. 2004. Finding scientific topics. Proc. Nat. Acad. Sci. 101 Suppl. 1,
5228–5235.
H OFFMAN, M., B LEI , D., AND B ACH , F. 2010. Online learning for latent Dirichlet allocation. In Proceedingsof the Conference on Advances in Neural Information Processing Systems (NIPS).
H OFMANN, T. 1999. Probabilistic latent semantic analysis. In Proceedings of the International Conference
on Uncertainty in Artificial Intelligence (UAI). 289–296.
H OFMANN, T. 2003. Collaborative filtering via Gaussian probabilistic latent semantic analysis. In Proceedings of the Annual ACM SIGIR Conference on Research and Development in Information Retrieval(SIGIR). 259–266.
I WATA , T., WATANABE , S., YAMADA , T., AND U EDA , N. 2009. Topic tracking model for analyzing consumer purchase behavior. In Proceedings of the International Joint Conference on Artificial Intelligence(IJCAI). 1427–1432.
I WATA , T., YAMADA , T., AND U EDA , N. 2008. Probabilistic latent semantic visualization: Topic model for
visualizing documents. In Proceedings of the International SIGKDD Conference on Knowledge Discovery
and Data Mining (KDD). 363–371.
M INKA , T. 2000. Estimating a Dirichlet distribution. Tech. rep., MIT.
N ALLAPATI , R., C OHEN, W., D ITMORE , S., L AFFERTY, J., AND U NG, K. 2007. Multiscale topic tomography. In Proceedings of the International SIGKDD Conference on Knowledge Discovery and Data Mining(KDD). 520–529.
PAPADIMITRIOU, S., B ROCKWELL , A., AND FALOUTSOS, C. 2003. Adaptive, hands-off stream mining. In
Proceedings of the International Conference on Very Large Databases (VLDB). 560–571.
PAPADIMITRIOU, S., S UN, J., AND FALOUTSOS, C. 2005. Streaming pattern discovery in multiple timeseries. In Proceedings of the International Conference on Very Large Databases (VLDB). 697–708.
R EN, L., D UNSON, D. B., AND C ARIN, L. 2008. The dynamic hierarchical Dirichlet process. In Proceedingsof the International Conference on Machine Learning (ICML). 824–831.
S AKURAI , Y., PAPADIMITRIOU, S., AND FALOUTSOS, C. 2005. Braid: Stream mining through group lag
correlations. In Proceedings of the ACM SIGMOD International Conference on Management of Data
(SIGMOD). 599–610.
S ATO, I., K URIHARA , K., AND N AKAGAWA , H. 2010. Deterministic single-pass algorithm for LDA. In
Proceedings of the Conference on Advances in Neural Information Processing Systems (NIPS).
S TEPHENS, M. 2000. Dealing with label switching in mixture models. J. Royal Statist. Society B 62, 795–809.
T EH , Y. W., J ORDAN, M. I., B EAL , M. J., AND B LEI , D. M. 2006. Hierarchical Dirichlet processes. J. Amer.
Statist. Assoc. 101, 476, 1566–1581.
WANG, C., B LEI , D. M., AND H ECKERMAN, D. 2008. Continuous time dynamic topic models. In Proceedingsof the International Conference on Uncertainty in Artificial Intelligence (UAI). 579–586.

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.

Sequential Modeling of Topic Dynamics with Multiple Timescales

19:27

WANG, X. AND M C C ALLUM , A. 2006. Topics over time: A non-Markov continuous-time model of topicaltrends. In Proceedings of the International SIGKDD Conference on Knowledge Discovery and Data Mining (KDD). 424–433.
WATANABE , S., I WATA , T., H ORI , T., S AKO, A., AND A RIKI , Y. 2011. Topic tracking language model for
speech recognition. Comput. Speech Lang. 25, 2, 440–461.
W EI , X., S UN, J., AND WANG, X. 2007. Dynamic mixture models for multiple time-series. In Proceedings ofthe International Joint Conference on Artificial Intelligence (IJCAI). 2909–2914.
Z HANG, J., S ONG, Y., Z HANG, C., AND L IU, S. 2010. Evolutionary hierarchical Dirichlet processes for multiple correlated time-varying corpora. In Proceedings of the 16th ACM SIGKDD International Conferenceon Knowledge Discovery and Data Mining (KDD). 1079–1088.
Received January 2011; revised June 2011; accepted September 2011

ACM Transactions on Knowledge Discovery from Data, Vol. 5, No. 4, Article 19, Publication date: February 2012.