測定誤差を考慮した隠れマルコフ劣化モデル
小林潔司1・貝戸清之2・林秀和3
1フェロー会員

2正会員

京都大学教授 経営管理大学院 (〒 606-8501 京都市左京区吉田本町)
E-mail:kkoba@psa.mbox.media.kyoto-u.ac.jp
大阪大学特任講師 大学院工学研究科グローバル若手研究者フロンティア研究拠点 (〒 565-0871 吹田市山田丘 2-1)
E-mail: kaito@ga.eng.osaka-u.ac.jp
3学生会員 京都大学大学院工学研究科都市社会工学専攻 (〒 615-8540 京都市西京区京都大学桂)E-mail:hhide@psa.mbox.media.kyoto-u.ac.jp

土木施設の劣化予測モデルとして，劣化過程を健全度間の推移確率を用いて表現したマルコフ劣化モデルが提案されている．土木施設の健全度の評価結果には，ランダム誤差や，多くの劣化事象の中から評価対象となる事象を選択する場合に生じるシステム誤差が存在する．本研究では，土木施設の劣化過程をマルコフ劣化モデル
として表現するとともに，健全度の測定結果に誤差が発生するメカニズムを隠れマルコフ劣化モデルを用いて表現する．その上で，真の健全度により定義されるマルコフ劣化モデルを，マルコフ連鎖モンテカルロ (MCMC)
シミュレーションによりベイズ推計する方法論を提案する．さらに，道路舗装を対象とした適用事例を通じて，
隠れマルコフ劣化モデルの有効性について考察する．

Key Words ：hidden Markov chain, measurement errors, representation matter, MCMC

1．はじめに

結果に誤差が存在する場合，測定結果によって得られる情報は「見かけの健全度」であり，「真の健全度」に関土木施設の統計的劣化予測モデルとして，マルコフする情報は獲得できない．見かけの健全度を用いて推劣化モデルが提案されている．マルコフ劣化モデルは，土木施設の劣化状態を複数の健全度を用いて評価し，土計したマルコフ推移確率が，真の推移確率に一致する保証はない．本研究では，現実に測定された「見かけの木施設の劣化による健全度間の時間的推移過程をマル健全度」と測定されない「真の健全度」の関係を混合分

コフ推移確率を用いて表現する．マルコフ劣化モデル

布モデルで表現するとともに，「真の健全度」を用いては操作性が高く，米国の標準的な橋梁マネジメントシス

定義されるマルコフ推移確率と，システム誤差を介在

テムの 1 つである PONTIS 等をはじめとして，アセッして測定される「見かけの健全度」の関係を隠れマル

トマネジメントにおける劣化予測モデルとして幅広く

コフ劣化モデルを用いて表現する．のちに，5．で考察適用されている1) ．

するように，隠れマルコフ劣化モデルは，通常の最尤マルコフ劣化モデルの適用事例では，複数の時刻で法を用いて推計することが困難である．そこで，本研実施された劣化状態の測定結果に基づいて，マルコフ

究では，隠れマルコフ劣化モデルをマルコフ連鎖モン

推移確率が推計されることになる．しかし，劣化状態のテカルロ (Markov Chain Monte Calro:以下，MCMC

測定結果には，誤差が含まれる場合が少なくない．測定と略す) シミュレーションを用いてベイズ推計する方法

誤差が発生する原因は多様であるが，本研究では，「対について考察する．

象とする道路区間に存在する劣化事象の中から，当該施設の劣化度をどの劣化事象で代表させるか」という以上の問題意識の下に，本研究では測定結果にシス代表値問題により発生するシステム誤差の問題に着目

テム誤差が介在するような土木施設の劣化過程を，隠する．特に，舗装のような土木施設では，ひび割れ，われマルコフ劣化モデルを用いて推計する方法論を提案

だち，非平坦性等といった劣化事象が膨大な数に及ぶする．さらに，道路舗装を対象として，本研究で提案

ため，路面性状調査において，常にもっとも損傷度がした方法論の有効性を検証する．以下，2．では，本研進んだ劣化事象が特定されるわけではない．このよう究の考え方を説明する．3．では，システム誤差が存在な代表値問題による測定誤差が発生する場合，路面性

しない劣化過程をマルコフ劣化モデルを用いて記述す

状調査における測定結果に基づいたマルコフ劣化モデ

る．4．では，システム誤差が存在するような劣化過程

ルの推計結果にシステム的なバイアスが発生する．

を隠れマルコフ劣化モデルを用いて定式化する．5．でマルコフ劣化モデルは，「真の健全度」間のマルコフ

は，隠れマルコフ劣化モデルを推計する方法論を提案

推移確率を用いて定義される．しかし，劣化状態の測定

する．6．では，道路舗装を対象とした実証分析につい1

て考察する．
50

6年間隔

45
事後測定わだち掘れ量（mm）

2．本研究の基本的な考え方(1) 従来の研究概要
土木施設の劣化予測に関しては，数多くの研究事例がある．土木施設の統計的劣化予測モデルとしてマル
コフ劣化モデルが提案されている．マルコフ推移確率

40
35
30
25
20
15
10
5

の推計方法として，1) 集計的推計方法と，2) 非集計的

0

推計方法が存在する．前者の方法は，ある一定の測定

0

期間の中で生起した健全度間の推移状態に関するデー
タに基づいて，マルコフ推移確率を直接推計すること

10
20
30
40
事前測定わだち掘れ量（mm）

50

推移状態に関する実データの数え上げにより，推移確注) 45 度線より下方に位置するサンプルは，時間の経過により健全度が改善された事例を表している．このようなサンプルは，理論的には存在しえない．

率を直接定義する方法2) である．これに対して，最尤法

図-1 道路舗装のシステム誤差を目的とする．もっとも単純な算定方法は，健全度間のにより，推移確率を推計する方法3) も提案されている．
マルコフ推移確率は，推移確率を定義する測定間隔に依存する．現実に測定される健全度データには，測定る階層型指数劣化ハザードモデル10) が提案されている．

間隔が異なる多様なデータが混在している場合が多い．

また，マルコフ推移確率の推計方法に関しては，測定

この場合，実データが測定された測定間隔の差異がもデータが非常に少ない段階で，技術者の経験情報と測たらす影響を補正することが必要となる．杉崎等4) は，定結果を結合してマルコフ推移確率を推計するベイズ

異なる測定間隔を有する目視点検データを用いて，マ推計モデル11) ; 12) ，予防補修により測定データが欠損す

ルコフ推移確率を集計的に推計する方法を提案していることにより発生する欠損バイアスを補正する方法13) ，る．しかし，このような集計的劣化予測方法では，個々ハザード率の異質性を考慮したランダム比例ワイブルの施設が置かれている使用環境や，施設が有する構造

劣化ハザードモデル14) が提案されている．以上の方法

的，機能的特性と推移確率との関係をモデル化できなは，いずれも測定結果に誤差が存在しない場合を想定いという限界がある．

したものである．しかし，本研究の適用事例としてとこれに対して，非集計的推計方法は，個々の土木施りあげる道路舗装をはじめとして，土木施設の健全度

設の劣化過程に関する情報に基づいて，その背後にあの測定結果に多大な誤差が含まれる場合が少なくない．

る劣化過程の統計的規則性を推計する方法である．こ

本研究では，システム誤差を考慮した劣化予測モデルのような非集計的推計方法として，貝戸等5) は，ニューを隠れマルコフ劣化モデルを用いてモデル化することヨーク市における橋梁の目視検査データを用いて，橋を目的とする．筆者等の知る限り，システム誤差を考

梁の劣化速度に着目した平均劣化曲線の算出方法を検慮してマルコフ推移確率の非集計的推計方法を開発し討している．また，劣化速度を確率変数と捉えて，過去た研究事例は他に見あたらない．

の検査履歴を反映したマルコフ推移確率の推計方法を(2) システム誤差と代表値問題提案した．その後，非集計的推計方法は，ハザードモデ
ルの適用により，飛躍的な発展を遂げている．その中本研究では，隠れマルコフ劣化モデルの適用事例とで，青木等は，ワイブルハザードモデルを用いて，トして舗装の劣化過程をとりあげる．図-1 は，実証分析

ンネル照明の寿命解析を行っている6) ．Mishalani andの対象としてとりあげた舗装の劣化過程に関するデー

Madanat7) は，2 つの隣接する健全度のみを対象として，
マルコフ推移確率を指数ハザードモデルを用いて表現タを示したものである．図中の各点は，同一の道路区間する方法を提案した．これとは独立に，津田等8) は，2測定されたわだち掘れ量の対応関係を示している．こつ以上の任意の健全度間における推移状態を表現する

こで，各サンプルにより，道路舗装の特性や，測定間多段階指数ハザードモデルを提案し，マルコフ推移確隔は多様に異なっていることを断っておく．ただし，同率を推計する一般的な方法論を提案した．その後，マ

図では，2 つの測定時刻の間に，舗装の予防的補修が実

ルコフ推移確率が過去の記憶を有する非斉次マルコフ

施されなかった道路区間に関するサンプルのみを取り推移確率を推計するための多段階ワイブル劣化ハザー上げている．道路舗装の予防的補修を実施しない限り，

ドモデル9) ，異なる劣化パターン間の推移過程を表現す

時間の進行に伴って道路舗装の健全度が，改善されるこに対して，2 つの測定時刻úA とúB (úA < úB ) において2

とはあり得ない．したがって，測定サンプルは，すべて図中の 45 度線より下方に位置しなければならない．しかし，同図に示すように，本来あり得ない 45 度線より

見かけの健全度

真の健全度

m(ú)

mÉ(ú)

上方の領域に位置するサンプルも相当程度存在し，舗n

装の健全度データには，かなりの程度の測定誤差が存m

在していることが推察できる．

i

このような測定誤差として，1) ランダム誤差，2) 代

j表値問題によるシステム誤差が考えられる．道路舗装の劣化状態を，たとえば路面性状測定車で計測する場
Z

合，必然的に機器測定に伴う測定誤差が介在してくる．
これらの測定誤差は，基本的にはランダム誤差と考えúAることができる．もちろん，測定機器特有の誤差や測定

Z

úB

ú

úA

úB

ú

対して，平均値等マクロ的な特性を比較検討すること

注) システム的誤差により，時間の経過により見かけの健全度が改善された事例を表している．システム的誤差が存在するため，見かけの健全度は，常に真の健全度より上方に位置している．

により，システム誤差が存在するかどうかを検討する

図-2 健全度とシステム誤差

日特有の誤差が介在する場合，システム誤差が発生す
る．異なる測定業務により測定された測定サンプルにことが必要であることは言うまでもない．しかし，舗装の劣化状態の測定結果には，別のシステム的な理由により発生する測定誤差が存在する．道路舗装の劣化も損傷が進んだ箇所の健全度を求めるとともに，これにより現れるわだち，ひび割れ，非平坦性等の個別のらの健全度の平均値を区間全体の健全度として定義す現象を劣化事象と呼ぼう．ある道路区間における舗装ることもある．このような方法で健全度を定義しても，の健全度は，その道路区間の中に存在する劣化事象に時間の経過とともに，健全度が回復するという事象は基づいて定義される．ところが，道路区間の中には，膨理論的には起こり得ない．しかし，路面性状調査によ

大な数の劣化事象が存在しており，路面性状調査により獲得されたデータは，劣化事象の中からサンプル抽

りすべての劣化事象に関する情報が得られるわけでは出されたものであり，常に損傷がもっとも進んだ劣化ない．その結果，道路舗装の健全度は，道路区間に存事象が抽出されるわけではない．したがって，ランダ在する劣化事象の中からランダムに抽出されたサンプム抽出された劣化事象で定義される健全度は，真の健ルに基づいて定義されることになる．この場合，今回全度より損傷の程度が少なくなるというシステム的なの測定時刻において，前回の測定時に健全度を定義しバイアスが発生することになる．本研究では，現実のたサンプルより，損傷の程度が小さい劣化事象がサン舗装の測定データに基づいて，システム誤差を考慮しプルされる (健全度が前回よりも改善する) 場合も起こた隠れマルコフ劣化モデルを作成するとともに，舗装

