隠れマルコフモデルを用いた軌跡データに対するカーネルの設計

DEIM Forum 2015 P2-1

隠れマルコフモデルを用いた軌跡データに対するカーネルの設計

岩館 洸太

^†

和佐 州洋

^†

有村 博紀

^†

†

北海道大学大学院情報科学研究科 〒

060–0814

北海道札幌市北区北

14

条西

9

丁目

E-mail: †{k-iwadate,wasa,arim}@ist.hokudai.ac.jp

あらまし 本稿では，道路ネットワーク上を移動する二つの移動体の軌跡データ間の新しい類似度指標を提案する．

そのために，道路ネットワーク上の物体の移動軌跡に関する生成モデルとして，二次元平面上の軌跡に確率を割り当 てる，疎軌跡隠れマルコフモデル

(疎軌跡HMM)

を導入する．疎軌跡

HMM

は，隠れ変数として道路ネットワーク内 で通過したリンク列をもち，リンク列を生成するリンク間の遷移確率と，物体の真の位置から

GPS

の観測点を生成す る観測確率の組である．ここでは，疎軌跡

HMM

を用いたリンク列の推定と，2 つの軌跡データにおける類似度の計 算について議論する．とくに，この疎軌跡

HMM

に周辺化カーネル

(Tsuda, Kin, Asai, Bioinformatics, Vol.18, 2002)

を適用することで，軌跡データのための類似度として，軌跡カーネルを与える．

キーワード 疎な軌跡データ，隠れマルコフモデル，GPS データ，系列カーネル

1. ^{は じ め に}

1. 1

^{背 景}

近年，センシング技術の発展にともない，膨大な量の機器か ら多種多様なデータを収集することが可能になった．特に，自 動車や人などの移動体から大域位置計測システム（

GPS, global positioning system

）などを用いて得られたセンサーデータは，

軌跡データと呼ばれる．軌跡データは，軌跡予測

[3, 11]

や，

交差点のターン予測

[4]

，旅行時間予測

[13, 14]

，移動形態予 測

[8, 12]

，生活習慣解析

[2]

など，様々な応用があり，注目さ れている．軌跡データの特徴として，その量の膨大さや，欠損，

誤差を含む実数値であることなどが挙げられ，そのままでは扱 いにくいデータである．そこで，このような大量でかつ高次元 なデータを，効率良く扱うための手法の開発が望まれている．

1. 2

研 究 目 的

本稿では，軌跡データの確率的生成モデルとして，疎軌跡隠 れマルコフモデル（疎軌跡

HMM

）を設計し，疎軌跡

HMM

に 基づいた軌跡データのための適切なカーネル関数の設計につい て考察する．一般に，実数値座標の列からなる軌跡データのよ うな高次元データに対して，単純に類似度を測るのは，次元の 呪いの影響を受けるため，望ましくない．そこで，この問題を 解決するために，連続データに対する離散カーネルに注目する．

カーネル法（

kernel method

）

[9]

は，機械学習やデータマイニ ングにおいて，高次元データを効率よくとりあつかうための技 術であり，高次元属性空間のデータ間の距離を，属性ベクトル 間の内積として表すことで，指数的な数の組合せ属性を効率よ く暗黙に扱うことができる．

交通情報処理で良く用いられる軌跡のネットワークモデル

[5]

においては，移動体は与えられた道路ネットワーク上を移動す ると仮定する．ネットワークは，リンクと呼ばれる短い区間に 一定間隔で区切られ，軌跡データを各観測点の列に対応するリ ンクの列で表現する．単純なカーネルである軌跡計数カーネル

は，どのリンクを何回通過したか，その頻度を取り，比較する ことで類似度を計算する．しかし，単にリンクの出現回数を比 較するだけでは，リンク列が少しずれるだけで，軌跡データ間 の類似度が非常に低くなってしまい，現実的ではない．

また，現実のデータには，値の欠損や観測誤差があり，観測 データの各座標から，単純に対応するリンクを割り当てるだけ では，軌跡データから移動体が通ったリンク列を推定すること は難しい．とくに，マップ照合の精度が低い場合に，信頼でき ないリンク列を確定的に割り当ててしまう．したがって，得ら れた軌跡データに対して，欠損や誤差を考慮して，背後にある リンク列を推定する必要がある．

