DEIM Forum 2016 G5-2

大規模疫病データのための非線形テンソル解析松原靖子†

櫻井 保志†

Willem G. van Panhuis††

Christos Faloutsos†††† 熊本大学自然科学研究科情報電気電子工学専攻†† Dept. of Epidemiology, University of Pittsburgh††† Dept. of Computer Science, Carnegie Mellon University
E-mail: †{yasuko,yasushi}@cs.kumamoto-u.ac.jp, ††wav10@pitt.edu, †††christos@cs.cmu.edu
あらまし 本論文では，大規模疫病データのための非線形モデル解析手法である F UNNEL について述べる．F UNNEL
は，(disease, location, time) の三つ組で構成された d の疫病と l の州を含む，長さ n のテンソルデータが与えられたときに，そのテンソルの中から非線形パターンや外れ値を統合的に解析，要約，表現する．提案手法は (a) 疫病の季節
性，ワクチンの出現，外部要因による疫病の発生等の重要なパターンを発見し，(b) パラメータのチューニングを必要とせず，(c) 計算量は入力データのサイズに対して線形である．さらに，本手法は (d) 麻疹やインフルエンザ等の疫病だけでなく，コンピュータウィルスの感染パターンについても同様に解析することができる．実データを用いた実験
では，F UNNEL が大規模疫病テンソルの中から (P1) 病気の季節性，(P2) ワクチン効果，(P3) 地域性，(P4) 外因性による疫病パターン，(P5) データ入力エラー等の重要パターンを正確に発見することを発見した．
キーワード

疫病，テンソル，時系列データ，特徴自動抽出

1. ま え が き

感受性保持者の潜在人口を示している．図 1 (c) は各疫病の季節
性の強さ（半径）とそのピークの時期（角度）を示したもので大規模な疫病データの中から有用なパターンを発見し，感染ある．例えば，1 月から 2 月にかけてインフルエンザの流行の症の性質や法則を見つけ出すことは，医療分野や公衆衛生の研ピークがあり，麻疹のような小児呼吸器疾患は春，ライム病の究，および政策や社会活動等において非常に重要な問題である．

ようなダニ媒介疾患は夏がピークとなっている．また，淋病の本研究では，大規模疫病データを対象とし，重要なパターンをような性感染症（STD）には周期性がないことが分かる．提案

自動抽出するための非線形モデル解析手法として F UNNEL を提手法は，事前知識やパラメータ設定を要せず，すなわち，ユー

案する [4]．（注 1） より具体的には，以下の問題を扱う．

ザの介入なしに，様々な疫病に関する重要な特徴を自動的に抽d 種の疫病，l 箇所の地域，n の期間で構成される大規模疫病
データ集合 X が与えられたとき，以下の重要なパターンを自動抽出する．

•

季節性や感染率等の各疫病の基本的な特徴•

突発的な流行の外部要因•

入力エラー等の外れ値

本研究では，疫病に関するデータとして Tycho [10] を扱う．

Tycho は，アメリカ合衆国における 125 年間に渡る 56 種の疫病の患者の報告件数で構成される．（注 2）
具体例. 図 1 は，F UNNEL による解析結果を示している．具体的には，図 1 (a) において灰色の丸印は 1928 年から 1982 年までの，麻疹の感染者数を示しており，赤の実線は提案手法の学
習結果を示している．シーケンスは年の周期性があり，そして数年毎に大規模な感染（1941 年，1958 年）と小さな波（1940
年と 1947 年）が現れており [6]，これはスキップ現象として知られている [9]．また，感染者数は 1965 年において急激に減少しており，これはワクチン接種が 1963 年に開始されたことによるものである．考案手法とモデルはこのようなパターンを高精度で表現することができる．図 1 (b) は各州の麻疹に関する出することができる．
本論文の貢献. F UNNEL は以下の特長がある．
（ 1 ） 提案モデルは，麻疹，インフルエンザ，天然痘をはじめとする様々な感染症の振る舞いを表現し，図 1 に挙げられる
ような季節性，ワクチン効果，外部要因による突発的流行等の重要な時系列パターンを自動抽出する．
（ 2 ） 計算コストは入力データのサイズに対して線形である．
（ 3 ） 提案モデルは，SIRS モデル等の既存の疫病モデルを一般化すると同時に，コンピュータウィルスの感染パターン等の様々な種類の非線形伝染過程を柔軟に表現することができる．

2. 提案モデル
本章では，提案モデルである F UNNEL について述べる．

2. 1 FUNNEL の概要
本研究で扱う Tycho データは，1888 年から現在にかけてのアメリカ合衆国 50 州における 56 種の疫病を患者の報告件数（計 87,950,807 件）で構成される．本データは (disease, location,

timestamp) の三つ組で表現され，それぞれ，d 種の疫病，l 箇所の地域/州，n の期間（一週間単位）から構成される．この疫病
データは，3 階のテンソル X ∈ Nd×l×n として表現することが（注 1）：http://www.cs.kumamoto-u.ac.jp/∼yasuko/software.html

でき，X の要素 xij (t) は時刻 t において i 番目の疫病 (disease)（注 2）：Project Tycho at University of Pittsburgh: http://www.tycho.pitt.edu/

が j 番目の地域/州 (location) が出現した頻度を示す．例えば，

［知見 2］（疫病減少効果） 疫病の多くは，ワクチン，抗生物質，

4

Count

6

x 10

Vaccination

公衆衛生管理により，大幅に減少あるいは撲滅している．

4
Original
I(t)

2
0

1930

1940

1950

1960

1970

Count (log)

Year
5

10

1930

1940

1950

1960

1970

次に，各地域における傾向について述べる．一般に，各州の感染者数は高い相関を持つが，同時に，患者数の割合は異なる．

1980

これは，各州での感受性保持者の潜在人口が異なることに起因

Original
S(t)
I(t)
V(t)

する．例えば，麻疹の主な感染対象は小児であるため，小児が

1980

Year

