普遍クリギングを用いた最短旅行時間予測三浦英俊

1. はじめに本研究は 2 地点間の最短旅行時間予測にクリギングPBを適用する．図 1 のように相互の距離が近い平面上の 3

80

点 PA ，PB ，PC ならびに PX ，PY ，PZ があるとする．

PC から PZ まで自動車で移動するときの最短旅行時間を予測したい．PA から PX までおよび PB から PY まPX

60
PC?
PAでの最短旅行時間データがあるならば ，最短旅行時間PZ
PYが直線距離におおよそ比例することに加えて，これら4 地点と PC ，PZ 相互の距離が推測の手がかりとなる．

図 1: 未知の旅行時間の推測

クリギングとはこのような地理データに基づいた推測を統計学の枠組みで展開した手法である．距離が近速度は道路の等級や幅員あるいは通過交通量に応じていほど 似た値を示す性質をもつ空間データに対して，既定められることが多いが，機械的に定めると都市部な知観測データの空間分布に基づいて未知デ ータを推定

どで実態にそぐ わない予測値が出ることがある．これすることができる．クリギングはもともと地下の鉱床を避けるためには 詳細な速度調査に基づく設定が必要の空間分布を推定する鉱床学の手法とし て発展し ，近になる．

年では空間分布を取り扱う分野で広く応用されている

本研究で述べるクリギングを利用した最短旅行時間

11]．都市地域分野では地価の予測や細街路の自動車速度の推定など の既往研究がある 1]3]5]7] ．

予測は，道路ごとに 速度を設定した道路網データの代わりに，複数の地点間の最短旅行時間観測データを使

本研究では，出発地と目的地の位置を表すそれぞれ用する．従って予測精度を高めるためには信頼性の高の座標によって構成される 4 次元空間においてクリギい観測データを数多く用意しなければならない．もしングの一手法である普遍クリギングを適用する．PA →観測データに最短とは言えない旅行時間データが多数

PX ，PB → PY ，PC → PZ それぞれの最短旅行時間経路は，経験的に PA ，PB ，PC および PX ，PY ，PZ が含まれているならば，予測値は最短旅行時間とし て信用できないからである．

相互に近いほど 経路の重複部分が増える．迂回の程度

先に少し 触れたようにクリギングの適用範囲は幅広

や移動速度も類似し ，従って両経路の旅行時間は似かいが，旅行時間の予測に適用した既往研究に限ればそよった値になる．このことは，平面上の 2 地点間の移れほど 多いとは言えない．Braxmeier, H., Schmidt, V.

動を出発地と目的地の位置座標で構成される 4 次元空

and Spodarev, E. 1] および Du, J. and Aultman-Hall,
L. 3] は，GPS によって計測された道路の自動車走行間における点で表すならば，点間の距離が近いほど 似た旅行時間になると言い換えられる．それゆえ，この速度を観測データとして使用し ，データが取得できな

4 次元空間でクリギングの利用による最短旅行時間予い道路の速度をクリギングによって予測する方法を提案している．ただしこれらは 4 次元空間でクリギング測が可能であると言える．

2 地点間の最短旅行時間を予測するためには，道路を適用した研究ではない．距離の近い道路は似たよう網ネットワークデータにダ イクストラ法など の経路探な走行速度になることを仮定して，平面上で推定が行索を適用する方法が一般的である．精度の良い予測をわれている．

行うためには各道路に付与される速度が重要である．

2. 普遍クリギングによる旅行時間の予測
本節では Cressie 2] と間瀬，武田 11] に依拠しつつ，みうらひでとし 明海大学不動産学部〒 279-8550 千葉県浦安市明海 1-2

普遍クリギングによる旅行時間予測値の導出方法につhmiura@meikai.ac.jp1

Pn  = 1 を満たすとする (ただし  は正値とは限i=1 iいて述べる．

らない)．パラメータベクトル

対象地域を平面領域と仮定する．この地域の出発地・目的地ペアの座標によって構成される 4 次元ユークリッ

二乗予測誤差ド 空間の部分集合を D とする．任意の出発地・目的地(1)

小化すべき目的関数を(  m) =C(0) + 0 ; 2 0 c ; 2m(2)(3)とする．((Xni=1)

i f(si ) ; f(s)(7)

 m) を  m に関して偏微分して得られる正規方程式を解いたならば ，E ((s)) = 定数，および (3) が成り立つとき，(s) は二次定常性を満たすという．これは任意の s 2D について (s) の平均が一定であり，相対的位置が h 離れたときの (s) と (s + h) の共分散がコバリオグラムと呼ばれる関数 C(h) で記述できることを意味する．

=