1. 3

主 結 果

前述の問題を解決するために，本稿では前述の軌跡計数カー ネルを，

Tsuda

らの周辺化カーネル

[10]

を用いて拡張し，

2

次 元平面上の軌跡データから確率モデルで周辺化した軌跡計数 カーネルを計算する手法を与えた．この周辺化によるソフト化 の効果で，データの背後にある確率的構造を取り込み，多少の ずれによる類似度の低下を防ぐことができる．

また，疎な軌跡データを生み出すリンク列の生成モデルとし て，疎軌跡隠れマルコフモデル（疎軌跡

HMM

）を導入した．

ここに，疎な軌跡データとは，欠損や誤差を含む軌跡データで

ある．このモデルによるマップ照合法として，軌跡データの生

成確率を最大化させるリンク列を計算する

Viterbi

アルゴリズ

ムを与えた．以下では，疎軌跡

HMM

に対する

Viterbi

アルゴ

リズムを，疎軌跡

HMM

マップ照合法と呼ぶ．さらに，周辺化

軌跡計数カーネルを計算するために，軌跡データの背後にある

リンク列の事後確率を与える

Forward-Backward

アルゴリズム

を提案した．計算機実験では，提案した疎軌跡

HMM

マップ照

合法の，疎な軌跡データに対する性能評価実験を行い，その有

用性を確認した．これらの結果は，筆者らの知る限り，実数座

標列としての軌跡データから，直接，離散カーネルを計算する

初めての結果である．

1. 4

関 連 研 究

本稿であつかった概念は，先行研究

[5, 13, 14]

に負うところ が大きい．井手と加藤ら

[14]

と，井手

[13]

は，軌跡データに 対する回帰問題に関する先駆的な研究を行なっている．彼ら

は

[13, 14]

，リンク列として与えられた，ネットワーク上の軌跡

のコスト（移動時間）の予測問題に対して，カーネル回帰の枠 組みを与えている．

[14]

では，系列カーネルに基づく軌跡カー ネルを提案しており，

[13]

では，正則化のための二つのリンク の類似度として，ネットワーク上の最短距離に基づく類似度を 導入している．しかし，いずれも軌跡の前処理としてマップ照 合を仮定しており，平面上の点列としての軌跡データは扱わな い点で，本稿のカーネルとは異なる．

本稿でも考察しているネットワークモデル上の隠れマルコ フモデルを用いたマップ照合は，

Krumm [7]

らによって提案 され，

GIS

分野で広く使われている．彼らのアルゴリズムは，

マップ照合における精度としては，問題ないが，疎な軌跡デー タを扱う際に，とびとびの観測点の間を，最短経路で近似する という経験則的工夫がなされており，厳密には，隠れマルコフ モデルとはなっていない．

Bernstein

と

Kornhauser [1]

は，軌 跡データの各座標を，幾何的に近いリンクに置き換えることで，

リンク列を推定するマップ照合アルゴリズムについて提案して いる．彼らのアルゴリズムは，一般道の上を高速道路が通過し てるケースや，丁字路で交差点付近で誤差があるケースなどで，

マッピングがずれるといった欠点がある．

2. ^{準 備}

本節では，本稿で扱う基本的な用語の定義を述べる．まず，

軌跡データについての定義を与える．次に，道路ネットワーク についての定義を与える．

2. 1

軌跡データ

本稿では，空間の領域として，

2

次元平面

D=R²

^{を考える．}

D

の任意の要素

p= (x, y)

を位置（

location

）または点

point

といい，

p.x=x, p.y=y

を

x

座標と

y

座標とする．今，

2

次 元平面

R²

^{上を動き回る移動体（}

moving object

）

o

を考える．

この移動体の位置をある時間間隔ごとに，観測して得られる点 列

x= (x1, . . . , xm)∈D^m

を，移動体

o

の軌跡（

trajectory

） とよぶ．このとき，軌跡

x

の長さを

|x|=m

と定義する．各 系列の長さ

m

は異なっていて良いものとする．

以下では、

D

上の点

q∈D

を長さ

1

の点列とみなし、空列

ε

を 長さ

0

の点列とみなす。このとき、列

y= (y1, . . . , yn)∈D∪{ε}