(a) F UNNEL F IT の学習結果 (麻疹)多い地域 (州) には，より多くの麻疹患者が見込まれる．図 1 (b)では，ニューヨーク (NY), ペンシルベニア (PA), カリフォルニ
ア (CA), テキサス (TX) に，多くの感染対象者 (つまり小児) がいることがわかる．
［知見 3］（地域性） 各州における感染者数の推移は正の相関を持つが，一方で，感受性保持者の潜在人口が異なる．
最後に，外因性による傾向について述べる． 図 1 (a) の 1941
年と 1958 年等において，麻疹の感染者数が大幅に増加して(b) 各州における感受性保持者の潜在人口 (麻疹)
0.5

在する．

March (3)Measles 0.4
Rubella0.3

実データの中には，予期しないエラー値や入力ミスによる突

Mumps

June (6)

February (2)Smallpox

0.2

Chickenpox
Influenza

July (7)January (1)Streptococcal sore throat
Lymedisease

発的な外れ値が存在する．
［知見 5］（入力エラー） 疫病の感染パターンに合致しない外れ

0.1
Rocky mountain spotted fever［知見 4］（外部要因による突発的流行） 疫病の多くには，数年に一度現れる，外部要因に由来する大規模な流行パターンが存

April (4)May (5)いる．

Gonorrhea

December (12)Typhoidfever
Typhus fever
Cryptosporidiosis
November (11)September (9)August (8)値や入力エラーが存在する．
本研究では，大規模疫病データを対象とし，上記の重要な要素を表現する統合モデルとして，F UNNEL を提案する．具体的には，次の 5 つの要素を全て表現する．

•

October (10)(c) 各疫病の季節性: 半径は強さ，角度はピークを示す．
図1

F UNNEL F IT を用いたパターン発見と出力結果．
Fig. 1Modeling power of F UNNEL F IT.

(‘measles’, ‘PA’, ‘April 1-7, 1931’; 4740) の場合，麻疹 (measles)(P1): 季節性•

(P2): 疫病減少効果•(P3): 地域性•

(P4): 外部要因による突発的流行•(P5): 外れ値，入力エラー

がペンシルベニア州 (PA) で 1931 年 4 月 1 日から 4 月 7 日の次章では，提案モデルの詳細を示す．まず，(a) 単一の疫病

間に 4740 件報告されたことを表す．ここで，i 番目の疫病の j

シーケンスの場合（例えば，ニューヨーク州 (NY) の麻疹の感番目の州におけるシーケンス xij = {xij (t)}nt=1 を，ローカル

染者数など）についてのモデルを述べ，次に発展として，(b) 複（州）の疫病シーケンスと呼ぶ．同様に，各州のシーケンスの数シーケンスの場合（つまり，d 種の疫病，l 箇所の地域に対す

合計 x̄i = {x̄i (t)}nt=1 を，グローバル（国）の疫病シーケンス

る個々のパターン）のモデルについて述べる．

と呼び，x̄i (t) は i 番目の疫病の時刻 t における出現頻度の国内
Pの総数（つまり，x̄i (t) = lj=1 xij (t)）とする．

2. 2 FUNNEL - 単一の疫病シーケンスの場合

疫病データの考察. ここでは，疫病データのモデル化にあたり

まず最も簡単な場合として，単一の疫病シーケンスのモデル化について述べる．

重要となる要素について事前考察を行う．図 1 (c) で示したよう

2. 2. 1 基本モデル - F UNNEL - BASEに，疫病には季節性を伴う場合が多く，例えば麻疹やおたふく提案モデルは以下の 3 つの状態から構成される．

かぜに代表される小児呼吸器疾患は春に，ライム病のようにダ
ニ媒介疾患は夏がピークとなっている． 同様に，インフルエン
ザは気温と強い負の相関があり，気温の低い季節に活発になる．
［知見 1］（疫病の季節性） 疫病の多くは年単位の周期性を伴い，気温や季節と強い相関を持つ．
ワクチンの導入，抗生物質の使用や公衆衛生の管理により，•Susceptible （感受性保持者）: 現時点では疫病に感染していないが，近隣の感染者から伝染する可能性のある状態•

Infected （感染者）: 既に疫病に感染しており，他者へ伝染させる可能性のある状態•

Vigilant （免疫保持者）: 疫病から回復した状態，あるいはワクチン接種等により疫病に感染する可能性がない状態

多くの疫病が，数十年の間に大幅に減少している．例えば，図 1

図 2 (a) は，基本モデルである F UNNEL - BASE の状態遷移の(a) で示したように，麻疹の感染者数はワクチン接種が 1963 年

様子を示す．ここで，β(t) は感染者と感受性保持者の間の疫病に開始されて以来，大幅な減少傾向にある．

ε (t)
!!

β (t)

γ

部要因によって突発的に流行することがある．例えば，豚イン

β (t)

!!

γ

フルエンザの流行が発生した場合には，過去の状況よりも多くの感染者が見込まれることになる．

θ (t)

#!

"!

δ

率 ǫ(t) を導入する．具体的には，図 2 (b) に示すように，時刻 t(a) F UNNEL - BASE(b) F UNNEL -RE
図 2 F UNNEL の状態遷移の様子．

において外部ショックによるイベントが起きた場合，感受性保

#!
"!

δ

Fig. 2 F UNNEL diagrams.

の感染率を，δ は感染者の回復率を示す．γ は，免疫力の減少率を示す（注 3）．さらに，先述の要素 (P1) を表現するために，本この現象を表現するために，本研究では一時的な感受性変化

持者 S(t) の感染者数が一時的に変化する．これにより，大規模な感染拡大パターンを表現する．
［モデル 3］（F UNNEL -RE ） 提案モデルは次の式で構成される：S(t + 1) =S(t) − β(t)ǫ(t)S(t)I(t) + γV (t) − θ(t)S(t)て表現する．