りえる．しかし，補修が実施さらない限り，前回にサの測定結果にシステム的な測定バイアスが存在するかンプルされた劣化事象は，少なくとも損傷が改善されどうかを実証的に分析する．なお，舗装劣化予測の実ずに残存しているわけであり，健全度が改善されるは務では，図-1 の 45 度線以上に位置するサンプルをデー

ずはない．この場合，当該道路区間の健全度を定義す

タベースから除去し，45 度線より下方に位置するサンる際に，「どの劣化事象が代表としてとりあげられていプルのみを用いて劣化曲線を推計する場合が多い．しるのか」ということが問題となる．このように，土木

かし，この場合でも，代表値問題によるシステム的な施設の健全度が異なる劣化事象に基づいて定義される

推計バイアスの問題は解決していない．

ことによるシステム誤差を代表値問題 (representation(3) システム誤差と隠れマルコフ劣化過程

matter) と呼ぶこととしよう．
土木施設の劣化状態を評価する場合，対象とする施健全度の測定にシステム誤差が存在する場合に生起す

設に存在する劣化事象の中で，もっとも損傷度が進行しる問題を図-2 を用いて説明しよう．土木施設の真の劣化

た劣化事象に基づいて健全度が定義されることが多い．

状態が I 個のレーティング指標 (本研究では，健全度と呼例えば，実証分析の対象事例では，同一道路区間に対ぶ) i (i = 1; ÅÅÅ; I) で表現されると考えよう．i = 1 の場

して 3 箇所の舗装状態を測定し，その中からもっとも損合，健全度がもっとも優れた状態にあり，i = I はもっと傷度が進んだ箇所の健全度を用いて，当該道路区間のも悪い状態と対応している．いま，ある時刻úA で健全度

健全度を代表させている．あるいは，道路区間をさらを測定した結果，健全度が m(úA ) = m (m = 1; ÅÅÅ; I)に複数区間に細分化し，細分化された各区間のもっとであると判定されたとしよう．しかし，測定の結果に3

は，測定誤差が含まれる可能性がある．時刻úA における
真の (システム誤差がない) 健全度を mÉ(úA ) = i (i =健全度1; ÅÅÅ; I) と表そう．健全度調査の結果にシステム誤差が存在する場合，測定によって得られた「見かけの健全度」m(úA ) = m が，「真の健全度」mÉ(úA ) = i に一致

12

する保証はなく，真の健全度より小さい値が観測され3る可能性がある．ここで，健全調査による測定結果が，ô

離散的確率分布 fi (mjãi ) に従って分布すると考えよう．

i

ただし，確率分布 fi (mjãi ) は，「真の健全度」が i であ

i +1

るときに，健全度が m と判定される条件付き確率分布(以下，条件付システム誤差分布と呼ぶ) である．また，

ú0

ú1

úA

ú2

úB

úi

ú

ãi は，確率分布を特徴づけるパラメータベクトルである．つぎに，時刻úA から，ある一定の時間 z が経過した

図-3 劣化過程

時刻úB = úA + z に，再び測定を実施したとしよう．健
全調査による測定結果を m(úB ) = n と，時刻úB における「真の健全度」を mÉ(úB ) = j と表そう．この 2 つの3．マルコフ劣化モデル

健全度も互いに一致する保証はない．以上の測定結果(1) モデル化の前提条件

より，期間 [úA ; úB ) の間に生起する「見かけの健全度」の変化パターンは m ! n となる．一方，「真の健全度」

土木施設の劣化過程をモデル化するために，ひとまの推移パターンは i ! j である．システム誤差が存在すずシステム誤差が存在せず，常に真の健全度が測定さる場合，それぞれの時刻で測定された健全度が「真のれる場合を考えよう．システム誤差を考慮した劣化過健全度」であるかどうかは判らない．したがって，健全程に関しては，4．で議論する．土木施設の劣化予測モ

度調査によって測定されるのは「見かけの推移パター

デルを推計するためには，施設の劣化状態に関する時ン」m ! n であり，「真の推移パターン」i ! j は測定

系列データを蓄積することが必要となる．いま，ある土

されないことになる．
「見かけの推移パターン」に関す

木施設の劣化に関する履歴が図-3 に示すように与えらるデータに基づいて，健全度間の推移確率ômn を推計れたとする．同図は，施設が補修されずに放置された

しても，それが「真の推移パターン」に基づいて定義

時に，劣化がどのように進展するかを表したものであされる推移確率ôij に一致する保証はない．このように，る．現実には，施設の劣化過程には不確実性が含まれ，測定結果に誤差が介在する場合，「真の健全度」間で定しかも劣化状態は時間軸上の限られた時刻で実施され義されるマルコフ劣化過程が，システム誤差を含んだる測定を通じてのみ知ることができる．図中，時刻úは「見かけの健全度」の間で定義されるマルコフ連鎖の背カレンダー上の実時刻 (以下，時刻と呼ぶ) を表す．時

後に隠れてしまうという現象が生じる．

刻ú0 で土木施設の使用が開始された直後から劣化が始
まる．2． (2) で議論したように，施設の劣化状態は I

本研究では，このようにシステム誤差を含んだ見か個の健全度で記述される．施設の健全度を表す状態変けのデータの中に，真の劣化過程を表すマルコフ連鎖数を i (i = 1; ÅÅÅ; I) で表現しよう．施設がもっとも健が隠れている状況を隠れマルコフ劣化モデルを用いて全な (劣化が進展していない) 状態を i = 1 で表し，状モデル化する点に特徴がある．隠れマルコフ劣化モデル

態変数 i の値が大きくなるほど，劣化が進展しているこでは，測定結果に誤差が含まれるメカニズムを表す確率とを表す．i = I の場合，当該の施設が使用限界に到達

fi (mjãi ) と，「真の健全度」に対して定義される推移確率ôij を同時推計することが課題となる．なお，正確には，「真の健全度」i が測定できないため，「見かけの健全していることを表す．図中の時刻úi (i = 1; ÅÅÅ; I Ä 1)

において，健全度が i から i + 1 に移行している．しか度」m が，異なる条件付システム誤差分布 fi (mjãi ) (i =される．したがって，測定者が獲得できる情報は，測1; ÅÅÅ; I) の中から，どの条件付システム誤差分布により発生したデータなのかも判らない．したがって，「見
かけの健全度」に関するデータを発生した条件付シス
テム誤差分布を推計する方法論が必要となる．このため，隠れマルコフ劣化モデルの推計問題は複雑となるが，推計方法に関しては 5．で考察する．

定時刻における健全度のみであり，健全度が推移した

し，測定は同図における 2 つの時刻úA ; úB において実施

時刻úi (i = 1; ÅÅÅ; I Ä 1) に関する情報は得られない．

(2) マルコフ劣化モデル
土木施設の劣化過程をマルコフ劣化モデルを用いて表現しよう．いま，2 つの時刻間における真の健全度間

4の推移状態をマルコフ推移確率で表現する．時刻úA にと表せる．すなわち，ハザード関数ïi (yi ) は，初期時おける真の健全度を状態変数 mÉ(úA ) を用いて表そう．

刻から時間 yi が経過するまで真の健全度 i の状態が継続

時刻úA における真の健全度が i (i = 1; ÅÅÅ; I) であればしたという条件の下で，期間 [yi ; yi + Åyi ) 中に真の健m (úA ) = i と表せる．マルコフ推移確率は，時刻úA で

全度 i + 1 に進展する条件付確率である．ハザード関数

測定された真の健全度 mÉ(úA ) = i を与件とし，将来時がサンプル時間軸上の時刻 yi に依存せず，常に一定値

É

É

刻 (たとえばúB ) において真の健全度 m (úB ) = j が生

íi > 0(i = 1; ÅÅÅ; I Ä 1) をとる場合，指数ハザード関数

起する条件付推移確率として定義される．すなわち，
É

É

Prob[m (úB ) = jjm (úA ) = i] = ôijïi (yi ) = íi

(5)

(1)

が成立する．指数ハザード関数を用いることにより，劣と表せる．このような推移確率を真の健全度ペア (i; j)

化過程が過去の履歴に依存しないというマルコフ性をに対して求めれば，マルコフ推移確率行列10
ô11 ÅÅÅ ô1I
B .
.. C
..
C
.
Ö =B
.
. A
@ .
0 ÅÅÅ ôII

表現できる．さらに，指数ハザード関数を用いれば，真の健全度 i の寿命が yi 以上となる確率F~i (yi ) は，

F~i (yi ) = exp(Äíi yi )(2)と表現できる．

さらに，サンプル時間軸上のúA で，真の健全度が i でを定義できる．マルコフ推移確率 (1) は所与の 2 つの時あり，かつ時刻úA から追加的に z (ï 0) 以上にわたって

刻úA ，úB の間において生じる真の健全度間の推移確率

真の健全度 i が継続する確率F~i (úA + zjêi ï úA ) は，を示したものであり，当然のことながら，対象とする

F~i (úA + zjêi ï úA ) = Probfêi ï úA + zjêi ï úA gexpfÄíi (úA + z)g(7)== exp(Äíi z)exp(Äíi úA )測定間隔が異なれば推移確率の値は異なる．補修がな
い限り常に劣化が進行するので，ôij = 0 (i > j) が成PI立する．また，推移確率の定義より j=i ôij = 1 が成立する．すなわち，マルコフ推移確率に関して9ôij ï 0 (i; j = 1; ÅÅÅ; I) >=ôij = 0 (i > j の時)>PI;
j=i ôij = 1(6)と表される．すなわち，測定時刻úA において真の健全度が i と判定され，次の測定時刻úB = úA + z においても真の健全度が i と判定される確率は，(3)Prob[mÉ(úB ) = ijmÉ(úA ) = i] = exp(Äíi z) (8)が成立しなければならない．状態 I は，補修のない限りとなる．ただし，z は 2 つの測定時刻の間隔を表す．確マルコフ連鎖における吸収状態であり，ôII = 1 が成立

率 Prob[mÉ(úB ) = ijmÉ(úA ) = i] はマルコフ推移確率

すると考える．なお，マルコフ推移確率は過去の劣化ôii (z) にほかならない．指数ハザードを用いた場合，推移確率ôii (z) はハザード関数íi と測定間隔 z のみに依存
し，時刻úA ，úB に関する情報を用いなくとも推移確率を推計することが可能となる．以上の議論を拡張し，指数ハザード関数を用いて，測定時刻úA とúB = úA + z の履歴とは独立して定義される．マルコフ推移確率モデ
ルでは，真の健全度が i Ä 1 から i に推移した時刻に関わらず，測定時刻úA から測定時刻úB の間に推移する確率は時刻úA における真の健全度のみに依存するという
性質 (マルコフ性) を満足する．

間で真の健全度が i から j (> i) に推移するマルコフ推
移確率ôij (z) (i = 1; ÅÅÅ; I Ä 1; j = i; ÅÅÅ; I) は，(3) 多段階指数ハザードモデル

ôij (z) = Prob[mÉ(úB ) = jjmÉ(úA ) = i]

マルコフ推移確率は，多段階指数ハザードモデルを用いて推計できる．本研究では，津田等8) が開発した多=k=i

段階指数ハザードモデルを用いるが，読者の便宜を図Yímímexp(Äík z)
í
Ä
í
í
k
m+1 Ä ík
m=i m
jÄ 1

m=k

(i = 1; ÅÅÅ; I Ä 1; j = i + 1; ÅÅÅ; I)