に対して、

⊕

y=y1· · ·yn∈D^∗

で

y

中のすべての点列を連 結して得られる点列を表す。

2. 2

道路ネットワーク

本稿では，軌跡データの生成について，移動体の動きは次の ように定義される道路ネットワークモデル上に制約されると仮 定する．ただし，実際の観測点は，観測誤差などから必ずしも 道路ネットワーク上にある必要はなく，

D

上の任意の場所をと りうる．

OSM

（

Open Street Map

）

^（注1）

などを一般化し，道路ネット ワークデータの定義を以下のように与える．

道路ネットワーク（

road network

）は，各頂点が

D

上の位 置に対応する平面埋め込みグラフ

N = (V,L, loc)

であり，次 を満足するものとする：

•

頂点集合

V={1, . . . ,|V|}

^．各頂点

v∈ V

^{は交差点や道} 路の接続点等の結節点を表す．

•

^{リンク集合}

L={1, . . . ,|L|}⊂=V²

．各リンク

ℓ∈ L

^に対 して，その両端点の頂点を

ℓ.s

と

ℓ.t

で表す．

•

位置関数

loc:V →R²

は，ネットワークの各節点

v∈ V

に対して，平面上の点

loc(v)∈D

をその位置として割り当てる．

ネットワークのサイズ（

size

）をそのリンク数

L=|L|

^と定 める．以後，簡便のため，ネットワーク

N

^{とリンク集合}

L

^を 同一視し，リンク

ℓ

に対して，

ℓ∈ N

のように書く．さらに，

ある点

µ∈R²

^がリンク

ℓ

の表す線分上にあることを

µ∈ℓ

と 書く．

3. 疎な軌跡データのための確率的生成モデル

本節では，疎な軌跡データの確率的な生成モデルである疎軌 跡隠れマルコフモデル（

Sparse Trajectory HMM

） （または 疎軌跡

HMM

）を導入する．これは，交通情報処理（

ITS

）分 野や地理情報処理（

GIS

）分野であつかわれるネットワークの リンクを状態とする隠れマルコフモデルに，観測のオンオフを 表す隠れ変数を追加したものである．

ランダムウォークの観点からは，このモデルは，移動体がネッ トワーク上のランダムウォーク（乱歩）を行いながら，ときど き，その位置を二次元平面上の点として，確率的に出力するモ デルとみなすこともできる．

ネットワークリンクを状態とし，実数平面上の出力分布をも つ隠れマルコフモデルにおいて，出力状態と未出力状態を設け ることで同様の生成モデルを構築できるが，ここでは明示的に 観測のオンオフを表すを表す隠れ変数を設ける選択をしている．

これにより，次節で考察する動的計画法に基づく推定アルゴリ ズムの構築が容易になる．

今，任意の正整数

m, n(1<=m <=n)

を考える．任意の添え 字

i∈[1, n]

を時刻という．とびとび

HMM

^{の観測変数は，長} さ

m

の点列である入力点列

x= (x1, . . . , xm)∈D^m

である．

各

xi∈D(i= 1, . . . , n)

は平面

D

中の位置であり，

i

番目に観 測された位置を表す

.

疎軌跡

HMM

は、次の隠れ変数をもつ．

•

^{長 さ}

n

の リ ン ク 列（ 移 動 リ ン ク 列 ま た は 移 動 経 路 ）

h = (h1, . . . , hn) ∈ Lⁿ.

移動体が動きまわるリンクの時系 列を表す．

•

長さ

n

のコイン投げの列（発火列）

b= (b1, . . . , bn)∈ {0,1}ⁿ.

各

i

に対して，

bi∈ {0,1}

^{は表が出る確率が}

β∈[0,1]

のコインである．

•

^長さ

n

の出力の列である隠れ出力系列

y= (y1, . . . , yn)

． 各出力

yi∈D∪ {ε}

^は，時刻

i

における観測を表し、

（注

1）：http://openstreetmap.com/

– bi= 1

のとき，位置

yi=q∈D

，または

– bi= 0

のとき，空列

yi=ε

のいずれかを値としてとる．最終的な観測点列は、

y

が含む出 力を時間順に連結して得られる点列

x=⊕

y=y1· · ·yn∈D^∗