I(t + 1) =I(t) + β(t)ǫ(t)S(t)I(t) − δI(t)［モデル 1］（F UNNEL - BASE ） S(t), I(t), V (t) を時刻 t におけV (t + 1) =V (t) + δI(t) − γV (t) + θ(t)S(t)研究では，感染率 β(t) は時刻 t における周期性を伴う関数としる感受性保持者，感染者，免疫保持者の総数とする．F UNNEL -ここで，ǫ(t) は，一時的な感受性変化率を表す：

BASE は，次の式で構成される．

S(t + 1) =S(t) − β(t)S(t)I(t) + γV (t)I(t + 1) =I(t) + β(t)S(t)I(t) − δI(t)V (t + 1) =V (t) + δI(t) − γV (t)
(1)

2π
ここで，β(t) = β0 · 1 + Pa · cos Pp (t + Ps ) , Pp = 52 とし，ǫ(t) = 1+

k
X

(T )

f (t; ei), f (t; e(T ) ) =

i=1

(

ǫ0

(tµ − tσ < t < tµ + tσ )0(else)

k は外部ショックの回数を示す（k = 0 の場合には ǫ(t) = 1）．
各外部ショックは次のパラメータで構成される：e(T ) ={tµ , tσ , ǫ0 }.
•tµ : 外部ショックイベントの中心時刻．

1, V (1) = 0 とする．

•

tσ : イベントの長さ．

まとめると，F UNNEL - BASE は次のパラメータ集合で構成され•

ǫ0 : 外部ショック効果の強さ．

（注 4）

さらに，N = S(t) + I(t) + V (t), S(1) = N − 1, I(1) =る：b = {N, β0 , δ, γ, Pa , Ps }.

•N : 疫病の潜在的人口．

•

β0 : 感染率の年間平均値．

•

δ: 疫病の回復率．

•

γ: 免疫力の減少率．

•

Pa : 季節性の強さ，振幅．

•

Ps : 季節性の位相．

2. 2. 2 疫病減少効果- F UNNEL -R
2 つ目の要素 (P2) のために，モデルを次のように拡張する．
［モデル 2］（F UNNEL -R ） θ(t) は疫病減少効果の強さを示す．

(3)

2. 3 FUNNEL- 複数シーケンスの場合
前節では，単一の疫病シーケンスについての F UNNEL の振る舞いについて述べた．ここからは，疫病テンソル X が与えられた上での d 種の疫病，l 箇所の地域の個々の振る舞いを表現するモデルについて考える．まず，最も簡易的な解決法として挙げられるのが，テンソル X を (d × l) 個の長さ n のシーケン
ス集合：{xij }d,li,j=1 （つまり，ローカル（州）の疫病シーケンス）とみなし，各シーケンスに対し個別にモデルパラメータ集
合：{b, r, e(T ) } を推定する場合である．しかしながら，これらのローカルの疫病シーケンスのペア (disease, state) の中には，
例えば，アラスカ州におけるライム病のように，感染者数が非S(t + 1) =S(t) − β(t)S(t)I(t) + γV (t) − θ(t)S(t)

常に少ないスパースな時系列データを含む場合が多く，直接的

I(t + 1) =I(t) + β(t)S(t)I(t) − δI(t)に個々のモデルパラメータを推定することが困難である．さらV (t + 1) = V (t) + δI(t) − γV (t) + θ(t)S(t)(2)
ここで，疫病減少効果の開始時刻をtとし，次のように定義すθ(に，上記の方法では，グローバル（国）の疫病の振る舞いをモ

このモデルは，F UNNEL - BASE に疫病減少効果の要素 θ(t) を加FUNNEL - パラメータ集合. 本研究の目的は，大規模疫病デーえたものである．これは，図 2 (b) に示すように，感受性保持タ X ∈ Nd×l×n の中から重要なトレンドや特徴を抽出するこ者が，ワクチン接種や抗生物質投与等により，直接免疫保持者とである．図 3 は提案モデルの概要を示す．疫病テンソル Xの状態に移行する．基本パラメータ b に加え，F UNNEL -R は次が与えられたとき，提案手法は 5 つの重要な要素を発見する：のパラメータを要する：r = {tθ , θ0 }.

(P1) B: 季節性等の疫病の基本的な振る舞い，(P2) R: 疫病減•tθ : 疫病減少効果の開始時刻．

少効果，(P3) N: 地域性，(P4) E : 外部ショックイベント集合，•

θ0 : 疫病減少効果の強さ，拡散率．

る：θ(t) =0

θ0(t < tθ )(t >= tθ )

デル化することができない．この問題を解決するために，本研
究では，l 箇所の全ての州で，一部のパラメータを共有する手法を用いる．

(P5) M : 入力エラー．ここで，(P1) と (P2) は，グローバル

2. 2. 3 外部ショック- F UNNEL -RE（国）の振る舞いを表現するパラメータ集合，(P3) はローカル

次に，3 つ目の要素：(P4) について述べる．疫病の多くは外（州）の個々のパターン，そして，(P4) と (P5) は X の外部要

（注 3）
：γ は，出生率および死亡率も表現する．

因に関するパラメータ集合である．

（注 4）：1 年間は 52 週で構成される．

［定義 1］（FUNNEL のパラメータ集合） F を疫病テンソル X

(P1, P2): global/country view! (P3): local/state view! (P4, P5): extra view - E: shocks & M: mistakes!

diseases!

a
st

t

!
es

l

l
=
X

d

N
B

d
time!

n

l

l

θ 0 , tθN, β 0 , δ, γ , Pa , Ps

,

d

R

E

d

d

,

,

,
n

2

6

M

d

n

(a) F UNNEL のモデル構造 (F = {B, R, N, E , M })

(ek)!

l
d

l

e1(S )

=

E

(e2)!

n

ek(S )e2(T )+ ... +ek(T )
=

n

d
(e1)!

e2(S )e1(T )+e1( D)

(e1)!

e2(D)

(e2)!

(ek)!

ek( D)kl

E(S ) t , t , εµe(D) E(D) 1

σ

0

E(T )

k

3

k

(b) 外部ショックテンソル (E = {E(D) , E(T) , E(S) })図3

F UNNEL の概要．提案手法は疫病テンソル X の中から重要なパターンを自動抽出する．

Fig. 3 llustration of F UNNEL structure: we extract the important behavior of epidemics from X .

次に，各パラメータの詳細を述べる．

(P1), (P2) グローバルパラメータ集合. 次のパラメータ集合は，

15000

5000
0

Count (log)［定義 2］（基本行列 B (d × 6)） B を d 個の疫病に関する基本パラメータ集合 B = {b1 , . . . , bd } とする．ここで， bi は i 番
目の疫病のパラメータを表す．
フロリダ州で同じである．同様に，麻疹のワクチンが導入され［定義 3］（減少行列 R (d × 2)） R を疫病減少効果に関するパラメータ集合 R = {r1 , . . . , rd } とする．ここで，ri は，i 番目の疫病のパラメータを表す．

(P3) ローカルパラメータ集合. 次の課題は，テンソル X の中の地域固有のパターンの表現方法である．ニューヨーク州とフロリダ州それぞれの麻疹の感染者数の推移に変化をもたらす原因
は一体なんであろうか．各州における麻疹のウィルスに地域性や大きな個体差はない．唯一の違いは，感受性保持者の潜在人口の差である．そこで本研究では，l 箇所の州でグローバルパ
ラメータ集合を共有し，疫病の潜在的人口 Nij のみを i 番目の疫病，j 番目の州毎に導入する．具体的には，モデル 3 において，次の式を満たす：Nij = Sij (t) + Iij (t) + Vij (t). 例えば，
小児呼吸器疾患である麻疹について考えると，ニューヨーク州にはフロリダ州よりも多くの小児がいるため，（つまり，潜在的な麻疹感染対象者が多いため）より多くの感染者が見込まれる．
［定義 4］（地域行列 N (d × l)） N を d 種の疫病，l 箇所の地域における潜在的人口を表すパラメータ：N = {Nij }d,li,j=1 とする．ここで，Nij は，i 番目の疫病の j 番目の州における潜在的人口を表す．

(P4) 外部ショック. 本研究では，突発的に発生する疫病の流行
イベントを 3 つの要素：(disease, state, time) として表現する．
例えば，1946 年に国家規模で起きたインフルエンザの歴史的大流行の場合は，(e1) “influenza, country-wide, 1946” として表現される．同様にして，2007 年にユタ州で起きた地域規模のOriginal
2007 2008 2009 2010 2011
2012
S(t)
Year
I(t)
V(t)

0

10

0

5

10

2007 2008 2009 2010 2011 2012
Year
Original
S(t)
I(t)
V(t)

0

2007 2008 2009 2010 2011 2012(a) 外部ショックとしての学習結果

102007 2008 2009 2010 2011 2012(b) 入力エラーしての学習結果

図 4 ジアルジア症に対する学習結果の比較．

る時期は，どの州においても同じであり，導入されたワクチンの効果が出始める時期も一致する．

OriginalI(t)5000

グローバルな振る舞いを表現し，l 箇所の州で共有する．

例えば，麻疹ウィルスの感染率及び回復率はニューヨーク州と10000
Original
I(t)

10000

Count (log)を表現する全パラメータ集合 F = {B, R, N, E , M }) とする．

Fig. 4External shock vs. mistake for giardiasis.

クリプトスポリジウム症の流行は，(e2) “cryptosporidiosis, Utah,2007” として表現される．これらの複数のイベントを表現する
ために，新たにパラメータ集合として外部ショックテンソル Eを導入する．図 3 (b) に示すように，E は，k 個のイベント集
合で構成される．外部ショックテンソル E は，3 つの要素行
列 {E(D) , E(S) , E(T) } として圧縮することもできる．さらに，
単一の外部ショックイベントは，三つ組のベクトル集合 {e(D) ,e(S) , e(T ) } として次のように表現される．
•

疫病ベクトル e(D) : 疫病の ID•

地域ベクトル e(S) : 各地域のショックの強さ•

時間ベクトル e(T ) : ショックの発生期間

ここで，2. 2. 3 節で示したとおり，e(T ) はグローバルなパラメー(S)

タであり，e(S) はローカルなパラメータ：e(S) = {ej }lj=1 である．これは，モデル 3 におけるパラメータ ǫ0 の l 箇所の州毎の変化率を表す．より具体的には，j 番目の州での外部ショッ(S)

クの強さを ǫ0 · ejとする．

［定義 5］（外部ショックテンソル E ） E を k 個の外部ショック
イベントで構成される 3 階のテンソル：E = {E(D) , E(S) , E(T) }とする．ここで，3 つ組の行列は，疫病，州，時間の各要素を表現するパラメータである．

(P5) 入力エラー. 実データには，多くのエラー値や報告の誤りが含まれる．提案モデルは，これらの入力エラーを外れ値として取り除く必要がある．
［定義 6］（入力エラーテンソル M ） M を入力エラーに関する

3 階のテンソルとする．M の各要素 mij (t) は，i 番目の疫病，

入力エラーテンソル M は以下で構成される．

j 番目の州，時刻 t の値を表す．

•

テンソル M 内の非ゼロ要素の個数：log∗ (|M |).

ここで，M は非常にスパースであり，多くの場合において，•

非ゼロ要素の位置とその値 mij (t)：log(d), log(l), log(n),log∗ (mij (t)).

mij (t) = 0 である．
図 4 は，ジアルジア症に対する 2 種類の学習結果の様子を示したがって，CostM (M ) = log∗ (|M |) +∗P |M |

mij (t)>0 (log(d) +している．灰色の丸印はオリジナルデータ，赤線はフィッティン

log(l) + log(n) + log (mij (t))). ここで |M | は M 内の非ゼログ結果を示す．図に示すように，ジアルジア症は 2006 年におい要素の個数を示す．

て大きな外れ値（およそ 10, 000）がある．図 (a) は外部ショッ

疫病データの符号化コスト. ハフマン符号を用いた情報圧縮

クとして学習した結果であり，図 (b) は入力エラーとして学習

では，モデルパラメータ集合 F が与えられた際の X の符号した結果である．(a) において，近隣の点が外れ値に大きく影化コストを負の対数尤度を用いて次のように表現することがP−1できる [2]：CostC (X |F ) = d,l,ni,j,t=1 log2 pGauss(µ,σ 2 ) (xij (t) −響されている一方で，(b) では外れ値が独立した点として認識

mij (t) − Iij (t))．ここで，符号化にはガウス分布を用い，xij (t),され，より正しく学習できていることがわかる．

mij (t) は，それぞれ，疫病テンソル X ，入力エラーテンソルの3. 最適化アルゴリズム

要素の値を示す．Iij (t) は，モデル 3 で示したように，感染者

本章では，モデルの学習アルゴリズムである F UNNEL F IT について述べる．提案アルゴリズムの目的は，与えられた疫病テ数の推定値を示す．さらに，µ, σ 2 はそれぞれ，オリジナルデー
タの値と，推定値の間の平均，分散を示す．（注 7）

ンソル X に対し，重要なパターンを自動抽出することである．

符号化コスト関数. まとめると，モデルパラメータ集合 F が与［ 問題 1］ (disease, state, time) の三つ組で構成される疫病テンえられたときの X の符号長は次のように表現される．

ソル X が与えられたとき，X を表現する最適なモデルパラメー
タ集合 F = {B, R, N, E , M } を発見する．

CostT (X ; F ) = log∗ (d) + log∗ (l) + log∗ (n)

ここでの重要な課題は，(a) テンソル X の中の重要なパターン+CostM (B) + CostM (R) + CostM (N)を表現するための最適なパラメータをどのように推定するか，+CostM (E ) + CostM (M ) + CostC (X |F )(4)(b) 外部ショックの最適な個数 k をどのように決定するか，そして (c) データ X の中に含まれる 入力エラーをどのように取りしたがって本論文の次の目標は，上記のコスト関数を最小化す

除くかである．

ることである．

3. 1 モデル化とデータ圧縮

3. 2 多階層最適化アルゴリズム

本研究では，大規模疫病テンソルデータ X を表現するためここまでは，パラメータ集合の候補解 F が与えられた上でテに最小記述長 (MDL:minimum description length) に基づく新た

ンソル X を表現するための方法について述べた．次の課題は，な符号化スキームを導入する．直感的に，データがより圧縮でどのようにして全てのパラメータ集合 F = {B, R, N, E , M }きれば，よりよいモデルであるとみなす．

を最適化するかである.

モデル表現コスト. F UNNEL のパラメータ表現コストは以下の要素から構成される．

•∗

疫病の総数 d，州の総数 l，および，時系列の長さ n：∗

図 2. 3 で示したように，F UNNEL は，複数のパラメータ集合
から構成され，それらは，疫病テンソル X の中のグローバル∗

log (d) + log (l) + log (n) ビット.

（国）かローカル（州）のいずれかのパターンを表現する．例
えば，基本行列 B と減少行列 R は，各疫病のグローバル（国）（注 5）

基本行列 (B), 減少行列 (R), 地域行列 (N) にそれぞレベルの振る舞いを表現し，一方で，地域行列 N は，ローカ

れ d × 6, d × 2, d × l のパラメータを要する．まとめると，

ル（州）なトレンドを表現する．また，外部ショックテンソル

CostM (B) + CostM (R) + CostM (N) = cF · d(6 + 2 + l). こ

E = {E(D) , E(S) , E(T) }) は，グローバル (E(D) , E(T) )，ロー