るためモデルの概要を説明しておく．いま，真の健全度i (i = 1; ÅÅÅ; I Ä 1) の寿命を確率変数êi で表す．真の健
全度 i の寿命が，確率密度関数 fi (êi )，分布関数 Fi (êi )に従うと仮定する．時刻úA における真の健全度 i であり，
そこから時間 yi が経過した時刻で真の健全度 i + 1 に到kÄ1Y現する．この時，ハザード関数は，供用時間 yi まで真の健全度が i のまま継続する生存確率F~i (yi ) を用いて，

kÄ1
Y

=
5

ím
exp(Äík z)ím Ä íkYímím
exp(Äík z)í Ä íkím+1 Ä íkm=i m

m=i;6=k(4)(9)と表すことができる8) ．ただし，表記上の規則として，( QkÄ 1ím(k = i の時)ím Äík = 1
Qm=i
jÄ 1
ím
m=k ím+1 Ä ík = 1 (k = j の時)が成立すると考える．さらに，表記の便宜上，達する確率密度をハザード関数15) ; 16) ïi (yi ) を用いて表fi (yi )Åyiïi (yi )Åyi =F~i (yi )

j kÄ
X
Y1

jÄ 1

m=k

と簡略化する．また，ôiI に関しては，マルコフ推移確

ちに 4． (3) で言及するような困難性が存在するため，率の条件より次式で表せる．
IÄ1X
ôiI (z) = 1 Äôij (z) (i = 1; ÅÅÅ; I Ä 1) (10)

十分な研究が蓄積されていない．そのため，隠れマルコ(4) 期待劣化曲線

研究が進展した．隠れマルコフ連鎖モデルは，動的混

土木施設の劣化特性は，施設の環境条件や構造・機合分布モデルの中で，レジーム遷移 (regime switching)

能的な特性に依存して変化する．多段階指数ハザード

モデルと類似の確率構造を有している．このため，レモデルを用いれば，土木施設の特性別にマルコフ推移

ジーム推移モデルの推計方法を，隠れマルコフ連鎖モ

確率を求めることができる．いま，多段階指数ハザー

デルの推計にも適用できることが判明した．

フ連鎖モデルの適用事例もそれほど多くない．近年になり，時系列モデルに関する研究の発展とともに，動的j=i

混合分布モデルと呼ばれる時系列統計モデルに関するレジーム遷移モデルは，Hamilton18) によって提案さドモデルのハザード率íi (i = 1; ÅÅÅ; I Ä 1) を施設特性

x = (x1 ; ÅÅÅ; xM ) を用いて次式により表そう．
íi = íi (x) = xå0i

た だ し ，åi=れ，景気分析や金融計量経済学の分野で応用研究が急速に進展しつつある．レジーム遷移モデルでは，時系(11)

列データの構造変化を，レジーム間の推移現象として(åi1 ; ÅÅÅ; åiM ) は 未 知 パ ラ メ ー タ

把握し，レジーム間推移確率を非線形時系列モデルと

åim (m = 1; ÅÅÅ; M ) の行ベクトルである．記号 0 は転置操作を表す．劣化推移確率はデータが観察されたして定式化する19) ．たとえば，Hamilton18) は，レジーム間の推移確率過程を Markov 過程で表現したような

測定間隔 z に依存する．したがって，推移確率ôij を目

マルコフ遷移モデルを提案している20) ．マルコフ遷移

視測定による実測データ (z; x) と未知パラメータå =モデルは，時系列データの推移過程をマルコフ連鎖を(å1 ; ÅÅÅ; åIÄ 1 ) の関数としてôij (z; x : å) と表すことができる．さらに，ハザード率íi (x) を用いた生存関数を
F~i (yi jíi (x)) と表記する．生存関数F~i (yi jíi (x)) が式 (6)用いて表現するモデルである．マルコフ遷移モデルの推計方法の発展により，隠れマルコフ連鎖モデルの推計上の困難性も克服されることが判明した．とりわけ，で表されることに留意すれば，特性xを有する土木施設

MCMC (Markov Chain Monte Carlo: マルコフ連鎖モンテカルロ) 法の発展が，マルコフ遷移モデルの推計に大きく貢献している．すでに，MCMC 法を用いた隠れマルコフ連鎖モデルの推計方法に関しても，いくつかの研究事例がある21) ．しかし，これらの既往研究では，が，健全度 i にはじめて到達した時刻から，劣化が進展して次の健全度に進むまでの期待期間長 (以下，健全度

期待寿命と呼ぶ) RM Di (x) (i = 1; ÅÅÅ; I Ä 1) は，
Z 1
F~i (yi jíi (x))dyi
RM Di (x) =
0
Z 1
1
=
expfÄíi (x)yi gdyi =
(12)
xå0i
0
と表される．土木施設の補修直後の健全度を i = 1 としマルコフ推移確率を直接求めることを目的としており，

2． (1) で言及した集計的手法に属する．本研究では，
非集計的方法でマルコフ推移確率を求めることを目的

よう．この時，初期時刻の健全度 i = 1 の状態から劣化としており，このような立場から隠れマルコフ連鎖モが進行し，健全度が j (> 1) に推移するまでに要する平デルを推計する研究事例は他に見あたらない．

均的経過時間 ETj (x) (j = 2; ÅÅÅ; I) は

ETj (x) =

j
X
1
xå0i
i=1

(2) 混合分布モデル(13)

土木施設の劣化過程は不確実であり，マルコフ劣化モ
デル (2) に従うと考えよう．時間軸上の時刻úA においてと表される．ETj (x) (j = 2; ÅÅÅ; I) は，直近の補修時

刻から健全度 j に到達するまでの平均的な経過時間を表測定が実施され，対象施設の見かけの健全度 m(úA ) =す．本研究では，健全度 j (j = 1; ÅÅÅ; I) と平均的経過

m (m = 1; ÅÅÅ; I) が測定されたとする．ただし，測定結果に誤差が含まれるため，測定された見かけの健全度 m が，真の健全度 mÉ(úA ) = i (i = 1; ÅÅÅ; I) に一

時間 ETj (x) の関係を期待劣化曲線と呼ぶ．

致しているかどうかは判らない．いま，時刻úA におけ

4．隠れマルコフ劣化モデルの定式化

る真の健全度が mÉ(úA ) = i であると仮定しよう．こ(1) 隠れマルコフ劣化モデルの時，見かけの健全度 m(úA ) = m が測定される確率

隠れマルコフ連鎖モデルは，測定される状態変数に分布を fi (mjãi ) と表そう．ãi は，確率分布 fi (mjãi ) を

誤差が介在するようなマルコフ連鎖モデルである．隠れ

特徴づけるパラメータベクトルである．以下，この確マルコフ連鎖に関する理論的研究に関しては研究の蓄率分布を条件付システム誤差分布と呼ぶ．条件付シス積があり，いくつかの成書17) で詳細に紹介されている．

テム誤差分布は，真の健全度ごとに定義される．条件

しかし，隠れマルコフ連鎖の推計方法に関しては，の付システム誤差分布 fi (mjãi ) は，時間に依存しないと6仮定する．ここで，測定者にとって，真の健全度 i が不度の生起確率は，「時刻úA における見かけの健全度がど確実であることに留意しよう．すなわち，測定されたの確率密度から生成されたか」に依存しているため，時見かけの健全度が，いずれの条件付システム誤差分布

間を通じて測定される見かけの健全度の時系列データ

fi (mjãi ) (i = 1; ÅÅÅ; I) から発生されたサンプルなのかは互いに独立にはならない．

が測定される尤度 `(m(úA ) = m) を

(3) 初期値問題が不確実である．ここで，見かけの健全度 m(úA ) = m`(m(úA ) = m) =IXôi (úA )fi (mjãi )尤度 (17) には，真の健全度に関する初期分布ôi (úA )，(14)

マル コ フ推 移 確率ôij (z)，条 件付 シ ス テム 誤 差分 布

i=1

fi (mjãi ) という 3 種類の未知確率が存在する．このうち，初期分布ôi (úA ) に関する先験的情報が存在しないと表そう．ôi (úA ) は見かけの健全度 m(úA ) = m がfi (mjãi ) から生成される確率である．この確率は，時刻úA で真の健全度が m である確率に等しい．言い換えという問題 (初期値問題と呼ぶ) がある．初期値問題をれば，互いに排反な I 個の原因のうち，原因 i によって克服する 1 つの方法は，初期分布としてノンパラメト

事象 m(úA ) = m が発生する確率がôi (úA ) であるとみ

リック分布を与える方法がある．しかし，多様なサンプなせる．式 (14) は，複数の条件付システム誤差分布を

ルごとに初期分布が異なるため，膨大な数のパラメー

加重平均した確率分布を表しており，混合分布モデル

タが必要となる．本研究では，更新，補修時刻で，健全度が i = 1 に回復すると考えよう．すなわち，更新・(mixture distribution model)22) と呼ばれる．
つぎに，2 つの測定時刻úA ; úB (úA < úB ) に着目し，
時刻úA に見かけの健全度 m(úA ) = m が測定され，かつô(ú0 ) = fô1 (ú0 ); ÅÅÅ; ôI (ú0 )g

時刻úB = úA + z に，見かけの健全度 m(úB ) = n が同

= (1; 0; ÅÅÅ; 0)時に測定される尤度 `(m(úA ) = m; m(úB ) = n) を再帰

時刻における初期値分布を設定することにより，初期ことが判明しているときに，時刻úB で見かけの健全度値問題を克服することが可能になる．

m(úB ) = n が測定される条件付尤度 `i (m(úB ) = n) は`i (m(úB ) = n) =ôij (z)fj (njãj )いま，直近の補修・更新時刻ú0 以降，健全度が測定された時刻をú1 ; ÅÅÅ; úT と表そう．T は直近の更新・補修以(15)降の測定回数を表す．時刻út に測定された見かけの健全j=iと表される．しかし，時刻úA における真の健全度も観測度を m(út ) = mt (t = 1; ÅÅÅ; T ) と表そう．さらに，第できない．いま，時刻úA で真の健全度が i である確率をt Ä 1 回目の測定と t 回目の測定間隔を zt (t = 1; ÅÅÅ; T )と表そう．この時，測定間隔ベクトルz = (z1 ; ÅÅÅ; zT )を与件とし，健全度m = (m1 ; ÅÅÅ; mT ) が測定される尤

ôi (úA ) とすれば，時刻úB に，見かけの健全度 m(úB ) = nが測定される尤度 `(m(úB ) = n) は`(m(úB ) = n) ==IXôi (úA )

i=i

I
X

I
X

度関数 L(ã; m; z) を再帰的に定義すれば

ôi (úA )`i (m(úB ) = n)i=iôij (z)fj (njãj )L(ã; m; z) =

(16)

j=i

`h (t) =と表される．さらに，時刻úA において，真の健全度が i

I
X

j=h

の時に，見かけの健全度 m(úA ) = m が測定される確率

I
X
j=1

ô1j (z1 )fj (m1 jãj )`j (1) (19a)ôhj (zt )fj (mt jãj )`j (t + 1)

(19b)

(1 î t î T Ä 1)

が fi (mjãi ) で表されることより，時刻úA に見かけの健`h (T ) =全度 m(úA ) = m が測定され，かつ時刻úB に，見かけの

I
X

j=h

健全度 m(úB ) = n が同時に測定される尤度 `(m(úA ) =ôhj (zT )fj (mT jãj )(19c)と表される．津田ら8) は，多段階指数ハザードモデル

m; m(úB ) = n) は`(m(úA ) = m; m(úB ) = n)
0
1
I
I
X
X
ôi (úA )fi (mjãi ) @ôij (z)fj (njãj )A (17)
=
i=1

(18)

が成立すると考える．このように，更新・補修直後の的に定義しよう．まず，時刻úA で真の健全度が i である
I
X