c + X (X 0 ;1 X );1 (f(s) ; X 0 ;1c) 0 ;1m = ;(f(s) ; X 0 ;1 c)0 (X 0;1X );1(8)(9)

コバリオグラムの関数型にはいくつかの定型が提案0を得る．ただし X = ff(s1 ) : : :  f(sn)g とする．z =されている．本研究では，コバリオグラムが等方的 (距fz(s1) : : :  z(sn)g0 として (8) を (5) に代入すると予測

離 jjhjj のみに依存する) であることを仮定し ，指数型

^ s) の計算式を得る．
値 Z(と呼ばれる次の関数形を使用する．

n^ s) = X i z(si) = zZ((C(h) = 0 +;jjh1jj =2 jjhjj = 0(4)
1 e jjhjj > 0．
ただし 0 ，
1 ，2 はパラ メータである．0  0，1  0，
2 > 0 とする．既往研究 (例えば 5]8]) では sphericali=1= ^gls f(s) + c0 ;1 z ; ^gls X ．(10)ただし^gls = (X 0;1X );1 X 0 ;1 z型と呼ばれる型が使用されることが多いが，これは対象空間の次元が 3 を超えると適用できない．以上の事は柄がクリギング適用に必要な仮定である．

(11)の一般化最小二乗推定量である．予測値の計算に必要な (1) のパラメータ

さて， n 組の出発地・目的地ペア間の旅行時間観測と (4) のパラメータ 0 ，
1 ，2は制限付き最尤推定法によって推定する 2]9]．詳しくデータが得られているとしよう．第 i 番目のデ ータのは 3.2 節で述べる．

出発地・目的地ベクトルを si 2D，旅行時間を z(si ) と^ s) の 95 ％信頼区間を得るには (8) を (6)
予測値 Z(する (i = 1 : : :  n)．クリギングによる任意の s に関す^ s) は観測データ z(si ) の重み付き
る旅行時間予測値 Z(に代入すればよい．すなわち

2 (s) =C(0) ; c0 ;1 c+ (f(s) ; X 0 ;1c)0 (X 0 ;1 X );1 (f(s) ; X 0 ;1 c)(12)

平均i=1(6)と置いた．(6) にラグランジュ乗数 m を導入して，最すと仮定する．

n

 = fC(si ; sj )gij ，c = (C(s ; s1) : : :  C(s ; sn ))0

2 地点間の距離， はパラ メータである．2 地点間の距離は必ずし も直線距離でなくてもよい．E() を平均，Cov( ) を共分散，h を任意の 4 次元空間ベクトル，コバリオグラムを C(h) として，(s) は次の性質を満た^ s) = X i z(si )Z(を最小にする値とし て求められる．ここでという構造を持つと仮定する．ただし f(s) は平面上のE ((s)) = 0Cov ((s + h) (s)) = C(h):= (1  ::: n)0 は平均^ s))22 (s) = E (z(s) ; Z(= C(0) + 0 ; 2 0 c

ペア s 2 D について，2 地点間の旅行時間 z(s) が，z(s) = f(s) + (s)i(5)によって得られる．ここで重みパラメータ i は2から 95 ％信頼区間h^ s) ; 1:96(s) Z(^ s) + 1:96(s)Z(iが求められる．

3. 適用例：千葉県内の自動車旅行時間予測
3.1. 使用データ
千葉県を対象地域とし て任意の 2 地点間の最短旅行時間予測を行う．図 2 の点は 2007 年 7 月現在千葉県にある市役所 (36ヶ所) を示す．北部は人口が多く道路密度も高いが，南部は平地が少なく人口は海岸地域に集中している．従って南部地域では道路と市役所は海の近くに分布しており，半島をよぎ る道路は比較的少ない．
図 2: 千葉県の 36 市分布図

図 2 の 36 市を出発地・目的地とする全てのペアを取り出す．これらを si (i = 1 : : :  1260(= 36 × 35)) とし ，それらの最短旅行時間を観測データ z(si ) とする．

z(si)

z(si)

f(si)

対象地域は平面ではないから座標と距離の扱いに少々注意を払う必要がある．市役所の緯度経度は国土数値
情報ダウンロード サービス 6] と 2 万 5 千分１地形図から取得した．si の出発地と目的地の合わせて 4 つの座

標値は，緯度と経度を平面直角座標系 (第ＩＸ系) に変換したものを使用する．出発地から目的地までの距離f(si)

f(si ) は，ユークリッド 距離ではなく，地球を 6370kmの球体と仮定して経度と緯度から計算される大円距離

図 3: 大円距離と最短旅行時間の散布図を使用する．
最短旅行時間 z(si ) はナビ タイムの車ルートサービ

布に従うと仮定する．まず制限付き最尤推定法によっ

ス (有料道路優先) によって取得した 10]．ナビ タイムて 0  1  2 を推定し (2]9])，次に (11) によって得ら