こで，cF は浮動小数点のコストを示す．（注 6）

カル (E(S) ) のパラメータとして分解できる．同様にして，入同様に，外部ショックテンソル E = {E(D) , E(S) , E(T) }) は次

力エラーテンソル M は，{M(D) , M(S) , M(T) } として分解での要素で構成される．

きる．そこで本研究では，上記の性質を用いて，より効率的か

•

•
•
•

つ効果的にパラメータの推定を行う手法として，多階層最適化∗

外部ショックの個数 k ：log (k) ビット．
外部ショック (疫病) 行列 E

(D)

外部ショック (時間) 行列 E

アルゴリズムを提案する．より具体的には，全パラメータ集合：k log(d).

(T)の中の各要素 e(T )=FL に分解しそれらのパラメータ集合を個別に学習し，高速に{tµ , tσ , ǫ0 }：log(n), log(n), cF .
•外部ショック (地域) 行列 E

F を，グローバル，ローカルの二つの部分パラメータ集合 FG ,(S)

：cF · kl.

まとめると，外部ショックテンソル E のモデルコストは次のようになる．CostM (E ) = log∗ (k)+k log(d)+2 log(n)+cF ·(1+l) .

モデル推定を行う．
提案アルゴリズムは次の 2 つのステップで構成される．