として定まる。得られた列

x

はちょうど

m

個の点を含むもの とする．

•

^長さ

m

の時間列である発火時間列

t= (t1, . . . , tm)∈ [0,∞)^m

．各

i= 1, . . . , n

に対して，

i

番目の発火時間，すなわち，

位置

xi∈D

を出力した時刻は

ti= min{t∈[1, n]|∑_i

j=1bj>= i}

で定められる．これは発火列から決まる．

上記の値が与えられた場合に，同時経路対（

joint path

）は，

対

z=⟨h,b⟩

である．数値パラメータ

m, n

と入力列

x

が固 定された場合，同時経路対

z

は正確に上記の隠れパラメータす べてを一意に定める。非負整数

i

に対して，

z1:i=⟨h1:i,b1:i⟩

を表す．

各確率変数の確率分布を次のように与える．始点と終点の事 前確率は，それぞれ，

p(s|Θ)

と

p(t|Θ)

で与えられる．遷移確 率（

transition probability

）は式任意のリンク対

k, ℓ∈ L

^に対 して，リンク

k

からリンク

ℓ

への

p(hi=ℓ|h_i−1 =k, A) =Akℓ

で定義する．ここに，遷移行列

A= (Akℓ)

は，サイズ

|L| × |L|

の実数行列である．行列の各行は

1× |L|

^{の確率ベクトルであ} り，

∑

ℓAkℓ= 1

を満たすものとする．リンク

h

と

h^′

が隣接し ているとき，そのときのみ，

Ahh′ >0

とする．

発火列

b

は，

p(b|m, n,Θ)

は，ちょうど

m

個の

1

を含むコ イン投げ全体上の一様分布に従うとする．

移動体の現在の位置

µ∈D

に関して，その観測位置

x∈D

は，中心

µ

と分散

σ

の

2

次元正規分布

p(x|µ, σ) = 1 2π|Σ|^1/2 exp

{

−1

2(x−µ)^TΣ⁻¹(x−µ) }

に従うものとする．

(3.)

をリンク

hi

上で積分すると，リンク

hi

上での観測点

xi=x

の累積生成確率は

p(xi=x|hi=h, σ) =

∫

µ∈h

p(x|µ, σℓ)dµ

となる．しかしこの積分計算は計算を要するので，簡便のため 本稿では，リンク

hi

の中点

ci

に対して，中心

µ=ci

と分散

σ

をもつ正規分布

N(c, σ)

を考えて，これで

p(xi=x|hi=h, σ)

を近似する．

以上の準備の下で，パラメータ

x,b,h

の同時確率分布は次 の式で与えられる

:

p(x,b,h|m, n,Θ)

= p(x|b,h, m, n,Θ)p(b|m, n,Θ)p(h|n,Θ)

= p(s=h1.s|Θ)p(t=hn.t|Θ)

∏n i=1

p(hi|hi−1, A)p(bi|Θ) (β p(xi|hi, σ))^bi(1−β)^(1−bi).

入力列

x

の周辺確率分布は次の式で与えられる．

p(x|m, n,Θ)

= ∑

h

∑

b

p(x,h,b|m, n,Θ)

= ∑

h

p(h|n,Θ)∑

b

p(b|m, n,Θ)p(x|h,b, m, n,Θ)

次節と次々節では，疎な軌跡隠れマルコフモデルに対して，

各種の推定アルゴリズムを与える。

4. ^疎軌跡 HMM のための最尤経路推定アルゴ リズム

本節では，入力系列から最尤経路推定を行うための，

HMM

でいうところのいわゆる

Viterbi

アルゴリズムに相当するアル ゴリズムを与える．

とびとび

HMM

から観測された点列

x

に対して，事後確率

p(x,h,b|Θ)

が最大となる

h

と

b

の組を求めるための

Viterbi

アルゴリズムと，求めた

h

を用いて

p(h|x)

を求めるための

Forward-Backward

アルゴリズムを提案する．以下では，ユー ザーが与える移動体の走行時間の最大を

n

^{とおく．観測点数}

m=|x|

^は，

m <=n

であることに注意されたい．

略記法として，出力コイン

b∈ {0,1}

^に対する

k

番目の観測 点の

xk

の出力確率を次のようにおく．