車ルートサービスとはインターネットで出発地から目れる

gls をの推定値とする．

0  1  2 を推定するための対数尤度関数は以下のと的地までの自動車による最短旅行時間と経路を検索できるサービスである．道路の設定速度は不明であるが，おりである．

結果の出力から等級・幅員・上り下り別に設定されていL(0  1 2 ) =; n ; 1 log 2 + 1 log jX Xj ; 1 log jj22211;
1; 2 log jX  Xj ; 2 Z Z ．
(13);1
ただし ， = ;1 ; ;1X (X ;1 X );1 X  とると推測される．観測デ ータの大円距離と最短旅行時間の散布図を図 3 に示す．図 3 において通常の最小二0乗法によって大円距離を説明変数とする切片が 0 の回0帰直線を求めたところ回帰係数は 2.147 であった．仮00にこの回帰直線によって最短旅行時間を予測したなら0ば，予測区間の 95 ％信頼区間は，f(s) の大きさにもよする．

3.2. セミバリオグラムの推定112，gls = 1:894 である．これらを当てはめたコバリオグラムを図 4 に示す．

得られたパラメータ推定値は 0 = 0 1 = 861 2 =るが 2.147f(s) を中心とした約 40 分の範囲である．
クリギングによって予測値を求めるためには，(1) のパラメータ

3.3.と (4) のコバリオグラム C(h) のパラメー

任意の 2 地点間の最短旅行時間予測

千葉県内の任意の 2 地点間の最短旅行時間を予測すタ 0  1 2 を導出しなければならない．(s) が正規分る．千葉県内の 3538 個の郵便番号から 2 つを 1 組とし