•

G LOBAL F IT: グ ロ ー バ ル（ 国 ）の 疫 病 シ ー ケ ン ス（注 5）：ここで，log∗ は整数のユニバーサル符号長を表す．

（注 7）：ここで， µ, σ 2 は 2cF ビットを要するが，これらのコストは定数である（注 6）：本論文では 4 × 8 ビットとする．

ため，モデル推定の際には除外することができる．

{x̄i }di=1 に 対 し ，最 適 な モ デ ル パ ラ メ ー タ FG = {B, R,

ことがわかる．表 6 は，ワクチンがライセンス化された年と，} を推定する．

F UNNEL において自動発見された疫病減少効果の開始年を比較

L OCAL F IT: ローカル（州）の疫病シーケンス {xij }d,li,j=1に対し，最適なモデルパラメータ FL = {N, E(S) , M(S) } を推したものである．基本的に，ライセンス化されて 2，3 年後に減少効果が見られることがわかる．これは，ワクチンプログラ

定する．

ムの普及にラグが生じる影響と考えられる．一方で，インフル

ここで，i 番目の疫病に対するグローバル（国）のシーケンス

エンザに対する疫病減少効果は発見されなかった．これは，イx̄i は l 箇所すべての州のローカルシーケンスの合計値とする：
Px̄i (t) = lj=1 xij (t). アルゴリズム 1 は F UNNEL F IT の概要を