補修直後の時刻ú0 における初期値分布に関して(9) を最尤法を用いて推計する方法を提案している．しかし，隠れマルコフ劣化モデルの尤度関数 (19a)-(19c)は最尤法に適さない性質を持っていることが知られている23) ．特に，尤度関数がôij (z) に関して高次の非線j=iと表される．式 (17) に示すように，見かけの健全度形多項式となっており，1 階の最適化条件が (複素数解m(úA ) = m と m(úB ) = n が抽出される確率 fi (mjãi ),fj (njãj ) は，互いに，マルコフ推移確率ôij (z) を通じて相関がある．すなわち，時刻úB における見かけの健全を含めて) 非常の多くの解を有している点にある．さらに，最適化条件の (複数個の) 実数解が 0 と 1 の間にある保証がない点である．最尤法の代わりにベイズ推計

7法を用いれば，高次の非線形多項式を解く問題を回避fstÄ1 ;ÅÅÅ; st+1 g) の全条件付事後確率は，ベイズの法則よりできる．しかし，尤度関数 (19a)-(19c) が，極めて多くの項を含んでおり，計算量が膨大になってしまう欠点

Probfst = ijsÄ t ; ã; òg~ i ; ã; m; z)L(s= Pst+1 Ä t i~i=stÄ 1 L(sÄ t ; ã; m; z)がある24) Ä 27) ．このような最尤法の難点を克服するた
めに，尤度関数の完備化操作が必要となる．

!it fi (mt jãi )= Pst+1
j=stÄ 1 !jt fj (mt jãj )(4) 完備化操作
ある土木施設に対して，測定時刻út (t = 1; ÅÅÅ; T ) にと表される．ただし，
8>< ô1j ôjs2!jt =ôstÄ 1 j ôjst+1>:
ôsT Ä 1 j

おいて見かけの健全度m = (m1 ; ÅÅÅ; mT ) が測定できたと考えよう．さらに，隠れマルコフ劣化モデルを推計するために，見かけの健全度mが，どの真の健全度を持つ条件付システム誤差分布から生成されたかを示
す潜在変数ベクトルs = (s0 ; ÅÅÅ; sT ) を導入しよう．た

(24)

t=T

と 表 さ れ る ．す な わ ち ，真 の 健 全 度 間 の 推 移 確 率

だし，劣化過程の性質より，施設が補修されない限り，

s0 = 1 î s1 î ÅÅÅî sT î I

t=1
2îtîT

(23)

ôij (z) (i = 1; ÅÅÅ; I; j = i; ÅÅÅ; I) と条件付システム誤差分布 fi (mjãi ) (i = 1; ÅÅÅ; I) が求まれば，sÄ t を与件(20)とした時刻 t の真の健全度 st 2 fstÄ 1 ; ÅÅÅ; st+1 g の全条を満足する．真の健全度は本来測定不可能であり，潜件付事後確率を求めることができる．完備化尤度関数

在変数sは，本来測定不可能な変数である．議論の便宜(22) では，潜在変数sは確定的である．ただし，条件付

上，ひとまず潜在変数が仮に測定できたと考えよう．さシステム誤差分布, 真の健全度間の推移確率には未知パ

らに，潜在変数sの測定結果に基づいて，ダミー変数(
1 st = i(21)éti =0 st 6= i

ラメータã; åが含まれており，潜在変数に関する全条
件付事後確率を先験的に求めることができない．全条
件付事後確率 (23) を用いた MCMC 法を用いて，反復的(t = 1; ÅÅÅ; T ; i = 1; ÅÅÅ; I)に潜在変数sをランダム発生させ，パラメータã; åをベを導入しよう．見かけの健全度m，潜在変数s，測定間

イズ推計することになる．これらの未知パラメータと隔z を与件とした尤度関数 (19a)-(19c) は

全条件付き事後確率を求める方法についても，5． (4)で改めてとりあげる．このような手続きにより，完備~ ã; m; z)
L(s;
I n
Y
ô1i (z1 )é1i fi (m1 jãi )é1i
=
T Y
I
Y

化尤度関数を用いて求めたパラメータのベイズ推計値
が，真の尤度関数を用いて求めたパラメータの最尤推

i=1

T n
Y

t=2 j=i

=

T
Y

t=1

=

ôij (zt )étÄ 1i étj fj (mt jãj )étjt=1

ôstÄ 1 st (zt )fst (mt jãst )ôstÄ 1 st (zt )

T
Y

t=1

o

fst (mt jãst )o

計値に収束することが証明されている22) ．

(5) 条件付システム誤差分布
条件付システム誤差分布 fi (mjãi ) は，真の健全度がÉ

m (út ) = i の時に，見かけの健全度 m(út ) = m が測定される確率分布を表している．ただし，条件付システム
誤差分布の関数形は，測定時刻út (t = 1; ÅÅÅ; T ) に依存(22)

せず，時間を通じて一定であると仮定する．いま，対象と表現できる26) ．以上の操作を完備化 (completion) ととする施設の真の健全度が i (i = 1; ÅÅÅ; I) であるとしよう．しかし，測定機器により測定される健全度が，真言う．完備化された尤度関数 (以下，完備化尤度関数と呼ぶ) (22) は，通常の尤度関数 (19a)-(19c) より大幅に簡略の健全度に一致する保証はない．健全度の測定結果にシ化されていることが理解できる．ただし，尤度関数 (22)ステム誤差が介在し，見かけの健全度 m (m = 1; ÅÅÅ; i)が測定されることになる．

の中に含まれる潜在変数sは，測定できない変数である．
そこで，完備化尤度関数を用いて，潜在変数の確率分布

測定結果にシステム誤差が介在する場合，真の健全を推計することを考える．完備化尤度関数を展開すれば，度 i と見かけの健全度 m の間に m î i という関係が存

潜在変数sに関する全条件付事後分布 (full conditional

在する．例えば，本研究の実証分析では，道路舗装のposterior distribution) を導出できる．劣化過程の特性により，補修が実施されない限り，条件 (20) が成立する．ここで，sÄ t = (s1 ; ÅÅÅ; stÄ 1 ; st+1 ; ÅÅÅ; sT ); siÄt =劣化問題をとりあげるが，そこでは単位道路区間内の

3 つの代表地点における健全度を測定し，その中でもっとも損傷が大きい地点の健全度により，当該区間の健全度として定義している．しかし，当該路区間におけ(s1 ;ÅÅÅ;stÄ1 ; i; st+1 ; ÅÅÅ; sT ) とすれば，st = i (i 28るすべての地点における劣化事象を測定しているわけのことが，ベイズ統計学をハザード解析へ適用する際ではないため，ランダムサンプリングにより測定したに，大きな障害になっていた．しかし，近年，MCMC

劣化事象が，当該道路区間においてもっとも損傷が進法21) がベイズ統計学の分野に導入され，多重数値積分

んだ劣化事象である保証はない．ここで，対象とするにより基準化定数を求めなくても，効率的に事後分布

土木施設において，もっとも損傷が進んだ劣化状態でを求めることが可能となった．その結果，ベイズ推計法

定義される真の健全度 i に対して，見かけの健全度 m がの適用範囲は大幅に拡大したと考えることができる．す選択される確率を定義しよう．本研究では，このようでに，MCMC 法を用いたベイズ推計法に関して，いくなシステム誤差が右切断されたノンパラメトリックなつかの研究が蓄積されている21) ．代表的な MCMC 法

離散的確率分布 gi (mjãi ) (m = 1; ÅÅÅ; i) に従うと考えとして，ギブスサンプリング (Gibbs sampling) 法，メ

度が m (m î i) となる確率 gi (mjãi ) は(
0m > i の時
gi (mjãi ) =ãim m î i の時

MH と略す) 法等が提案されている21) ．この内，ギブス
サンプリングはもともと画像復元のアルゴリズム31) とる．すなわち，真の健全度が i の場合に，見かけの健全

トロポリス・ヘイスティングス (Metropolis-Hastings：(25)

して知られていたが，ベイズ推計法における事後分布の推計に応用された32) ．ギブスサンプリング法，MHiと表される．ただし，ãmは定数であり，

法は，いずれも事後確率密度関数を直接求めることが0 î ãim î 1

(26a)

i
ãm
=1

(26b)

i
X

難しい場合に，各パラメータの条件付き事後確率密度関数を用いて，反復的にパラメータåのサンプルを乱数発生させることにより，事後分布からパラメータサン

m=1を満足する．誤差発生メカニズム (25) にノンパラメト

プルを獲得する方法である．すでに，筆者等は MCMC

リックな条件付システム誤差分布を仮定しているため，法を用いて，マルコフ推移確率を効率的にベイズ推計

サンプル数が十分であれば，条件付システム誤差分布することを明らかにしている．本研究では，筆者等がの特定化誤差を回避することができる．隠れマルコフ

提案したワイブル劣化モデルのベイズ推計法11) を拡張

劣化モデルを推計する場合，システム誤差パラメータ

し，MCMC 法を用いてシステム誤差を考慮した隠れマiãm
を直接推計することとなる．

ルコフ劣化モデルを推計する方法を提案する．

なお，測定誤差には，システム誤差だけでなく，ラ隠れマルコフ劣化モデルを含む混合分布モデルの推ンダム誤差も存在する．しかし，識別性問題が存在す

計では，前述したように尤度関数が特殊な形をしているため，両者を分離計測することは困難である．両者るため，通常の最尤法やベイズ推計法を用いることがの測定誤差を分離計測する 1 つの方法は，例えば，確率困難である23) ; 24) ．このようなことから，混合分布モデ

的フロンティア理論で採用されているように 2 種類のルの推計方法として，通常の尤度関数ではなく，完備

誤差項の間に構造的関係を想定することが必要となる

化尤度関数を定義するとともに，MCMC 法を用いて混28) ; 29) ．そのためには，両者の誤差分布をパラメトリッ

合分布モデルを推計する方法が提案されている22) Ä 25) ．

クな条件付誤差分布で表現し，両者の関係構造を明示

しかし，既往の隠れマルコフ劣化モデルでは，マルコ的にモデル化することが必要となる．しかし，現時刻

フ推移確率が定数で与えられ，これらの定数パラメーにおいては，誤差発生メカニズムに関する先験的情報

タを集計的に推計するに留まっている．しかし，本研が存在せず，本研究ではノンパラメトリックな条件付

究では，3．で言及したように，多段階指数劣化ハザー

システム誤差分布を用いている．測定誤差の分離計測

ドモデルを用いて，マルコフ推移確率を推計する点にの問題に関しては，誤差発生メカニズムに関する知見特徴がある．このような隠れマルコフ劣化モデルを推の蓄積を待たざるを得ない．

計するためには，既往の隠れマルコフ劣化モデルを推計するための MCMC 法の中に，マルコフ推移確率のベ
イズ推計アルゴリズムを内包したような MCMC アルゴ

5．推計方法リズムを開発することが必要になる．

(1) ＭＣＭＣ法(2) 完備化尤度関数の定式化

伝統的なベイズ統計学では，共役な事前・事後分布を用いて，パラメータを推計する方法が採用される21) ．

土木施設 k (k = 1; ÅÅÅ; K) に対して，直近の更新時刻

しかし，ハザードモデルの場合，簡単な指数ハザードモ

以降，それぞれ合計 T k 回にわたる測定結果が得られた

デルを用いても，共役事前確率分布が存在しないことがとしよう．測定時刻útk (t = 1; ÅÅÅ; T k ) において測定さ知られている30) ．共役事前確率分布が存在しない場合，
数値解析により多重積分を求めることが必要となる．これた見かけの健全度の測定値をm(úñ tk ) と表す．ここに，PKñ
記号 は，測定値であることを表す．また，
Tk
k=1

9