3θ0+θ1^ )Z(siθ0=0,θ1=861,θ2=112=θ1/θ2

図 4: 推定されたコバリオグラムの図示

z(si)て 100 組をランダムに選び出し 出発地・目的地ペアの

図 5: サンプルの最短旅行時間とクリギング予測値のサンプル ti (i = 1 : : :  100) とする．サンプルの郵便番号に対応する住所間を移動する最短旅行時間 z(ti) を散布図ナビタイムによって求める．z(ti ) とクリギング予測値^ ti) の散布図は図 5 のとおりである．図 6 は比較のたZ(め縦軸を 最短旅行時間 = 2:147 ×大円距離 による予測値に置き換えたものである．両図を比較すればクリギングによって予測の精度が上がっていることが見て取れる．決定係数は図 6 の 0.8138 から図 5 の 0.9354 へ上昇した．(12) によって予測値の 95 ％信頼区間を求めると，信頼区間の幅の半分 1:96(ti) の分布は図 7 に示すとおりである．その平均は 16.1 分であり，大円距離z(si)だけで予測した場合の約 40 分と比べて信頼区間の幅はかなり短くなる．

図 6: サンプルの最短旅行時間と回帰直線による予測

なお (5) から明らかなように，観測データである 36

市役所間の最短旅行時間 z(si ) はクリギングによる最短値の散布図

3.4. 観測データ数と信頼区間の関係

値 (0 = 0 1 = 861 2 = 112，gls = 1:894) を共通し

^ si ) と完全に一致する．
旅行時間予測値 Z(て使用した．

予測値の信頼区間を小さくするためには計算に使用

間や手間を考慮すれば 観測データをむやみに増やすこ

1:96(ti) の算術平均は k が大きくなるに従って小さくなることが見て取れる．k が 600 を超えるあたりかとはできない．

ら低下の程度は緩やかになる．

する観測データ数はなるべく多い方がよいが，計算時

3.1 節で用いた 1260 個の観測データからランダムにこの結果から，対象地域が都道府県程度の大きさであれば観測データは 600 程度使用すればよく，それをk 個の観測デ ータを取り出す．取り出された観測デ ータを使用して，3.3 節の 100 サンプル ti についてあらためてクリギング予測値を計算し，さらに 1:96(ti) の算術平均を求める．これを k=10, 50, 100, 200, 400,600, 800, 1000 について各 10 回繰り返し た結果を図8 に表す．また，すべての観測データを使用したとき(k = n = 1260) の結果も示す．ただしコバリオグラム
パラメータは 3.2 節で導出した n = 1260 のときの推定

超えて観測データを用意しても予測信頼区間は大幅に小さくすることはできないことが分かる．

3.5.

予測誤差の検討

3.3 節の結果について予測誤差を吟味する．3(ti) <^ ti) ; z(ti )j となる予測誤差の大きい 16 サンプルをjZ(取り出し，予測誤差が生じた原因について考察する．図

9 はこれら 16 サンプルの出発地と目的地を結ぶ線分を4(%)σ(t )[]

図 7: 1:96(ti) の分布
σ(t )[]

図 9: 予測誤差の大きいサンプルした事例計算において信頼区間は予測値を中心とする

16.1 分の範囲であったが ，ナビ タイムデ ータには現実の走行時間データに想定される計測誤差が含まれていないことを割り引いておく必要がある．
最短旅行時間の推定にクリギングを適用する利点のk

一つは，時間帯や曜日によって道路の速度が大きく変化する地域を通過する旅行時間の推定が容易にできる

図 8: 観測データ数と 1:96(ti) の算術平均ことである．当該地域のさまざ まな時間帯および 曜日の最短旅行時間観測データを十分に取得することがで図示したものである．最短時間経路はいずれも山岳地

きれば，クリギング 計算に使用する観測データを出発域など 道路の限られた地域を経由しているため経路はする時間帯に対応して交換することによって精度の高直線から大きく離れており，クリギング予測時間はナ

い最短旅行時間推定が可能となる．

ビタイムによる旅行時間よりもかなり小さい．南部地

予測精度向上のために道路ネットワークデータの利域における観測データは数が少ないうえに，市役所が用を検討すべきであろう．例えば 出発地から 目的地までの距離 f(s) とし て使用した大円距離を，最短距離経海岸に近い地域に偏っていることが誤差の原因である．
これらは市役所を代表点とした観測データに基づい路の道路距離に取り替えれば，もっと信頼区間を小さてクリギングを行う際に当然想定すべき誤差であった．

くできるかもしれない．

南部の内陸地域に観測データ端点を設定しておけば，予一方でクリギングによる予測では旅行時間に対応し測誤差はもっと小さくなったであろう．クリギングにた経路が分からないという大きな欠点がある．さしあ

よって最短旅行時間を予測するためには ，観測データたっては最短距離経路など で代用するにせよ，予測旅の数だけでなくそれらの地理分布にも注意すべきこと

行時間から経路を推定する方法を考えるなど 対策を講を示している．

じなければならない．
本研究の地図作製には MANDARA および Illustra-4. おわりにtor を使用した．計算は主として Mathematica を用いた．

本研究では最短旅行時間の予測にクリギングの応用を試みた．さまざまな制約から適用例における観測デー
タにはナビ タイムデータを使用した．千葉県を対象と

＜謝辞＞

5

本研究は科学研究費補助金「都市内・都市間交通網と公共施設配置に関する数理的ならびに実証的研究」(基盤研究 (A) 課題番号 17201037) の補助を受けた．本稿をまとめるにあたり，東京大学大学院工学系研究科の井上亮先生と筑波大学大学院システム情報工学研究科の堤盛人先生から貴重なご 示唆をいただいた．ここに謝意を表し ます．
参考文献

1] Braxmeier, H., Schmidt, V. and Spodarev, E.："Kriged road-trac maps"，In: J. Pilz (ed.) Interfacing Geostatistics, GIS and Spatial Data Bases,
Proceedings of the International Conference, StatGIS 2003, Springer-Series: Advances in SpatialScience(2005), 39{52.
2] Cressie, N.： Statistics for Spatial Data. Revised
Edition, John Wiley & Sons, 29{170, 1993.
3] Du, J. and Aultman-Hall, L.："Using SpatialAnalysis to Estimate Link Travel Times on Local Roads"，Transportation Research Board 2006
Annual Meeting CD-ROM，06-0676(2006).
4] Wackernagel, H 原著：『 地球統計学』，森北出版，
2003.
5] 井上 亮, 木越尚之, 清水英範：「 時空間クリギングの地価推定への適用可能性の検討」, 地理情報シス
テム学会講演論文集, 14(2005), 39-42．
6] 国土交通省：国土数値情報ダウンロード サービス，
http://nlftp.mlit.go.jp/ksj
7] 杉浦 淳徳，島川 陽一：「 Kriging 法を用いた土地価格の推定と空間分析」，『 高等専門学校の教育と研究: 日本高専学会誌』，Vol.10, No.4(2005)，42{49．
8] 瀬谷創，堤盛人，井上亮，石田東生，岡本直久：「 covariogram に 移動 平均モデ ル を用 いた 時 空 間Kriging 」，応用地域学会，12 月 9 日〜 10 日，2006，広島.
9] 丹後俊郎:『統計モデル入門』，朝倉書店, 2000．
10] ナビ タイムジャパン：ナビ タイム車ル ート 検索，
http://www.navitime.co.jp/drive/11] 間瀬 茂，武田 純：『空間データモデ リング 』，共立出版，135{151，1999．
6