ンフルエンザウィルスの変異が活発であり，ワクチンの効果が示す．疫病テンソル X が与えられたとき，提案アルゴリズムは(P3) 地域性. F UNNEL は，国家規模の疫病だけでなく，地域固モデルパラメータ集合 F を自動抽出する．

有のパターンを発見することができる．例えば，図 1 (b) で示

Algorithm 1 F UNNEL F IT (X )多くの麻疹の感受性保持者が存在する．図 7 は，実際の感染者

E

(D)

,E

(T)

,M

(D)

,M

(T)

•

現れにくいことが原因である．

した通り，ニューヨーク州 (NY) とペンシルベニア州 (PA) には
1: Input: Tensor X (d × l × n)
2: Output: Complete set of parameters, i.e., F = {B, R, N, E , M }数とその学習結果である．図に示す通り，F UNNEL は周期的な

3: /* Parameter fitting for global-level sequences */感染パターンや疫病減少効果を表現すると同時に，ローカルな

4: {FG } =G LOBAL F IT (X );外部ショックのパターンも柔軟に表現している．

5: /* Parameter fitting for local-level sequences */(P4) 外部ショックイベント. 図 5 に示すように，F UNNEL F IT は

6: {FL } =L OCAL F IT (X , FG );次に挙げられるような重要な外部ショックイベントを自動発見

7: return F = {FG , FL };した．

•

4. 評 価 実 験
本論文では F UNNEL の有効性を検証するため，実データを用いた実験を行なった．具体的には，本章では以下の項目について検証する．

Q1 疫病のパターン抽出に関する提案手法の有効性
Q2 提案アルゴリズムの精度の検証
Q3 パターン抽出に対する計算時間の検証

小児性疾患である (a) ジフテリア (Diphtheria) と (b) 猩紅熱 (scarlet fever) は，第二次世界大戦前後に複数年に渡る流行パ
ターンを持つ．

•(g) インフルエンザ (Influenza) は，1929 年，1941 年等に大規模な流行が確認できる [5]．

•(h) ワクチンの普及が不十分だったことにより，天然痘(smallpox) の流行が 1937-39 年にかけて発生した [3]．
•(f) 水媒性の疫病であるクリプトスポリジウム症 (Cryp-

tosporidiosis) が，2007 年におけるユタ州の公共プール汚染を原

4. 1 疫病データからの重要パターンの発見

因として局所的に流行した [1].

本節では，大規模疫病データに対する F UNNEL の情報抽出の(P5) 入力エラー. 提案手法である F UNNEL F IT の強みの一つと効果を検証する．
図 5 は，9 種の主要な疫病に対するモデル学習の結果を示し

して，ノイズに対する頑健性が挙げられる．図 5 (i) に示す通り，F UNNEL F IT は，腸チフス (typhoid fever) のデータに含まれている． 灰色の丸印はオリジナルデータを，赤の実線は提案手る，1940 年代における欠損値を補完すると同時に，1953 年に法の学習結果 I(t) を示している．各シーケンスはそれぞれ，線おける外れ値を取り除くことができる．

形スケール（上段），対数スケール（下段）で示しており，下段については，感染者の推定値の他に，感受性保持者 S(t) と免疫保持者 V (t) の推定値も表している．
以下では，先述の 5 つの重要な疫病の要素に関して，提案手

4. 2 提案手法の精度と学習時間
まず，提案モデルの精度を検証するため，既存手法である

SIRS モデルと SKIPS [9] との比較を行った．本研究ではさらに，疫病減少効果の効果を検証するために，外部ショックと入力エ法で得られた知見について議論する．

ラーに関するパラメータを取り除いた F UNNEL -R とも比較を(P1) 季節性. 図 1(c) において既に示したように，F UNNEL は，行った．図 8 (a) は，オリジナルデータと推定値との二乗平均

季節性における次の 4 つのカテゴリを発見した．

誤差 (RMSE: root mean square error) を示している．ここで，(a)•インフルエンザは 1 月から 2 月にかけてピークを持つ．

は，グローバル（国）の疫病シーケンス集合 {x̄i (t)}d,ni,t ，そし•

麻疹，おたふくかぜ，水疱瘡などの小児呼吸器疾患は強て (b) はローカル（州）のシーケンス集合 {xij (t)}d,l,ni,j,t についい周期性を持ち，春にかけて感染者が増加する [7].

てそれぞれ比較している．SIRS モデルは，季節性を伴うパターライム病のようなダニ媒介疾患やクリプトスポリジウムンを表現できない．SKIPS は周期性を発見できるが，疫病減少症のような水媒性の疫病は，夏に強いピークを持ち，保菌生物効果や外部ショックイベントを表現できない．図に示す通り，

やヒトの行動，あるいは気候等の要素と深い関連性を持つ [8]．

SIRS モデルと SKIPS が疫病の複雑なパターンの表現に失敗し淋病のような性感染症（STD）には周期性が見られない．

たのに比べ，提案手法は高い精度でのデータの学習に成功した．

(P2) 疫病減少効果. F UNNEL は，大規模疫病データの中から，次に，F UNNEL F IT の計算時間を検証する．図 9 はデータのサ

自動的に疫病減少効果を発見することができる．例えば，図 5

イズを変化させた上での提案手法の計算時間を示している．こ(a-e) において，多くの小児性疾患は疫病減少効果が見られる

こでは，それぞれ，(a) 疫病の総数 d，(b) 州の総数 l，(c) シー••

World War II

Er= 1.35e+02 (2.0e+04, 0.54, 0.45, 0.13) , P=(0.1, 12), (1509, 0.03), k=16

10000

1000

4000
2000

1930

1940

1950
Year

1960

5

0

1970
Original
S(t)
I(t)
V(t)

Count