ñK ) と表ñ = (ò
ñ1 ; ÅÅÅ; ò
個の測定情報に関するデータをòñk ; ÅÅÅ; òñk k ) である．また，ò
ñk = (òñk =そう．ただし，ò度関数を採用せざるを得ない．事前確率密度関数設定には，任意性が介在せざるを得ないが，サンプル数が1Ttñkt ) であり，m(útk ) = m(m
ñkt ; zñtk ; x
ñkt は時刻útk における施kk
設 k の見かけの健全度，zñtk = útk Ä útÄ
1 は測定時刻ú
t とk
ñkt は，時刻útl における土木施設 kútÄ
1 の間の測定間隔，x増加するにつれて事前確率密度関数の特定化の影響は次第に低下する．
まず，条件付システム誤差分布関数 (25) に含まれるパラメータãi = (ã1i ; ÅÅÅ; ãii ) は，式 (26a),(26b) を満足の特性ベクトルである．この時，隠れマルコフ劣化モ
デルの完備化尤度関数は K 個の土木施設から得られたする定数である．これらの定数の事前確率密度関数とデータの完備化同時生起確率

して，ディリクレ分布を仮定しよう．すなわち，パラ

ñ)
~
L(ã;
s,å; ò
K nY
T
Y
k

=

K hY
Tk
Y

k=1 t=1=

ôsk skt (ñ
ztk )
tÄ 1
sk
t
m
ñk
t

ã

X n
sk
t

oi
exp(Äílk zñtk )k=1 t=1

T
Y

l=sk
tÄ 1

t=1

メータãi の事前確率密度関数をok

ëi (ãi jói ) = âi (ói )fskt (mñkt jãskt )

lÄ
Y1

i=sk
;6=l
tÄ 1

âi (ói ) =

íik

i
X

íik Ä ílk
(27)

i
Y

i

i óm Ä 1
(ãm
)

(30)

m=1

Ä(ó1i + ÅÅÅ+ óii )Ä(ó1i ) ÅÅÅÄ(óii )

i
ãm
=1

m=1

と表現する．ただし，ói = (ó1i ; ÅÅÅ; óii ) は定数パラメーで表される．ただし，íik = exp(xk å0i ) と表される．この尤度関数を用いて，未知パラメータベクトルã =タベクトルである．ãi がディリクレ分布に従う場合，こ(ã1 ; ÅÅÅ; ãIÄ1 )，å = (å1 ; ÅÅÅ; åIÄ 1 ), 潜在変数ベク
トルs = (s1 ; ÅÅÅ; sK ) を最尤推計する問題に帰着されらのパラメータは自動的に式 (26a),(26b) を満足する．

れる．なお，パラメータベクトルã，åを与件とすれ密度関数として用いられる多次元正規分布に従うと仮ば，潜在変数 skt (t = 1; ÅÅÅ; T k ; k = 1; ÅÅÅ; K) の定しよう．すなわち，åi ò NM (êi ; Ü i ) である．ただつぎに，åi の事前確率密度関数が，標準的な事前確率

全条件付事 後分布を求めることができる．ここで，skÄ t = (sk1 ; ÅÅÅ; sktÄ 1 ; skt+1 ; ÅÅÅ; skT k ) とすれば，skt (skt 2fsktÄ 1 ; ÅÅÅ; skt+1 g) の全条件付事後確率は，

h(åi jêi ; Ü i ) =k!itfi (mkt jãi )

t+1
j=sk
tÄ 1

T kÄ 1(28)

k f (mk jã )
!jt
j
j
t

と表される．ただし，
8
>
ô ô
>
< 1j js2kôsktÄ 1 j ôjskt+1
!jt =
>
>
: ôsk jt=1
2 î t î Tk

1
p

o
1
(31)
Ä (åi Ä êi )Ü i Ä 1 (åi Ä êi )0
2となる．ただし，êi は NM (êi ; Ü i ) の事前期待値ベクト
ル，Ü i は事前分散共分散行列である．この時，完備化
事後確率密度関数ö(ã; åjs; ò) は，
ÅexpProbfskt = ijskÄ t ; ã; òg= Psk

し，M 次元正規分布 NM (êi ; Ü i ) の確率密度関数は，

n

(2ô)

M
2

~
ö(ã; åjs; ò) / L(ã;å; s; ò)oëi (ãi jói )

(29)

2 k
sk
K
T
t
n
Y
Y
X
4
/

t = Tkである．

t=1 l=sklÄY1

jÜ i jIÄY1ni=1h(åi jñi ; Ü i )

o
íik
k k
exp(Äíz)l tíik Ä ílk

IÄ
n 1
o
Y1
Å exp Ä (åi Ä êi )Ü i Ä1 (åi Ä êi )0
2
i=1
1†
0 k
!3
i
T
I Y
Y
Y
k
i
s
i óm
(ãm
) Ä1 5
Å@ ã t k A
k=1(3) ベイズ推計法
一般に，ベイズ推計法は，1) 事前の経験情報などに基づいて，パラメータã，åの事前確率密度関数を設定
する．2) 新しく獲得したデータòに基づいて尤度関数ñ) を定義する．さらに，3) ベイズの定理に基L(ã; å; ò
づいて事前確率密度関数を修正し，パラメータã，åに

t=1

tÄ 1

m
ñt

i=sk
;6=l
tÄ 1

(32)

i=1 m=1となる．

関する事後確率密度関数ö(ã; åjò) を得る，という手順を採用することになる．以上の手順を，本研究ではベ(4) ギブスサンプリング

イズ推計ルールと呼ぶ．最尤法と異なり，未知パラメー
タã; åの確率分布が，事後分布として求まる点にベイ

隠れマルコフ劣化モデルでは，事後確率密度関数

ズ推計法の特徴がある．前述したように，ハザードモ
デルでは，共役事前確率密度関数を見出すことは不可ö(ã; åjò) を直接解析的に求めることができない．そこで，代表的な MCMC 法であるギブスサンプリング法31)

能30) であり，事前確率密度関数は，非共役事前確率密を用いて，パラメータã，åの標本サンプルを事後確率

10

IÄ
Y1

ステップ１ 初 期 値 設 置

exp

i=1

(0)

・パ ラ メ ー タ 推 定 量 ： ã (0 ) , å・サ ン プ ル 数 ： n , n
・潜 在 変 数 : s (0 )ii
・事 前 分 布 の パ ラ メ ー タ : ó , ê , Ü ino1Ä (åi Ä êi )Ü i Ä 1 (åi Ä êi )0 (33b)2と表せる．また，潜在変数sの全条件付事後分布は式 (28)で表される．以上のギブスサンプラーを用いた MCMC法により，隠れマルコフ劣化モデルを推計することがã (n ) の 標 本 抽 出

ステップ２i;(n Ä 1) の 計 算
１） N m可能となる．図-4 に，隠れマルコフ劣化モデルをベイñms (n Ä 1 )

式 (35)

ズ推計するための具体的手順を整理している．読者の便宜を図るために，同図中には，推計法の詳細を説明)
2 ） ã (nの 抽 出 式 (34)i

する節番号や式番号を明記している．以下では，以上
3） ステップの 終了 判定

i= Iの各ステップの内容を，より詳細に説明する．
i = i+1

No

ス テ ッ プ ４ へY esａ) ステップ 1

初期値設定

事 前 分 布 (30),(31) の パ ラ メ ー タ ベ ク ト ル (行 列)ス テ ッ プ 3へi

ó (i = 1; ÅÅÅ; I ), êi ; Ü i (i = 1; ÅÅÅ; I Ä 1) の値を任意に

å(n ) の 標 本 抽 出

ステップ３

１） ö
ö(å
^( å1 1 jåÄ 1 1(n )(n Ä 1); s (n Ä 1) ; ò) か ら å1(n1 ) を 発 生 式 (36)(n )(n Ä 1); s (n Ä 1) ; ò) か ら å1(n2 ) を 発 生 式 (36)

2） ö
ö(å
^( å1 2 jåÄ 1 2

設定する．潜在変数の初期値s(0) = (s(1;0) ; ÅÅÅ; s(K;0) )を設定する．ただし，s(k;0) = (sk;0; sk;01 ; ÅÅÅ
T ) であり，

1 î sk;0î ÅÅÅî sk;01
T î I を満足する．さらに，パラメータ推計量の初期値ã(0) ，å(0) を任意に設定する．これらの初期値の影響は，MCMC 法によるシミュレーショ
ン回数が蓄積されるにつれ，次第に薄れていく．MCMCのサンプル標本回数 n を n = 1 とする．

..
.
(n )

(n Ä 1 )

(n )

3） ö
ö(å
^( åI Ä 1;M jåÄ (I Ä 1M ) ; s (n Ä 1 ) ; ò) か ら åI Ä 1 M を 発 生 式 (36)１）ｂ) ステップ 2

s (n ) の 標 本 抽 出

ス テップ4s kt ;nkの 決 定 ( t = 1 ; ÅÅÅ; T ) 式 (37)パラメータã(n) の標本抽出ステップ 2 では，潜在変数s(nÄ 1) を与件とし，条件付
(n)

k = K
ス テ ッ プ ４ へY es(n)に関するパラメータ標本を獲得する．ただし，ãi =i;n: m = 1; ÅÅÅ; i) と表記する．ステップ 2 で用いる
(ãmk = k+1

No

ス テ ッ プ 5へ

ギブスサンプラーは式 (33a) で表される完備化条件付事

後密度関数ö(ã(n) js(nÄ 1) ; ò) で与えられる．真の健全度s(nÄ 1) と測定データòを与件としたとき，完備化された(n)
ãi に関するギブスサンプラーをステップ 5 ア ル ゴリズ ム の 終 了 判 定・ å (n ) の 記 録

n = n+1

n ï n

No

ö
~(ãi js(nÄ 1) ; ò)
8
9(
)
K Y
i
Tk
<Y
Y
k;(nÄ 1) =i
sti;n óm
Ä1
/
ãm
(ãm )
ñk
t
:
;

Y es

(n)

終了

図-4

隠れマルコフ劣化モデルのベイズ推計法iY

k=1 t=1=

密度関数から抽出する．式 (32) において，ã，åは互い

i
I Y
nY

sk
t
m
ñk
t

ã

k=1 t=1ëi(ãim )óm Ä 1

K hY
T
nY
X
i=1 m=1kö(åjs; ò) /
h

lÄ
Y1

i=sk
;6=l
tÄ 1

sk
t

i;(nÄ 1)Ä1

(34)i;(nÄ 1)は見かけの健全ñと潜在変数s
度の測定値mを与件とした時，nok;(nÄ 1)i;(nÄ 1)Nm=# m
ñkt = m \ st=i(35)(nÄ 1)ok=1 t=1 l=sk

tÄ 1

íik

i

óm +Nm(ãi;nm )と表すことができる．ただし，Nm

後密度関数ö(ãjs; ò)，ö(åjs; ò) は
K Y
Tk
êY

m=1

m=1

に独立であり，これらのパラメータの完備化条件付事

ö(ãjs; ò) /(n)

システム誤差分布のパラメータã(n) = (ã1 ; ÅÅÅ; ãIÄ 1 )2） ステップの終了判定

exp(Äílk ztk )íik Ä ílkと定義される．ただし，#f(33a)

g は，括弧 f

g 内に含まれる定義式が成立するような測定サンプル数を表す．式i;(nÄ1)i(34) は，パラメータóm
+ Nm
Ä 1 を有するディリ
クレ分布に他ならない．更新されたディリクレ分布 (34)を用いて，ギブスサンプリングにより，条件付システ
(n)
i;(n)
i;(n)
ム誤差分布のパラメータ標本ãi = (ã1 ; ÅÅÅ; ãi )を標本抽出する．すべての i (i = 1; ÅÅÅ; I ) に対してパ(n)ラメータ標本ãi を求める．

io11ｃ) ステップ 3(n)パラメータå(n) の標本抽出(nÄ1)

個の条件付き事後確率密度関数ö
^(åem jåÄem ; s(nÄ1) ; ò)を求めることが必要となる．本研究では，式 (36) からステップ 3 では，真の健全度で定義される多段階指数ハザードモデルのパラメータ標本を抽出する．ステッ

事後分布のパラメータåの標本をサンプリングする手法

プ 3 のアルゴリズムを説明するために，åからåem を除として，適応的棄却サンプリング33) を用いる．