q(xk|h, b) =





(1−βh) (ifb= 0), p(xk|h, σ)βh (ifb= 1)

つまり，

q(xk|h,0)

はリンク

h

上で何も観測されない確率であ り、

q(xk|h,0)

は観測点

xk

が観測される確率である．

各時点

i= 1, . . . , n

と，全てのリンク

h∈ L

^{，可能な累積出力} 数

0<=k <= min{i, m}

^の組合せ

(i, h, k)∈[1..n]×[1..L]×[0..m]

の全てに対して，実現可能な最大同時確率

MP[i][h][k] = max

z1:i h_i=h,|b1:i|1=k

p(x,z1:i|Θ)

= max

h1:i h_i=h

max

b1:i

|b_1:i|1=k

p(x,h1:i,b1:i|Θ)

を計算することを考える．ここで，過去の履歴

h1:i ∈ Lⁱ

^は，

末尾がちょうど

hi=h

となるリンク列の集合全体をとし，履 歴

b1:i ∈ {0,1}ⁱ

は，これまでの出力回数

♯b1:i =∑i

j=1bj

が ちょうど

k

となるスイッチ列の集合全体とする．

MP

の各要素の計算は，次の基底ケースと帰納ケースを用い て，動的計画法で行う．

(i)

基底ケースは次のとおりである：任意の

h ∈ L, k ∈ {0, . . . , m}

^{に対して，}

MP[0][h][k] =





1 (k= 0), 0 (k >0).

(1)

(ii)

帰納ケースは，任意の時点

i >1

に対して，次のとおり である：

MP[i][h][k]

= max

h′ max

b p(h|h^′)q(x|h, b, k)MP[i−1][h^′][k−b].(2)

Algorithm 1

疎軌跡

HMM

に対する

Viterbi

アルゴリズム

1: procedureSleepWakeViterbi(x)

INPUT:

観測点列

x;

OUTPUT:x

に対応するリンク列

h;

2: MP

の初期化;

▷

式

(1)

3: fori= 1, . . . , ndo 4: forh= 1, . . . , Ldo 5: fork= 0, . . . , mdo

6: MP[i][h][k]

の計算;

▷

式

(2) 7: (h_∗,b_∗) = arg maxhmaxbMP[n][h][m];

8: return(h_∗,b_∗)

ここに，

max

は任意の

h^′ ∈ L,b∈ {0,1}

^{に関してとるものと} する．これらの更新式を利用した提案の

Viterbi

アルゴリズム の擬似コードを，

Algorithm 1

に与える．

求めた

MP

に対して，リンク列の長さ

n

を与えた場合に，最 大事後確率は，

M P= min

h∈LMP[n][h][m]

で求めることができる．

次に，

M P

の計算コストを考えよう．可能な組

(i, h, k) ∈ [1..n]×[1..L]×[0..m]

の総数は

n×L×(m+ 1)

通りあるので，

一見すると動的計画法の計算は

O(nmL)

^{時間を要すると思わ} れる．しかし，各

i, h

に対して，実際にとれる

k

の値は限られ る．そのため，小さい計算量で実行できる．具体的には，全体 で

n+ 1

回のステップ

i= 0, . . . , n

に対して，各ステップ

i

で

L×2

通りの計算を行い，最小値をとると，全体では

O(nL)

時 間と領域で

M P

は計算可能であることがわかる．以上の議論 から次の補題が得られる．

［補題

1

］

M P

の計算量は，

O(nL)

時間と

O(nL)

領域である．

MP

の各要素の値は，

n

が大きくなるにつれ，非常に小さく なる．そこで，前節の提案の

Viterbi

アルゴリズムで用いられ た

MP

の代わりに，その負対数

LP[i][h][k] =−log(MP[i][h][k])

に書き換えることで，計算機でのアンダーフローを回避する．

このとき，

LP[i][h][k]

を計算するための基底ケースと帰納ケー スは，それぞれ，次で与えられる：

(i)

基底ケース：任意の

h∈ L,k∈ {0, . . . , m}

^{に対して，}

LP[0][h][k] =





0 (k= 0),

∞ (k >0).

(ii)

帰納ケース：任意の時点

i >1

に対して，

LP[i][h][k]

= min

h′ min

b