10

1930

1935

1940

1945 1950
Year

1955

1960

0

1965

1940

1950
Year

1960

1940

1945 1950
Year

1955

1960

10

1965

Er= 2.95e+01 (3.4e+03, 0.51, 0.46, 0.12) , P=(0.2, 27), (4889, 0.00), k=1

5000

OriginalI(t)Outbreakin Utah600
Count

Original
I(t)

400
200

1975

1980

1985
Year

1990

1995

0

2000

1975

1980

1985

1990 1995
Year

2000

2005

5
5

10
Count

Count

Original
S(t)
I(t)
V(t)

0

2010
Original
S(t)
I(t)
V(t)

(Missing values)

0

1970

1975

1980

1985
Year

1990

1995

2000

1975

(d) 風疹 (Rubella) (P1), (P2), (P4)
5

x 10

1985

1990 1995
Year

2000

10

2010

Count

1500
1000

1945

0

1950

10

1935

1940
Year

1945

1930

1935

1940

(g) インフルエンザ (Influenza) (P1), (P4)
図5

Fig. 5

2009
Year

2010

2011

2012

OriginalI(t)Mistake

1000
500
0

1945

1930

1940

1950

1960
Year

1970

1980

5

10

5

10

10

1950

2012

Er= 4.26e+01 (1.2e+04, 0.59, 0.48, 0.11) , P=(0.2, 24), (2081, 0.02), k=16

OriginalI(t)Original
S(t)
I(t)
V(t)

Original
S(t)
I(t)
V(t)

0

0

1930

2008

Year

Original
S(t)
I(t)
V(t)

5

2007

Count

1940
Year

Count

1935

2011

1500

Smallpox epidemics
in 1937−39

500
1930

2010

(f) クリプトスポリジウム症 (Cryptosporidiosis) (P1), (P3), (P4)

2000

1

2009
Year

Original
S(t)
I(t)
V(t)

Er= 1.07e+02 (1.9e+04, 0.58, 0.44, 0.09) , P=(0.1, 47), (1041, 0.03), k=8

Original
I(t)

2008

0

2005

(e) 水疱瘡 (Chickenpox) (P1), (P2)

Er= 7.04e+03 (8.8e+05, 0.53, 0.49, 0.05) , P=(0.5, 9), (1821, 0.00), k=8

2

1980

Count

10

2007

Count

1970

10

Count

2000

(c) おたふくかぜ (Mumps) (P1), (P2), (P4)

10000

0

Count

1990

800

1000

0

1980
Year

Original
I(t)

2000

3

1970

Er= 9.27e+02 (8.8e+04, 0.64, 0.52, 0.12) , P=(0.2, 45), (4889, 0.03), k=0

Count

Count

1935

(b) 猩紅熱 (Scarlet fever) (P1), (P2), (P4)

Er= 1.30e+02 (3.3e+04, 0.55, 0.46, 0.09) , P=(0.2, 45), (3173, 0.02), k=83000

Original
S(t)
I(t)
V(t)

0

1930

(a) ジフテリア (Diphtheria) (P1), (P2), (P4)
4000

2000

10

Original
S(t)
I(t)
V(t)
10

1990

5

5

10

1970

1980
Year

0

1930

1970

Count

0

5000

Original
I(n)

6000
Count

Count

World war II (1945)

2000

Er= 2.23e+02 (7.0e+04, 0.52, 0.47, 0.13) , P=(0.1, 50), (3329, 0.01), k=88000

OriginalI(t)Original
I(t)

Count

Count

3000

10

Year

1960
Year

(h) 天然痘 (Smallpox)(P1), (P2), (P4)(i) 腸チフス (Typhoid fever) (P1), (P2), (P4), (P5)

1930

1935

1940

1945

1930

1940

1950

1970

1980

主要な 9 種の疫病に対する F UNNEL の学習結果．

Fitting results of F UNNEL for 9 diseases (global-level counts).

ケンスの長さ n を変化させている．図に示すように，提案アル体的には，(a) インフルエンザ (Influenza) と (b) クリプトスポリ

ゴリズムである F UNNEL F IT は，データの入力サイズに対し線ジウム症 (Cryptosporidiosis) について，AR が外部ショックや外形である．

れ値に強い影響を受けているのに比べ，提案モデルは長期的な

5. ディスカッション

予測に成功している．
一般性 - コンピュータウィルスの流行パターン. F UNNEL は，そアプリケーション - 疫病の将来予測. F UNNEL は，様々な種類の一般性と柔軟性により，疫病のみならず，コンピュータウィ

の疫病データを柔軟に表現可能なため，最も重要なアプリケー

ルス等の流行パターンも学習することができる．例えば，図 11

ションとして，将来予測が挙げられる．図 10 は，2 つの主要は，IPA（注 8） による日本国内の企業および教育機関におけるコな疫病に対する予測結果を示している．ここでは，オリジナル

ンピュータウィルス感染数の報告データに対する F UNNEL の学

データの 2/3 の長さ（黒線）を用いてモデルを学習し，その後習結果を示している．ここで，丸印はオリジナルの報告件数，の数年間を予測している（赤線）．ここで，y 軸は対数スケー

F UNNEL の学習結果は線で表現され，期間は 2000 年から 2010

ルで表示している．提案モデルの精度を検証するため，自己回年である．例えば，(a) “Badtrans”, “Klez” は Microsoft Outlook

帰モデル (AR: autoregressive model) との比較を行った．AR ののセキュリティホールを利用したウィルスで，2001 年から 2002

回帰係数はそれぞれ，r = 52 (1 年間), 26 (半年間), 8 (F UNNEL

年にかけて流行した．これらのウィルスは，強い感染力を持っての基本パラメータ数 (6)，および疫病減少効果パラメータ数 (2)いたが，アンチウィルスソフトウェアの普及により鎮火された．

の合計) とした．さらに，学習の発散を抑えるため，AR につい(b) “Netsky” は，メールの添付ファイルとして拡散されたウィルては，対数スケールに対し予測を行った．図 10 に示した通り，
提案手法は AR と比較し高精度の予測を達成している．より具（注 8）：IPA - IT security center: https://www.ipa.go.jp/security/english/index.html