いた未知パラメータベクトルをåÄ em と表そう．この時，ｄ) ステップ 4式 (33b) より，åÄ em を既知とした時のåem の条件付き

事後確率密度関数ö(åem jåÄ em ; s; ò) はö^(åem jåÄ em ; s; ò)

j
e Y
T njÄ
I Y
K Y
Y
Y1
tk
tk X
/
(ílk )éij Ä éiek;nの全条件付事後確率は，

Åg)k;(nÄ 1)

h=i

tk

IÄ
n 1
oi
Y1
Å exp Ä (åi Ä êi )Ü i Ä 1 (åi Ä êi )0
2
i=1

I Y
K Y
e Y
T k hjÄ
Y
Y1
tktk
fexp(åem xkm )géij Ä éie
/

tk
iéij
1
k k
exp(
Äí
z
)
h t
ík Ä íhk
h=i l=i;6=h l
n õmm
o
^m )2
exp Ä e (åem Ä êe2MX^m = êm +ê(åeh Ä êeh )õehm
eek;(nÄ 1)(sk;n
2 fsk;n; st+1ttÄ 1 ; ÅÅÅProbfskt = ijã; sÄ t; òg
8k;(nÄ 1)(n)k!itfi (mt jãi )>>k;(nÄ 1)>Ps>k;(nÄ 1)
(n)
2
>
!jt
fj (mk>t jãj )
j=1>>>>
(t = 1 の時)>>>
k;(nÄ 1)(n)>!
fi (mk
)
>
t jãi
it
>
< Psk;(nÄ 1)k;(nÄ 1)
(n)
t+1
!jt
fj (mkt jãj )
=
k;n
j=s
tÄ 1
>
>
>
>
(2 î t < T k の時)>>>
k;(nÄ 1)(n)>!
fi (mk)> P itt jãi>>Ik;(nÄ 1)k jã(n) )
>
!
f
(m
j
k;n
>
t
jt
j
j=s
>
tÄ 1
>
>
:
k
(t = T の時)

oéij
1
k k
exp(
Äí
z
)
h
t
ík Ä íhk
l=i;6=h l
hÄ
Y1

l=i

全条件付事後確率 (28) に基づいて，新しい潜在変数s(n) をランダムサンプリングする．いま，潜在変数ベクトk;(nÄ 1)k;(nÄ 1)k;(nÄ 1)ルsÄt= (sk;nÅÅ; sk;n; ÅÅÅ; sT k
)を1 ;ÅtÄ1 ; st+1
定義する．この時，stk

i=1 j=e k=1 t=1

潜在変数の更新jhÄX
Y1

i=1 j=e k=1 t=1 l=iと表される．ただし，
8>
ô1j ôjsk;(nÄ 1)>>
2<
k;(nÄ 1)
ôsk;n j ôjsk;(nÄ 1)
!jt
=
tÄ 1
t+1
>
>
>
: ô k;n(36)h=1;6=mstktkと表せる．ただし，éie，éijは(
1 sktÄ 1 = i = e の時tkéie=
0 それ以外の時(
1 sktÄ 1 = i; skt = j の時tkéij=
0 それ以外の時T kÄ 1

j

(37)

t=1

2 î t < T k (38)t = Tkである．すべての k (k = 1; ÅÅÅ; K ) に対して，t = 1 よk;n

り逐次，潜在変数 stｅ) ステップ 5(t = 1; ÅÅÅ; T k ) を求める．

アルゴリズムの終了判定以上で求めたパラメータ推計量の更新値ã(n) ,å(n) ，潜
在変数の更新値s(n) を記録する．n î nの場合，n = n+1

として，ステップ 2 へ戻る．そうでない場合，アルゴリとなるダミー変数である．êem は事前期待値ベクトルêe のズムを終了する．

第 m 要素であり，õehm は事前分散共分散行列Ü e Ä 1 の第P
(h; m) 要素である．また， M
h=1;6=m は 1 から M までの

なお，以上のアルゴリズムの初期段階においては，パラメータの初期値設定の影響が残存している．このた要素のうち m を除いた要素の総和を意味する．これらのめ，シミュレーション回数 n が十分大きなnの到達する

条件付き確率密度関数から標本を発生させ，その標本をまでのパラメータ標本を除去することが望ましい．こ

用いてパラメータåの事後分布に関する各種の統計量をのため，ギブスサンプリングで求めたã(n) ; å(n) (n =計算することができる．MCMC 法から得られた標本を

n + 1; n + 2; ÅÅÅ; n) を，事後確率密度関数ö(ã; åjò) からの標本と見なすことができる．したがって，これらの標本を用いて，パラメータベクトルã; åの事後分布

用いて，事後分布の各種統計量を求める方法については，
(n)

(n)

5． (5) で述べる．このとき，å(n) = (å11 ; ÅÅÅ; åIÄ 1M )を以下の手順でランダムサンプリングする．
(n) (nÄ1)(n)・ステップ 3-1 ö^(å11 jåÄ11 ; s(nÄ 1) ; ò) からå11 を乱数発生する．
(n) (nÄ1)(n)・ステップ 3-2 ö^(å12 jåÄ12 ; s(nÄ 1) ; ò) からå12 を乱に関する各種の統計量を計算することも可能となる．

(5) 事後分布に関する統計量

MCMC 法によって得られた標本に基づいて，パラメ
ータベクトルã; åに関する統計的性質を分析することができる．MCMC 法を用いた場合，パラメータの事後確率密度関数ö(ã; åjò) を解析的な関数として表現することはできない．得られた標本を用いてノンパラメトリック数発生する．
・ステップ 3-3

返す．・ステップ 3-4(n)以下，同様の手順を繰り(n)(nÄ 1)
ö^(åIÄ 1;M jåÄ (IÄ 1M ) ; s(nÄ 1) ; ò)からåIÄ 1M を乱数発生する．

なお，ギブスサンプリングを行うためには (I Ä 1) × M

12に分布関数や密度関数を推計することとなる．いま，ギ
表-1
ランク
1
2
3
4
5

ブスサンプリングから得られた標本を (ã(n) ; å(n) ) (n =1; ÅÅÅ; n) と表そう．この内，最初のn個の標本は収束過程からの標本と考え，標本集合から除去する．その上で，パラメータの標本添字集合を M = fn + 1; ÅÅÅ; n)と定義しよう．このとき，パラメータã及びåの同時確率分布関数 G(ã); G(å) は#(ã(n) î ã; n 2 M)
G(å) =(39a)nÄn#(å(n) î å; n 2 M)
G(å) =(39b)nÄn
と表すことができる．ただし，#(å(n) î å; n 2 M) は論理式å(n) î å; n 2 M が成立するサンプルの総数である．また，パラメータåi の事後分布の期待値ベクト~ (å )，分散・共分散行列Ü~ i (å ) は，それぞれ
ルê
i
i
i

表-2
事前
健全度
1
2
3
4
5

~(åi;1 ); ÅÅÅ; ê~ (å ) = (ê~(åi;M ))0
ê
i
i
n
n
ê X
X
åi;1 (n)åi;M (n) ë 0=; ÅÅÅ;
(40a)
nÄn
nÄn
n=n+1
n=n+1
0
1
õ
~ 2 (åi;1 )ÅÅÅ õ~ (åi;1 åi;M )B .
C. . ..
C (40b).
Ü~ i (åi ) = B
. .
@ .
A
2
~ (åi;M åi;1 ) ÅÅÅ õõ~ (åi;M )

5
0
47
44
67
84

予測を試みる．対象とするデータベースは，補修年度
に関するデータと，路面性状の実測データにより構成されている．路面性状調査は，3 年もしくは 6 年間隔で行われている．分析対象とする区間は，M 県に存在する一般国道のうち，M 県が管理する約 60km の区間であり，100m を単位区間として測定サンプルが蓄積されている．また，1986 年以降，現在に至るまで，高速測定(41a)

車両による舗装のわだち掘れ測定を定期的に実施している．本研究の適用事例で用いる管理データベースで~ (åi;m åi;l )
õ
=

1
331
573
66
50
2

サンプル数
事後健全度
2
3
4
339
32
5
1919 468 187
240
382 163
63
52
82
22
16
27

面性状調査結果に適用し，道路舗装のわだち掘れ劣化と表される．ただし，nX~(åi;m )g2
fåi;m (n) Ä êõ~ 2 (åi;m ) =
nÄn
n=n+1

健全度レーティング
わだち掘れ量5mm 未満5mm 以上 10 未満
10mm 以上 15mm 未満
15mm 以上 20mm 未満
20mm 以上nX~(åi;m )gfåi;l (n) Ä ê~(åi;l )gfåi;m (n) Ä ê
nÄn
n=n+1

は，1998 年から 2005 年に至る期間中に実施された路面性状調査データが記載されている．これらのデータの中で，適用事例では，アスファルト舗装のわだち掘れ(41b)

過程に着目する．さらに，補修履歴データベースには，である．また，ギブスサンプリングによる標本を用い

分析対象期間中に実施された道路舗装の補修に関するて，パラメータã，åの信頼区間を定義できる．たとえ記録が蓄積されている．対象期間中に道路舗装の補修

ば，パラメータåの 100(1 Ä 2")% 信頼区間は，標本順序が実施された路面箇所に関しては，その時点を初期時"
統計量 (å"i;m ; åi;m ) (i = 1; ÅÅÅ; I Ä 1;
å"i;m = arg maxÉåi;m (n )"öm = 1; ÅÅÅ; M )õ#(åi;m (n) î åi;m (nÉ); n 2 M)
î"nÄn点とするサンプル時間軸を定義し，隠れマルコフ劣化
ハザードモデル推計のためのデータベースを作成した．
したがって，測定サンプルは，初期時点とそれ以降の連続する 2 つの測定時点におけるわだち掘れ量の実測(42a)値，および測定間隔に関する情報により構成されていåi;m = arg minÉÉ
åi;m (n )öõ#(åi;m (n) ï åi;m (nÉÉ); n 2 M)
î " (42b)nÄn

る．データベースには，個々の測定サンプルに関わる交通情報や，道路条件も記載されている．このようにして作成した測定サンプル数は，合計 5,261 個である．

"
を用いてå"i;m < åi;m < åi;m と定義できる．

舗装のわだち掘れ量の実測データは連続値として記録されている．本研究で提案した劣化予測モデルを推計するためには，わだち掘れ量の連続値を離散的な健6．適用事例

全度に置き換える必要がある．本研究では，わだち掘(1) データベースの概要

れに関する健全度を，表-1 に示すような，5 段階のレー本研究で提案した隠れマルコフ劣化モデルを，M 県

ティングに分類した． いま，連続する 2 つの時点におが管理する一般国道において過去 9 年間に実施された路ける健全度の実測値を，事前健全度，事後健全度と呼13

表-3 パラメータの推計結果
健全度 定数項
交通量
åi1
åi2
1
0.279
0.418
(16.135)
0.033
0.187
2
(223.255)
0.107
3
0.112
4
-

真の健全度１

真の健全度２

真の健全度４

真の健全度５

真の健全度３

1

0.8

割合

0.6

0.4
0.2

注) 括弧内は尤度比検定統計量を表している. また，損傷度 5 はマルコフ連鎖の吸収状態であり，ハザード率が0 となる．

0
1

2

3

4

5

見かけの健全度

図-5 システム誤差分布

ぶこととする．多段階指数ハザードに用いる測定サン
プルを，事前健全度と事後健全度に着目して整理した
結果を表-2 に示している．同表の各行は，事前健全度 i

以上のデータに基づいて推計した隠れマルコフ劣化に，各列は事後健全度 j に対応している．ここでは，2モデルのパラメータ推計値を，一括して表-3 に示してつの測定時点の間に，舗装の補修が実施されていないいる．ここで，説明変数 xi1 = 1 は恒常的に値 1 をとり，サンプルのみを抽出している．しかし，同表に示すよåi1 は定数項を表す．本研究では，説明変数の選定に際して，各説明変数の説明力を尤度比検定統計量を用いて評価した．本来，尤度比検定は最尤推計法において用いられる検定法であるが，本研究では MCMC 法で発うに，事前健全度に対して，事後健全度が回復しているサンプルが多数存在している．特に，事前健全度が2 の場合に事後健全度が 1 となるサンプルが多い．2 つの時点間において，道路舗装の補修が実施されていないことを考慮すれば，分析で用いるデータには，数多くのシステム的誤差が介在していることが理解できる．

生したパラメータサンプルそれぞれに対して尤度を求めるとともに，その平均値を用いて尤度比検定を実施
した．表-3 には，モデルの推計精度を表す尤度比を併記している．尤度比検定の結果，本適用事例では，xi1 ，(2) 推計結果
以上のデータベースを用いて，舗装わだち掘れ量の劣化予測モデルを推計した．推計にあたっては，道路特性を表す説明変数として，車線区分，大型車交通量，
構造特性等を説明変数の候補としてとりあげた．これらの説明変数のうち，大型車交通量 xi2 のみが，最終的xi2 が説明変数として採用されている．同表に示すように，大型車交通量は，道路舗装のわだち掘れの発生と初期段階の進行に多大な影響を及ぼす．しかし，わだ
ち掘れが進行した段階では，大型車交通量と関係なく
わだち掘れが進行することがわかる．つぎに，推計したシステム誤差分布を図-5 に示している．同図は，真に有意な説明力を有する変数として採用された．なお，の健全度に対する見かけの健全度の確率分布 fi (mjãi )説明変数 xi2 は大型自動車交通量について区間内での最を示している．同図より，以下の事項が読み取れる．す大値を 1 として基準化したものである．大型車交通量以

なわち，真の健全度と見かけの健全度が一致する確率外の変数は，有意な説明力を持たないため，説明変数

がもっとも大きい．健全度 4，5 の場合，システム的誤として選ばれていない．本研究では，M 県の国道とい差の影響により，見かけの健全度が広範囲に分散してう単一の路線を対象としたわだち掘れ予測モデルを推いることが理解できる．

計している．このため，舗装材料，施工時期等の説明

つぎに，推計されたパラメータを用いて，式 (11) よ

変数は，モデル推計に用いる各サンプルを通じて同一り各劣化状態のハザード率の期待値と期待寿命を算出の値をとっている．したがって，これらの舗装特性は，した．その結果を，表-4 に示している．同表より，平モデルの説明変数として選ばれていない．当然のこと均的なわだち掘れの進行過程は，初期時点から 3 年未満

ながら，推計したわだち掘れ予測モデルは，分析対象の早期段階で健全度 2 に推移し，その後は 10 年から 15とした一般国道にのみ適用可能であることは言うまで年程度の間隔で，健全度が低下していく．

もない．今後，環境条件が異なる多数の路線にも適用(3) 分析結果の考察

可能なわだち掘れ予測モデルを開発するためには，本適用事例で採択されなかった変数を説明変数に加えた

隠れマルコフ劣化モデルを用いて，マルコフ推移確モデルを開発することが必要となる．

率行列を求めよう．本研究で提案したハザードモデル14

表-4 レーティング期待寿命k
劣化状態E[íil ] E[RM Dil
](年)
1
0.362
2.763
0.070
14.225
2
0.107
9.305
3
0.112
8.951
4

表-6

注) 損傷度 5 はマルコフ連鎖の吸収状態であり，ハザー
ド率，期待寿命は定義されない．
表-5 推移確率行列 (隠れマルコフ劣化モデル)
事前
事後健全度
健全度
1
2
3
4
5
1
0.696 0.293 0.011 0.000 0.000
0.0
0.932 0.064 0.003 0.000
2
0.0
0.0
0.898 0.096 0.006
3
0.0
0.0
0.0
0.894 0.106
4
0.0
0.0
0.0
0.0
1.0
5
注) 推移行列は，1 年間の間に生起する状態推移確率を示している．ここでは，平均操作を行い，当該路線の平
均的なひび割れ発生，推移確率を求めたものである．