( −logp(h|h^′) −logq(xk|h, b) +LP[i−1][h^′][k−b]

) .

5. ^疎軌跡 HMM のためのリンク訪問確率推定 アルゴリズム

本節では，固定された観測点列

x

に対して，ネットワーク上

の経路が指定されたリンク列

h

を含む事後確率を求める，いわ ゆる

Forward-Backward

アルゴリズムを与える．

5. 1

概 要

入力列

x

を固定した場合に，長さ

n

の同時経路

z=⟨h,b⟩

の全体を考えて、移動体が時点

i

で指定されたリンク

h

を訪問 する事後確率は

P roblink(hi=h|x,Θ)

= ∑

z h_i=h

p(z|x,Θ)

= ∑

b₁+b₂=m







∑

z1:i h_i=h

p(z1:i|x,Θ)













∑

zi:n h_i=h

p(zi:n|x,Θ)







ここで，後半の同時系列の名前を変えて，

z^′=⟨h,b⟩

^とす る．前半の

z

と後半の

z^′

は、時刻

i

で接続する必要があるの で，制約

zi=z^′_i

を付加する．

= ∑

b1 +b2 =m z_i=z^′_i







∑

z1:i hi=h

p(z1:i|x,Θ)

















∑

z′i:n h^′_i=h

p(z^′i:n|x,Θ)











ここに略記として，

b1=|b1:i|1

および

b2=|bi+1:n|1

と書いた．

すると，先の接続制約

zi=z^′i

を

hi=h^′i

と

bi=b^′i

に置き換 えることができるが，前者はさらに

h

に等しいので，前半と後 半を完全に分離することができる．

= ∑

1<=k<=i b∈{0,1}

















∑

z1:i hi=h

bi=b

|b_1:i|1=k

p(z1:i|x,Θ)







































∑

z′ h′i:n

i=h b′i=b

|bi+1:n|1=m−k

p(z^′i:n|x,Θ)























= ∑

1<=k<=i b∈{0,1}

P robf w(i, k, b|h,Θ)·P robbw(i, m−k, b|h,Θ)

= ∑

1<=k<=i b∈{0,1}

FP[i][k][b][h]·BP[i][m−k][b][h] (3)

となる．

ここに，任意のリンク

h

と，ビット

b

，時刻

1<=i <=n

に対 して，前進確率（

forward probability

）とは，リンク

h

で終わ り，ちょうど

k

回発火したようなネットワーク上の長さ

i

の同 時経路対

z1:i=⟨h1:i,b1:i⟩

すべてについて，その生成確率を 足し合わせたものであり，

P robf w(i, k, b|h,Θ) = ∑

z1:i hi=h,bi=b

|b_1:i|1=k

p(z1:i|x,Θ) (4)

で与えられる．前進確率を

FP[i][k][b][h]

と書く．

同様に，

k^′ =m−k

とおくと，後退確率（

backward prob- ability

^{） は，リンク}

h

^{で始まり，ちょうど}

k^′

^{回発火したよう}

なネットワーク上の長さ

n−i

の同時経路対

z1:i=⟨h1:i,b1:i⟩

すべてについて，その生成確率を足し合わせたものであり，

P robbw(i, k^′, b|h,x,Θ) = ∑

z′ i:n h′

i=h,b′ i=b

|b_i+1:n|1=m−k

p(z_i:n^′ |x,Θ) (5)

以後，後進確率を

BP[i][k^′][b][h]

^と書く．

5. 2

アルゴリズム

前進確率と後退確率を求める動的計画法を用いた効率良いア ルゴリズムを与える．初めに，前進確率

FP[i][k][b][h] = ∑

z1:i hi=h,bi=b

|b_1:i|1=k

p(z1:i|x,Θ)

の計算を考える．計算のために，

Q(xk, h, h^′, b) =p(h|h^′, A)q(xk|h, b)

とおくと，漸化式の基底ケースと帰納ケースは次のように書け る．任意の

1<=i <=n

^，

1<=k <=i

^，

b= 0,1

^{を考える．}

(i)

基底ケース：任意の

h∈ L

^{に対して，}

FP[0][k][b][h] = 1 (k= 0, b= 0),

FP[0][k][b][h] = 0 (

それ以外

). (6) (ii)