図 6 ワクチンのライセンス化 [10] と F UNNEL の結果．

GlobalFig. 6 The year of vaccine licensure [10] vs. detection.
Measles

1963

1965

Mumps

1967

1975

Whooping cough (pertussis)

1948

1951

Rubella

1969

1972

1000
500

30

SIRS

SKIPS Funnel−RFunnel(a) Global fitting
図8

Fig. 8

Original
I(t)

5000

35

25
SIRS

10000

Original
I(t)SKIPS Funnel−R Funnel(b) Local fitting

F UNNEL と既存手法の精度比較．

Fitting accuracy for the global and local sequences.

Original
S(t)
I(t)
V(t)

5

10

0

1950
Year

1970
Original
S(t)
I(t)
V(t)

5

10

1950
Year

10

1970

(a) ニューヨーク州 (NY)
図7

1930

1950
Year

Count (log)

10 20 30 40 50
Number of diseases (d)

2500
2000

Fig. 9

10

20
30
40
50
Number of states (l)(b) States (l)図9(b) ペンシルベニア州 (PA)

3000
2000
1000
0
2000
4000
6000
Number of time−ticks (n)(c) Duration (n)F UNNEL F IT の計算コスト．

F UNNEL F IT scales linearly: wall clock time vs. dataset size.

Local-level fittings for measles.

6. む す び

5

10

0

10

本論文では，大規模疫病データのための非線形モデル解析手

FunnelAR(52)
1930 AR(26)
1935
AR(8)

1940
Year

1945

法として F UNNEL について述べた．F UNNEL は，大規模な疫病

1950

データの中から，季節性，ワクチン効果，地域性，外部ショッ(a) インフルエンザ (Influenza)

クイベントや入力エラー等の重要な要素を自動抽出し，長期的

proposed=2.62e+01, AR(52,26,8)=1.19e+02, 8.49e+01, 9.03e+01

Count (log)

0

1970

proposed=8.69e+03, AR(52,26,8)=9.62e+03, 9.76e+03, 9.77e+03

な将来予測の能力を有する．様々な種類の疫病データを用いて2実験を行い，F UNNEL が最新の疫病解析手法と比べてより高い

10

0

10

FunnelAR(52)
AR(26)
2007
AR(8)精度と性能を持つことを示した．
2008

2009
Year

2010

2011

文

2012

(b) クリプトスポリジウム症 (Cryptosporidiosis)
図 10 2 つの疫病に対する将来予測の結果．
Fig. 10 Forecasting result: we train the model parameters using 2/3 ofeach sequence (i.e., solid black lines). We then start forecasting (atthe vertical dotted line).

3000

Badtrans

Netsky2000 Sircam

Mytob

Klez

1000

図 11

1000

3000

麻疹のローカル（州）の疫病に対する学習結果．

Fig. 70

2000(a) Diseases (d)

0

1930

3000

Wall clock time (s)

0
1930

1970

Wall clock time (s)

1950
Year

Wall clock time (s)

5000

0
1930

Count

RMSE

licensure detected

10000

10

1500
RMSE

Disease

Local
40

2001

2002

2003

2004 2005 2006
Year (per month)

2007

2008

2009

2010

コンピュータウィルスにおける F UNNEL の学習結果．

Fig. 11 F UNNEL is general: our model (solid lines) fits computer virus data(in circles) very well.

スで，2004 年に大規模な感染が報告された．このウィルスは，

10 年間にわたり，徐々に減少傾向にあるが，依然として感染報告がある．さらに，このウィルスには，週単位の周期性が確認され，週末の感染数が少ないという特徴がある．(c) “Mytob”
は，社内ネットワーク上（ローカル）で伝搬するウィルスである．さらに，最近になり猛威をふるっているのが，“Koobface”,“Fbphotofake” などに代表されるソーシャルネットワークサイト
上でのウィルスである．これらは，Facebook や Twitter 等で急速に伝搬し，潜在的な感染者数（ユーザ数）も多い．

献[1] Promotion of healthy swimming after a statewide outbreak of cryptosporidiosis associated with recreational water venues–utah, 20082009. MMWR Morb Mortal Wkly Rep, 61(19):348–52, 2012.
[2] C. Böhm, C. Faloutsos, J.-Y. Pan, and C. Plant. Ric: Parameter-freenoise-robust clustering. TKDD, 1(3), 2007.
[3] D. CC. Smallpox in the united states: It’s decline and geographic
distribution. Public Health Reports, 55(50):2303–2312, 1940.
[4] Y. Matsubara, Y. Sakurai, W. G. van Panhuis, and C. Faloutsos. FUNNEL: automatic mining of spatially coevolving epidemics. In KDD,
pages 105–114, 2014.
[5] F. NM, G. AP, and B. RM. Ecological and immunological determinants of influenza evolution. Nature, 422(6930):428–33, 2003.
[6] F. PE and C. JA. Measles in england and wales–i: An analysis offactors underlying seasonal patterns. Epidemiol, 11(1):5–14, 1982.
[7] D. SF. Seasonal variation in host susceptibility and cycles of certain
infectious diseases. Emerg Infect Dis., 7(3):369–74, 2001.
[8] M. SM, E. RJ, M. A, and M. P. Seasonality in six enterically transmitted diseases and ambient temperature. Am J Trop Med Hyg., 2014.
[9] L. Stone, R. Olinky, and A. Huppert. Seasonal dynamics of recurrent
epidemics. Nature, 446:533–536, March 2007.
[10] W. G. van Panhuis, J. Grefenstette, S. Y. Jung, N. S. Chok, A. Cross,
H. Eng, B. Y. Lee, V. Zadorozhny, S. Brown, D. Cummings, and
D. S. Burke. Contagious diseases in the united states from 1888 tothe present. NEJM, 369(22):2152–2158, 2013.