切断データベースを用いた推計結果
健全度 定数項
交通量
åi1
åi2
1
0.242
0.350
(65.412)
0.040
0.196
2
(183.223)
0.139
3
0.163
4
-

注) 括弧内は尤度比検定統計量を表している. また，損傷度 5 はマルコフ連鎖の吸収状態であり，ハザード率が0 となる．同表では，吸収状態のハザード率を記述していない．

て，大型車交通量の影響が大きいことがわかる．具体的には，健全度 2 の平均寿命は，交通量最大区間では約

4.5 年であるのに対して，交通量最小区間では約 30 年となっている．ただし，前述したように，交通量の影響
は供用開始時点から健全度 3 に至る過程にしか及ばず，

0

10

経過時間(年)2030
40

健全度 3 以降では交通量に関わらず同じ速度で劣化する
50

結果となっている．

60

1(4) システム誤差と推計バイアス

健全度

2

システム的誤差の存在が，多段階指数劣化ハザード

3

モデルの推計結果に及ぼす影響を分析する．いま，図-14に示すようなデータに対して 45 度線より下方に位置する（事後健全度が事前健全度より改善した）サンプルの

5
BMケース

交通量maxみを除去して，多段階指数ハザードモデルのベイズ推

交通量min計を試みた．多段階指数ハザードモデルのベイズ推計法の詳細は参考文献12) に譲る．このように，事後健全
図-6 供用性曲線度が改善したサンプルを除去して作成したデータベー
ス（以下，切断データベースと呼ぶ）に基づいて，多段階指数劣化ハザードモデルを用いて推計した結果をは，説明変数の組み合わせごとにハザード率を定義す

表-6 に示している．さらに，平均的ハザード率を用いることができる．言い換えれば，道路特性，車線特性別て作成したマルコフ推移確率を表-7 に示している．隠のマルコフ推移確率を推計することができる．ここでれマルコフ劣化モデルを用いて作成したマルコフ推移

は，すべての道路区間にわたってハザード率を平均化

確率（表-5）と，切断データベースを用いて推計した

したような平均的なハザード率を用いてマルコフ推移

マルコフ推移確率（表-7）を比較しよう．2 つの表を比確率を求めた．その結果を表-5 に示している．ただし，較すれば，健全度 1 の場合を除いて，その他の健全度でマルコフ推移確率行列は，1 年間隔で定義されている．

は隠れマルコフ劣化モデルを用いた方が，同一健全度

つぎに，式 (13) を用いて，供用性曲線を求めた．そのにとどまる確率が大きくなっている．図-7 には，平均結果を図-6 に示す．同図には，対象区間全体の平均ハハザード率を用いて作成した供用性曲線を示している．

ザード率を用いた供用性曲線 (BM ケースと呼ぶ) を示

マルコフ推移確率の比較結果からも類推できるように，している．さらに，大型車交通量がひび割れ進行に及ぼ

切断データベースを用いた場合，劣化速度を過大に評す影響を分析するために，観測された最大交通量と最価する危険性が存在することが理解できる． このこと

小交通量を用いてハザード率を算出し，それぞれのハは，事前健全度と事後健全度が逆転するようなサンプ

ザード率に対応する供用性曲線を求めた結果も併記しルにおいては，真の事後健全度が真の事前健全度より，ている．同図から，健全度 2 から健全度 3 の推移におい劣化がそれほど進展していない可能性が大きいことを15

いように思える．しかし，このことは本適用事例にの表-7 推移確率行列 (多段階指数劣化ハザードモデル)
事前
事後健全度
健全度
1
2
3
4
5
1
0.728 0.261 0.010 0.000 0.000
0.0
0.927 0.068 0.005 0.000
2
0.0
0.0
0.869 0.121 0.010
3
0.0
0.0
0.0
0.850 0.150
4
0.0
0.0
0.0
0.0
1.0
5
注) 推移行列は，1 年間の間に生起する状態推移確率を示している．ここでは，平均操作を行い，当該路線の平
均的なひび割れ発生，推移確率を求めたものである．なお，多段階指数ハザードモデルの推計には，健全度が回復しているサンプルを用いていない点に留意されたい．

み成立する事項であり，システム的誤差の存在が舗装
劣化予測に及ぼす影響に関して，今後経験的な知見を増やす努力が必要である．

(5) シミュレーションによる再現性の検証
隠れマルコフ劣化モデルは，見かけの健全度データを用いて，真の健全度の間に成立する推移関係を推計することを目的とする．モデルの推計に用いる健全度にシステム的誤差が存在しているため，実測値を用いて隠れマルコフ劣化モデルの推計精度を直接検討することができない．そこで，本節では，モンテカルロシ

0

10

経過時間(年)
20

ミュレーションにより人為的に作成したデータベース
30を用いて，隠れマルコフ劣化モデルの現象再現能力に40

1ついて分析することとする．
いま，表-6 に示す多段階指数劣化ハザードモデルが，健全度

2

真のモデルであると仮定しよう．隠れマルコフ劣化モ
デルの推計で用いたデータベースの各サンプルに対し

3て，説明変数（大型車交通量）を入力変数として，表-54に示す多段階指数劣化ハザードモデルを用いれば，式(8),(9) より各サンプルに対してマルコフ推移確率ôij (z )を求めることができる．サンプル k (k = 1; ÅÅÅ; K ) に対して，時点útk (t = 1; ÅÅÅ; T ) で測定した時に得られる仮想的な健全度を，以下の手順で作成した．すなわ
図-7 供用性曲線
ち，１）初期時点から時点ú1k の推移確率ô1j (z1k ) を用い^ (úk ) = ^iをランダムにて，時点úik における真の健全度h1k発生する．ただし，z1 は，初期時点から時点 t = 1 まで意味している．したがって，切断データベースを用い^ (úk ) = ^iに対して，の時間間隔である．２）真の健全度hた場合，2 つの測定時点の間で劣化の進展が遅いサンプ
1^ (úk ) = ^j をラン
推移確率ô^ij (z2k ) を用いて，真の健全度hルが，データベースから除去される可能性がある．そ2
ダムに発生する．以下，同様の手順により，真の健全度の結果，切断データベースを用いて劣化予測を行った^ (úk ) までランダムに発生する．３）システム誤差分布
場合，舗装の劣化速度を過大に評価する可能性がある． hT^ (úk ) (t = 1; ÅÅÅ; T )（図-5 参照）を用いて，真の健全度h
しかし，適用事例に関する限り，健全度 3 に至るまでのtに対して，見かけの健全度m^ (útk ) をランダムに発生す供用性曲線に関しては，切断データベースを用いて推る．このようにして作成した人工的なデータベースを計した供用性曲線とシステム的誤差を考慮した供用性
人為データベースと呼ぶ．さらに，人為データベース
曲線の間に大きな差異はない．劣化が進行した段階にに対して，事後健全度 m(útk ) (t = 1; ÅÅÅ; T ) が事前健全おいては，劣化速度の過大推計の程度が大きくなる．しk
度 m(útÄ
かし，健全度が進行した段階では，道路舗装が補修さ1 ) より改善しているようなサンプルを除去しれることによるサンプル欠損バイアスが発生するため， た人為切断データベースを作成する．
5

隠れマルコフ劣化モデル

多段階指数劣化ハザードモデル

そもそも供用性曲線の推計精度自体に問題が発生する