帰納ケース：任意の時点

i >1

に対して，

FP[i][k][0][h] = ∑

h′

Q(xk, h, h^′,0)FP[i−1][k][h^′], FP[i][k][1][h] = ∑

h^′

Q(xk, h, h^′,1)FP[i−1][k−1][h^′], FP[i][k][h] = FP[i][k][0][h] +FP[i][k][0][h]. (7)

次に，後退確率

BP[i][k^′][b][h] = ∑

z′ h′ i:n

i=h,b′

i=b

|b_i+1:n|1=m−k

p(z_i:n^′ |x,Θ)

の計算を考える．計算のために，

R(xk, h, h^′, b) = p^R(h|h^′, A)q(xk|h, b)

とおく．さらに，逆向きの遷移確率行列

B = (Bh^′,h)^L_h′,h=1

を 次のように計算する．

B_h′,h = p^R(h|h^′, A) = p(hi−1 =h|hi=h^′, A)

= p(hi−1 =h, hi=h^′|A) p(hi=h^′|A)

= Ah,h′

∑

hA_h,h′.

上記の式を用いることで，漸化式の基底ケースと帰納ケース は次のように書ける．任意の

1<=i <=n

，

1<=k <=i

，

b= 0,1

を考える．

(i)

基底ケース

:

BP[n+ 1][k][b][h] = 1 (k=m, b= 0),

BP[n+ 1][k][b][h] = 0 (

それ以外

). (8) (ii)

帰納ケース

:

Algorithm 2

疎軌跡

HMM

に対する

Forward-Backward

ア ルゴリズム

1: procedureForwardBackward(x,h) INPUT:

観測点列

x,

同定されたリンク列

h;

OUTPUT:h

の事後確率

p(h|x);

2: FP

と

BP

の初期化;

▷

式

(6), (8) 3: fori= 1, . . . , ndo

4: forh= 1, . . . , Ldo 5: fork= 0, . . . , mdo

6: FP[i][h][k]

の計算;

▷

式

(7) 7: fori=n, . . . ,1do

8: forh= 1, . . . , Ldo

9: fork=m, . . . ,0do

10: BP[i][h][k]

の計算;

▷

式

(9)

11: p(h|x)

の計算;

▷

式

(3)

12: returnp(h|x);

BP[i][k][0][h] =∑

h′

R(xk, h, h^′,0)BP[i

＋

1][k][h^′], BP[i][k][1][h] =∑

h′

R(xk, h, h^′,1)BP[i

＋

1][k+ 1][h^′], BP[i][k][h] =BP[i][k][0][h] +BP[i][k][1][h]. (9)

上記の前進確率と後退確率の計算式を

(3)

と合わせると，リ ンク訪問確率を計算できる．これを

Forward-Backward

アルゴ リズムと呼ぶ．擬似コードを

Algorithm 2

に与える．

［補題

2

］

Forward-Backward

アルゴリズムは，リンク訪問確 率

P roblink(hi =h|x,Θ)

を

O(nm|L|)

時間と 領域 で計算す る．ここに，

m

は観測点列の長さであり，

n

は隠れリンク系列 の長さ，

|L|

はネットワーク中のリンク数である．

6. 疎な軌跡データに対する系列カーネルの設計

本節では，疎な軌跡データのためのカーネルを提案する．

6. 1

軌跡データに対する同時カーネル

はじめに，以下で導入するリンク系列計数カーネルは，文字 列に対する計数カーネル（

count kernel

^）

[10]

^{を軌跡データに} 拡張したものである．

今，長さ

m

の観測点列

x∈D^m

に対して、それを生成した 長さ

n

のリンク列

h∈ Lⁿ

^{を考える。この}

h

は、マップ照合 アルゴリズムで決定的に計算されたとしても良いし、確率的に 推定されたと考えても良い。

このとき，リンク

ℓ

に関するリンク系列

h

の正規化リンク計 数（

normalized count

）を

cℓ(h) = 1 n

∑n i=1

I(hi=ℓ) (10)

とおき，リンク計数の特徴ベクトルを

ϕ(h) = (c1(h), . . . , c_|L|(h))∈R^|L| (11)

とおくと，リンク系列計数カーネル（

link sequence count ker- nel

^）は式