以上のように作成した人為データベースを用いて，隠可能性が存在する．舗装劣化予測の実務において，かれマルコフ劣化モデルによりハザードモデルを推計しねてより路面性状調査結果のシステム的誤差の存在がた（ケース A）．さらに，人為切断データベースを用い指摘されてきた．現実の劣化予測においては，2 つの測て，多段階指数劣化ハザードモデルを用いて，ハザード

定時点において健全度が改善されるようなデータを除モデルを推計した（ケース B ）．以上の推計結果を表-8去したような切断データベースを作成し，供用性曲線に併記している．さらに，図-8 には，ケース A とケーが作成されてきた．本適用事例に関する限り，健全度 3

B のそれぞれのケースに対して，供用性曲線を作成した結果を示している．表-5 と表-8 を比較することにより，隠れマルコフ劣化モデルを用いることにより，シス
テム誤差の影響を排除して，真のマルコフ劣化ハザー程度の比較的わだち掘れが軽微な段階では，切断デー
タベースを用いて，多段階指数劣化ハザードモデルと供用性曲線の推計を行っても，実務上大きな支障はな

16

装の劣化情報に無視できないシステム的誤差が存在す

表-8
健全度

1
2
3
4

ハザードモデルの再現結果
ケース A
ケース B
定数項
交通量
定数項
交通量
åi1
åi2
åi1
åi2
0.270
0.415
0.310
0.150
(107.673)(14.340)
0.039
0.166
0.064
0.251
(149.528)(149.095)
0.105
0.162
0.113
0.160
-

ることを明らかにした．その上で，隠れマルコフ劣化
モデルを用いて，システム誤差の背後に存在する劣化過程をマルコフ劣化モデルとして効果的に表現できる
ことが明らかになった．しかし，本研究で提案した隠れマルコフ劣化モデルに関して，今後に残された研究課題がある．第 1 に，本研究で提案した隠れマルコフ劣化
モデルは，道路舗装以外にも多様なタイプの土木施設の劣化予測に適用することが可能である．ただし，土木施設によりシステム誤差の発生メカニズムは多様に

注) 括弧内は尤度比検定統計量を表している. また，損傷度 5 はマルコフ連鎖の吸収状態であり，ハザード率が0 となる．

異なることが予想される．システム誤差の発生過程に関しては，システム誤差関数の特定化を含めて，今後の実証分析の蓄積が不可欠である．さらに，システム

0

10

経過時間(年)20

誤差の発生を可能な限り抑制しうるような測定方法に
30

40

関する分析が必要である．第 2 に，劣化予測モデルに介1在するシステム的誤差が発生する背景として，本研究で取り上げた代表値問題以外にも，多様な要因が考え

健全度2られる．測定機器や測定時刻に固有なシステム的な測
定バイアスも存在する．それ以外にも，初期時刻にお3ける更新・補修方法が，その後の劣化過程に及ぼす影響も無視できない．このような初期時刻の施工方法が

4

劣化速度に及ぼす影響を分析するためには，ハザード率が確率変動するような混合ハザードモデルを開発す

5
真の供用性曲線

ケースAることによりアプローチが可能である．第 3 に，M 県

ケースBの一般国道を対象とした実証分析では，劣化が進行するほど，切断データベースを用いた場合に，劣化速度

図-8 供用性曲線の再現結果の過大推計の程度が大きくなる．しかし，わだち掘れ
ドモデルを再現することに成功していることが理解でが軽微な段階では，切断データベースを用いても，シきる．その結果，図-8 に示すように，隠れマルコフ劣ステム的誤差による推計バイアスはそれほど大きいも化モデルを用いて，真の供用性曲線をほぼ再現できてのではないことが判明した．しかし，このことは適用いると判断できよう．しかし，人為切断データベースを事例についてのみ成立する事項であり，今後，システ用いて多段階指数ハザードモデルを推計した場合，前ム的誤差が舗装劣化予測に及ぼす影響に関する経験的述したように劣化速度を過大評価されることが，図-8

知見を蓄積することが必要である．
なお，本研究の一部は文部科学省「若手研究者の自からも読み取れる．

立的研究環境整備促進」事業によって大阪大学グローバル若手研究者フロンティア研究拠点にて実施された．

7．おわりに

参考文献

本研究では，土木施設の劣化過程をマルコフ劣化モ1) 小林潔司，上田孝行：インフラストラクチャのマネジメン
ト：研究展望，土木学会論文集，No.744/IV-61，pp.1527，2003．
2) 保田敬，小林潔司：BMS における点検結果と状態推移確率
が LCC に及ぼす影響，建設マネジメント論文集，Vol.11，
pp.111-122，2004.
3) Lee, T.C., Judge, G.G. and Zellner, A.: Estimating
the Parameters of the Markov Probability Model from
Aggregate Time Series Data, Amsterdam, NorthHolland, 1970.
4) 杉崎光一，貝戸清之，小林潔司：目視検査周期の不均一性を考慮した統計的劣化予測手法の構築，構造工学論文
集，Vol.52A，pp.781-790，2006.
5) 貝戸清之，阿部允，藤野陽三：実測データに基づく構造物

デルとして表現するとともに，健全度の測定結果に誤差が発生するメカニズムを隠れマルコフ劣化モデルを用いて表現できることを示した．その際，測定誤差として代表値問題により発生するシステム誤差に着目した．さらに，システム誤差が介在するような見かけの健全度に関する測定データに基づいて，真の健全度により定義されるマルコフ劣化モデルを，マルコフ連鎖
モンテカルロシミュレーションによりベイズ推計する
方法論を提案した．さらに，道路舗装を対象とした実証分析の結果，路面性状調査により測定される道路舗

17

の劣化予測，土木学会論文集，No.744/IV-61, pp.29-38,
2003.
6) 青木一也，山本浩司，小林潔司：劣化予測のためのハ
ザードモデルの推計，土木学会論文集，No.791/VI-67，
pp.111-124，2005.
7) Mishalani, R. and Madanat S.：Computation of infrastructure transition probabilities using stochastic duration models, ASCE Journal of Infrastructure Systems, Vol.8, No.4, 2002.
8) 津田尚胤，貝戸清之，青木一也，小林潔司：橋梁劣化予測のためのマルコフ推移確率の推定，土木学会論文集，
No.801/I-73，pp.68-82, 2005．
9) 青木一也，山本浩司，津田尚胤，小林潔司：多段階ワイブ
ル劣化ハザードモデル，土木学会論文集，No.798/VI-68，
pp.125-136, 2005．
10) 貝戸清之，熊田一彦，林秀和，小林潔司：階層型指数劣化ハザードモデルによる舗装ひび割れ過程のモデル化，土木学会論文集Ｆ, Vol. 63, No. 3, pp.386-402，2007.
11) 津田尚胤，貝戸清之，山本浩司，小林潔司：ワイブル
劣化ハザードモデルのベイズ推計法，土木学会論文集，
No.798/VI-68，pp.125-136, 2006．
12) 貝戸清之，小林潔司：マルコフ劣化ハザードモデルのベイ
ズ推定，土木学会論文集 A，Vol.63, No.2. pp.336-355,
2007.
13) 小林潔司，熊田一彦，佐藤正和，岩崎洋一郎，青木一也：
サンプル欠損を考慮した舗装劣化予測モデル，土木学会
論文集 F，Vol.63，No.1，pp.1-15，2007.
14) 貝戸清之，山本浩司，小濱健吾，岡田貢一，小林潔司：
ランダム比例ワイブル劣化ハザードモデル：大規模情報システムへの適用，土木学会論文集 F，2008.(登載決定)
15) Lancaster, T.: The Econometric Analysis of Transition Data, Cambridge University Press, 1990.
16) Gourieroux, C.: Econometrics of Qualitative Dependent Variables, Cambridge University Press, 2000.
17) MacDonald, I.L. and Zucchini, W.: Hidden Markov
and Other Models for Discrete-valued Time Series,
Chapman & Hall, 1997.
18) Hamilton, J.: A new approach to the economic analysis of nonstationary series and the business cycle,
Econometrica, Vol. 57(2), pp.357-384. 1989.
19) Diebold, F.X. and Inoue, A.: Long memory andregime switching, Journal of Econometrics, Vol.105,
pp.131-159, 2001.
20) Kim, C.-J. and Nelson, C.R.: State-Space Models with
Regime Switching : Classical and Gibbs-Sampling Approaches with Applications, MIT Press, 1999.
21) 和合肇：ベイズ計量経済分析，マルコフ連鎖モンテカル
ロ法とその応用，東洋経済新報社，2005.
22) Diebolt, J. and Robert, C.P.: Estimation of ånite

mixture distributions through Bayesian sampling,
Journal of the Royal Statistical Society, Series B,
Vol.56, pp.363-375, 1994.
23) Titterington, D.M., Smithe, A.F.M. and Makov,
U.E.: Statistical Analysis of Finite Mixture Distributions, John Wiley ＆ Sons., 1985.
24) Robert, C.P.: Mixtures of Distributions: Inference
and Estimation, in: Gillks, W.R., Richardson, S. andSpiegelhalter, D.J. (eds.): Markov Chain Monte Carloin Practice, Chapman ＆ Hall, 1996.
25) Robert, C.P., Rydìen, T. and Titterington, D.M.:
Bayesian inference in hidden Markov models throughthe reversible jump Markov chain Monte Carlo
method, Journal of the Royal Statistical Society, Series B, Vol.62, pp.57-75, 2000.
26) Dempster, A.P., Laird, N. M. and Rubin, D. B.: Maximum likelihood from incomplete data via the EM
Algorithm, Journal of the Royal Statistical Society,
Series B, Vol.39, pp.1-38, 1977.
27) Celeux, G., Hurn, M. and Robert, C.P.: Computational and inferential diéculties with mixture posterior distributions, Journal of the American Statistical
Association, Vol.95, pp.957-970, 2000.
28) Bauer, W.: Recent developments in econometric estimation of frontiers, Journal of Econometrics, Vol.46,
pp.39-56, 1990.
29) Kumbhakar, S. C. and Lovell, C.A.K.: Stochastic
Frontier Analysis, Cambridge University Press, 2000.
30) Ibrahim, J.G., Ming-Hui, C. and Sinha, D.: BayesianSurvival Analysis, Springer Series in Statics, 2001.
31) Geman, S. and Geman, D.: Stochastic relaxation,
Gibbs distributions and the Bayesian restoration ofimages, Transactions on Pattern Analysis and Machine Intelligence, Vol.6, pp.721-741, 1984.
32) Gelfand, A.E. and Smith, A.F.M.: Sampling-based
approaches to calculating marginal densities, Journal of the American Statistical Association, Vol.85,
pp.398-409, 1990.
33) Gilks, W.R. and Wild, P.: Adaptive rejection sampling for Gibbs sampling, Applied Statistics, Vol.41,
pp.337-348, 1992.
34) Geweke, J.: Evaluating the Accuracy of Samplingbased Approaches to the Calculation of Posterior Moments, in: Bernardo, J.M., Berger, J.M., Dawid,
A.P., and Smith, A.F.M. (eds.) :Bayesian Statistics
4, pp.169-193, Oxford University Press, 1996.
(平成 20 年 3 月 10 日 受付)

A HIDDEN MARKOV DETERIORATION MODEL WITH MEASUREMENTERRORS
Kiyoshi KOBAYASHI, Kiyoyuki KAITO, and Hidekazu HAYASHIThe Markov chain models have been applied in many deterioration forecasting practices. The inspection
data may include measurement errors, categorized by random errors and system one by the so-called
representation matters. In this paper, the deterioration processes with measurement errors are formulated
as a hidden Markov chain model. The Bayesian estimation of the hidden Markov chain models can bemade by use of Markov Chain Monte Carlo simulation technique. The applicability of the methodology
presented in this paper is examined against the real world data concerning the national road pavement.

18