統計数理（2015）
第 63 巻 第 1 号 83–104
c 2015 統計数理研究所
特集「地震予測と統計モデル」［研究ノート］

活断層で繰り返される地震の点過程モデルとその長期確率予測
野村 俊一†（受付 2015 年 1 月 8 日；改訂 3 月 13 日；採択 3 月 17 日）

要旨

日本には非常に多数の活断層が分布しており，それぞれが大地震を将来引き起こす危険を孕んでいる．各々の活断層はある程度周期的に活動を繰り返しており，その活動周期や前活動時期も異なることから，現時点における地震危険は活断層によりまちまちとなる．本稿では，この “繰り返し地震” の危険性を長期確率予測として評価するための統計モデルおよび予測手法を解説する．大地震の予測においては，その活動周期の長さゆえにデータの不足および信頼性の問題がつきまとう．これらは統計モデルの適用に障害となるだけでなく，予測誤差として重大な影響をもたらしうる．上記問題に起因するパラメータの推定誤差やデータの不確定性に対処し，安定した予測性能を確保する手段として，ベイズ統計を用いた予測手法を紹介する．さらに近年では，海底下で発生する小規模の繰り返し地震が地震予測の一つの鍵として注目を集めている．2011 年東北地方太平洋沖地震のような海溝型巨大地震の前ぶれとして，震源域の周辺でゆっくり滑りと呼ばれる先行活動がしばしば起こり，それが小規模の繰り返し地震活動の変調として観測されるためである．ここでも，非定常な更新過程モデルが繰り返し地震活動の変化の解析に重要な役割を果たしている．
キーワード：繰り返し地震，長期確率予測，Brownian Passage Time 分布，更新過程，
ベイズ予測．
1.

はじめに

地震とは断層の急激なずれ滑りによって起こる自然現象である．逆に言えば，断層は過去の地震活動により形成されたものであり，なかでも今後も地震活動を続けると考えられるものは活断層と呼ばれる．活断層の規模は大小多様であるが，一般的に大きい活断層は大きい地震を，
小さい活断層は小さい地震を引き起こす．それゆえ，既に知られている大きな活断層があれば，その活断層で次にいつ地震が起きるかが防災上の重要な関心事となる．本稿では，地震の発生源を特定の既知の活断層（あるいは活断層群）に絞り，そこで繰り返される地震を予測するための統計モデルと予測手法について近年の研究成果も踏まえて解説していきたい．
ここで，繰り返し地震とは，先に述べた活断層において現在まで繰り返されている地震系列のことを指す．そもそも，同じ活断層で地震が起こり続けているのは，地震を起こそうとする力（断層面に掛かる応力）
がプレート運動によって供給され続けているためである．異なるプレー
ト同士の境界では互いのぶつかり合いに伴うひずみが発生しており，活断層はひずみを正常な†東京工業大学 情報理工学研究科：〒 152–8552 東京都目黒区大岡山 2–12–1

84

統計数理

第 63 巻

第 1 号 2015

状態に戻すべくずれ滑りを繰り返している．プレートは概ね年数 cm という速度で恒常的に移動を続けており，プレート内に存在する活断層におけるひずみによる応力もそれに伴って蓄積されていく．活断層は境界面の岩石同士の静的摩擦力により普段は静止しているが，蓄積された剪断応力が摩擦力（岩石破壊強度）を超えるまでに至ると急激なずれ滑りを起こし，それが地震として観測されることとなる．日本においては，日本列島がユーラシアプレートと太平洋プ
レート，フィリピン海プレート，北米プレートとの境界に位置していることから，非常に多くの活断層が存在している（例えば，活断層研究会，1991）．活断層は大小多数あり全てを網羅す
ることは不可能であるが，これまで地質学者らの精力的な調査により，日本の陸域および海域で相当な数の活断層が発見され，繰り返し地震予測の対象とされている．
なお，日本における繰り返し地震は，海溝型地震と内陸活断層の地震とに大きく分類される．
海溝型地震は主にプレート境界である海溝付近で発生する地震であり，国内観測史上最大の地（東日本大震災）はこちらに該当する．内陸活断層の地震である 2011 年東北地方太平洋沖地震
震は，プレート境界から離れた陸域の比較的浅い場所で起こる地震であり，1995 年兵庫県南部
地震（阪神・淡路大震災）
はこちらに含まれる．両者の繰り返し地震はその規模や発生周期に大きな差異があり，海溝型地震は巨大地震とされるマグニチュード 8 クラスの地震が同じ領域で百〜数百年の周期で繰り返されているのに対し，内陸活断層は最大でもマグニチュード 7 クラスの地震がそれぞれ千〜数万年以上の非常に長い周期で繰り返されている．
本稿では，まず第 2 章にて導入として日本の政府機関による繰り返し地震予測の実施内容とその予測手法を解説する．続く第 3 章は，繰り返し地震予測を実践するにあたって考慮すべき
様々な不確定性を説明し，その統計的扱いについて具体例を挙げながら幾つかの研究成果と著者の知見を示していく．第 4 章では，近年盛んになっている中小規模の繰り返し地震の研究を取り上げて，主に地震予測への繋がりの観点から研究動向とその成果を紹介する．
2.

国内政府機関による繰り返し地震予測の現状

2.1 地震予測の取り組みと予測結果の公表
現在，日本では，政府の特別の機関である地震調査研究推進本部（以下，推進本部と呼ぶ）の地
震調査委員会により，国内の主要な活断層に関する地震予測が行われており，予測結果が当該機
関のウェブサイト上に公表されている（http://www.jishin.go.jp/main/p hyoka02.htm）．その予測内容は，将来に大規模な地震を起こすと考えられる活断層について，将来起こりうる地震のマに地震が起こる確率
グニチュード（規模）と，将来の一定期間（30 年間，50 年間，100 年間など）を評価するというものである．マグニチュードの予測は主に，活断層に残された地震の痕跡から判断された，地震によるずれ滑りを起こした活断層の区間の長さに基づいている．また，地震発生確率を評価する期間が複数設けられているのは，0%や 100%付近ではない意味のある確率値が得られる予測期間が活断層により異なるためである．前述したように活断層ごとに地震が繰り返される周期は大きく異なるため，特に海溝型地震では今後数十年間の地震発生確率が90%を超えるものがあるのに対し，活断層地震は将来百年間でも地震発生確率がほぼ 0%となるものが多数存在している．そのほか，予測の元にもなっている活断層の位置，ずれ滑りの向きと速度，過去の地震活動履歴などの情報も一緒に掲載されている．地震調査委員会による予測手法
は推進本部のウェブサイト（http://www.jishin.go.jp/main/choukihyoka/01b/chouki020326.pdf）に詳しく説明されているが，以降ではその統計モデルと地震発生確率の算出方法に焦点を当てて，モデルや手法の背景なども踏まえて概略を解説していきたい．
2.2 繰り返し地震の統計モデル
地震調査委員会の繰り返し地震予測に用いられている統計モデルは，すべて時間軸上の点過活断層で繰り返される地震の点過程モデルとその長期確率予測

85程となっている．点過程とは，確率的に配置される点の数理モデルであり，地震や事故のような事象の生起あるいは人や粒子などの位置を表現するのに用いられる．点過程には平面や時空間などといった多次元空間上のものも存在するが，繰り返し地震は空間的には同じ震源位置（活断層）
でしか起こらないため，活断層ごとに地震の発生時刻を時間軸上の一次元の点過程とみなすことができる．地震調査委員会は活断層から得られる地震活動履歴などの情報量に基づいて，予測に適用する統計モデルを 3 種類用意して使い分けている．ここでは，それらの統計モデル
を単純なものから順に紹介していきたい．

1 ポアソン過程
ある活断層あるいは地域で起こる地震を確率的に予測するとき，最も単純なアイデアは次の例のようなものとなろう．
「過去 1000 年間の間に 10 回地震が起こっている．したがって，平均
して 100 年に一度地震が起きていて，1 年あたりに地震が起こる確率は 1%程度であろう」．こ
の考え方はとても単純ではあるが，現在も多くの地震予測に用いられている．点過程において，
このようなモデルはポアソン過程と呼ばれ，最も単純な形の点過程として知られている．ポア
ソン過程では，将来のある期間に地震が起こる回数はポアソン分布に従い，同じ長さの期間であれば地震が起こる確率はいつも一定となる．
一方，前章で述べた断層で地震が繰り返される物理過程と照らし合わせたとき，上記の考え方にはひとつの疑念が浮かぶ．すなわち，前回の地震が起きた直後と相当期間経過後とで，その活断層が再び地震を起こす確率ははたして同じであろうかという疑問が生まれる．断層の応力がその破壊強度（静止摩擦力）を超えて地震すべりを起こすと，地震すべりで解放されたひず
みの量だけ応力は一旦降下するため，次の地震が起こるまでにはまた同程度の応力が蓄えられなければならない．しかし，ひずみとそれに伴う応力は主にプレート運動に沿ってゆっくりと蓄積されていくため，前回の地震の直後に同じ断層で再び地震が繰り返されるとは考えづらい．
それゆえ，特定の活断層のみを対象とした繰り返し地震予測に対してポアソン過程を適用するのはあまり適切とは言えない．

2 更新過程上記の物理的解釈から，繰り返し地震の発生確率は前回起こった地震からの経過時間に依存して変化すると考えられるため，ポアソン過程から一歩踏み込んだ繰り返し地震の統計モデル
として更新過程と呼ばれる点過程が考えられてきた．活断層では，地震が起こるとそれまでに蓄えられたひずみがある程度リセットされて，そこからプレート運動によるひずみの蓄積が再開するという一連のサイクルによって地震が繰り返されている．したがって，この 1 サイクルにかかる期間すなわち地震活動間隔を確率変数とみなして，それぞれが独立に同一の確率分布に従うと仮定したモデルが更新過程である．
地震活動間隔の確率分布にはこれまで様々な種類のものが考えられ適用されてきた．例えば，
一般に故障確率を求めるのに利用されるワイブル分布（Hagiwara, 1974）や，岩石破壊実験によなどは古くから用いられており，また実際の繰り返し地り導かれた二重指数分布（Mogi, 1962）
震データへの当てはまりのよさの観点から，ガンマ分布（Utsu, 1984）や対数正規分布（Nishenkoなども適用されてきた．ちなみに，地震活動間隔の確率分布を指数分布としand Buland, 1987）
たとき，更新過程はポアソン過程へと帰着される．
地震調査委員会の繰り返し地震予測では，これらの分布と比較検討を行ったうえで，Brownian
Passage Time 分布を更新過程における地震活動間隔の確率分布として採用している．Brownian（以下「BPT 分布」という．
）
は，統計学では逆ガウス分布やワルド分布と呼ば
Passage Time 分布れるものであるが，地震学においては背景となる物理モデル（ブラウン緩和振動過程）との関連（応からこの名前で呼ばれている．Matthews et al.（2002）は，活断層で地震を起こそうとする力
力）の時間的推移を，プレート運動による定常的な応力蓄積に周辺地震活動による応力かく乱を

86

統計数理

第 63 巻

第 1 号 2015

図 1．ブラウン緩和振動過程における応力の時間推移と地震発生時間との関係．
Fig. 1. Temporal evolution of stress and earthquake recurrence times in the Brownian
relaxation oscillator.

加えたブラウン運動と見なして，応力が一定の破壊強度に達したときに地震が起こり，その時点で応力は一定の水準まで下がるとする，一連の物理的プロセスを仮定し，そのときの地震活
動間隔が BPT 分布に従うことから BPT 分布の適用を提案している．このような応力の時間発
展はブラウン緩和振動過程と呼ばれ，数式（確率微分方程式）では以下のように定義される．
(2.1)

dS(t) = λ · dt + σ · dW (t).

ここで，S(t) は時刻 t における応力状態を表す変数である．上式で応力の時間変分 dS(t) は，定常的な応力蓄積に対応するドリフト項 λ · dt と，周辺地震活動等による応力かく乱に対応する（W (t) は標準ブラウン運動）の和として定義されており，応力 S(t) はこの式に
拡散項 σ · dW (t)従うブラウン運動として推移することとなる．なお，ドリフト項の係数 λ はドリフト係数，拡散項の係数の二乗 σ 2 は拡散係数と呼ばれ，それぞれ単位時間あたりの応力増分の平均と分散になっている．さらに，応力状態 S(t) が断層の破壊強度 Sf に達すると地震が発生し，初期状態
S0 まで降下することを仮定すると，応力状態 S(t) の時間推移は図 1 のように上昇と降下を繰り返す．このときの地震活動の時間間隔は Brownian Passage Time と呼ばれ，それぞれが独立にBPT 分布に従う確率変数となる．BPT 分布は以下のような形の確率密度関数をもっている．
μ(x − μ)2(2.2)
f (x|μ, α) =exp−
(x > 0).
2πα2 x32μα2 x
ここで，x は地震活動間隔であり，μ と α はその平均および変動係数（標準偏差／平均）を表す
パラメータとなっている．これらのパラメータは，上述したブラウン緩和振動過程におけるパラメータと以下のとおり対応している．
(2.3)

μ=

(2.4)

α= 

Sf − S 0,λσλ(Sf − S0 ).

上式から分かるように，地震活動間隔の平均 μ は応力の蓄積率 λ と反比例し，変動係数 α は拡散項の係数 σ に比例した関係となっている．

3 時間予測モデル
1 ，
2 のモデルは，活断層の各地震活動間隔が同じ確率分布に従うこと
以上で紹介してきたを仮定しているが，この仮定は必ずしも正しくない．その主な理由のひとつに，地震ごとに断層

活断層で繰り返される地震の点過程モデルとその長期確率予測87

図 2．時間予測モデルにおける断層のすべり量と地震発生時間との関係．
Fig. 2. Fault slip and recurrence times in time-predictable model.

すべりの量が異なるという点がある．地震によるすべり量は，すなわち断層に溜まっていたひずみの解放量であり，応力の降下量とも概ね比例した関係にあると考えられる．それゆえ応力の蓄積率が一定であるとすると，次に応力が破壊強度に達して地震が起こるまでの時間は，前回の地震によるすべり量と比例することになる．時間予測モデル（Shimazaki and Nakata, 1980）
は，この考え方に従い，期待される地震活動間隔が前回の地震のすべり量に比例して変化することを仮定したモデルである．図 2 のように，長期的に見た活断層の単位時間あたりすべり量（すべり速度）を V ，前回の地震によるすべり量 U とおくと，次に起こる地震までの発生時間間
隔 T を T = U/V として予測できるため，‘時間予測’ モデルという名称が付けられている．
2 に述べた BPT 分布の更新過程にこの考え方を取り入れた時間予測モ
地震調査委員会は，
デルを第 3 の統計モデルとして採用した．すなわち，地震活動間隔が従う BPT 分布の平均パラメータ μ を上記のすべり速度 V と最新活動のすべり量 U から μ = U/V とおいて，次の地震
活動時期を予測するのがここでの時間予測モデルである．各地震のすべり量に応じて対応する
各地震活動間隔の分布が異なるという点が更新過程との違いとなる．BPT 分布の物理的背景であるブラウン緩和振動過程との関係で言えば，地震を起こす破壊強度 Sf は毎回一定であるが，地震活動時点の応力降下量 Sf − S0 が断層のすべり量に比例し，ゆえに降下後の応力 S0 が地震毎に異なると仮定したものが時間予測モデルとなる．
地震調査委員会が繰り返し地震の地震発生確率評価に適用しているモデルは，以上で紹介したポアソン過程，BPT 分布更新過程，時間予測モデルの 3 種類となっている．その使い分け方
は，主に過去の地震活動に関する情報量に依存している．南海トラフの海溝型地震は，過去複数回の地震活動についてすべり量および発生年代が得られており，各地震のすべり量とその後の活動間隔が概ね比例していることから，時間予測モデルが採用されている．それ以外で，最新活動の発生年代と平均的な地震活動間隔が判明している内陸活断層や海溝型地震に対しては，更新過程が適用されている．平均的な地震活動間隔しか分からない場合には，前回の地震活動
からの経過時間に依存しないポアソン過程が唯一の選択肢となる．また，海溝沿いの領域の中には，同じ震源域で繰り返し発生している地震がほとんど知られていないことから，領域内で起こる地震をまとめてポアソン過程で予測するところもある．
2.3 モデルパラメータの推定適用すべきモデルが定まれば，次にモデルに利用されているパラメータの推定が問題となる．
時間予測モデルで使われる長期的なすべり速度 V と最新の地震のすべり量 U は，それぞれ断層面を挟んだ地層のずれの大きさと地層に残る地震痕跡から地質学者の判断により推定される．
そのため，ここでは BPT 分布更新過程の適用に必要となる地震活動間隔の平均パラメータ μ と

88

統計数理

第 63 巻

第 1 号 2015

変動係数パラメータ α に焦点を当てて，それぞれがどのように推定されているかを解説する．
まず，地震活動間隔の平均パラメータ μ については 2 通りの推定方法が存在する．ひとつは直接法と呼ばれ，過去の地震活動年代から求まる地震活動間隔の平均値を μ の推定量とする方法である．もうひとつは間接法と呼ばれ，長期的なすべり速度 V と，地震 1 回あたりの平均的なすべり量 Ū が判明している場合に μ = Ū /V を推定量とする方法である．ここで，Ū は時間予測モデルで使われるすべり量とは異なり，過去のすべり量の平均あるいは活断層の地震活動
区間の長さとの経験的関係などにより推定される．地震調査委員会による評価では上記 2 つの手法について，直接法により平均活動間隔が推定できる場合は直接法を採用し，直接法が適用できないあるいは信頼度が低い場合は間接法を採用するという使い分け方をしている．
次に，地震活動間隔の変動係数 α については，活断層ごとに個別に推定するためには，最低でも 2 つ以上の地震活動間隔のデータが必要となる．ところが，現在その条件を満たす活断層
はむしろ少なく，多くの活断層では地震活動間隔のデータが 1 つないし全く得られていない状
況にある．そこで，地震調査委員会は最も多数の地震活動間隔のデータが得られている陸域の4 つの活断層から，ある程度信頼できる α の推定値を得ることを試みた．その結果，4 つの活断層について，変動係数 α をそれぞれ異なるものとして最尤推定量を当てはめるより，共通の変動係数 α をもつと仮定して最尤推定量を当てはめた方が情報量規準 AIC の点で優れていたため，地震調査委員会はここで求めた共通の変動係数の推定値 α = 0.24 を全国の活断層に一律で適用することとした．この方針により，過去の地震活動間隔が得られていない活断層であっても，平均的活動間隔 μ の推定値と最新の地震活動年代さえ判明していれば，変動係数を α = 0.24とした BPT 分布更新過程を適用して予測することができる．
2.4 地震発生確率の算出と伝達以上のプロセスにより将来地震確率の算出に必要な材料が揃い，最後に検討されるのが予測
結果を情報発信する際の「表現方法」となる．地震危険情報を伝達するための表現方法というのは，科学的観点よりむしろ社会的観点から検討されるものであるが，一般社会への地震に対する危機意識の浸透のために大変重要な点である．例えば，各活断層において次に地震が起こる
年代の期待値あるいは予測区間を算出し公開するとする．すると多くの内陸活断層で数百年〜数千年後に次の地震が起こるという予測結果になるため，そうした地域では当面地震の心配はないという感想を持たれかねない．一般社会の関心事は主に我々が生きている間の地震危険であるため，百年以上先の未来の予測が出されても当事者意識は持たれづらいと思われる．地震の危険を間近に感じるためには予測期間が短いほど良いが，日本の陸域の活断層における地震
活動周期は短くとも千年程度に及ぶため，1 ヶ月，1 年という単位で地震確率を見積もっても非常に小さな確率しか得られず，やはり地震への危機意識を薄めてしまう結果となる．これらのことが勘案されて，地震調査委員会では複数の将来スパン（30 年，50 年，100 年など）を用意して，活断層ごとにそれぞれの期間で地震が起こる確率を算出し公開している．将来地震確率は
陸域の活断層では 0〜数十パーセントまで分布しており，より発生頻度の高い海溝型地震はところにより 99%以上となるため，複数のスパンで地震発生確率を与えることでそのいずれかで意味のある予測結果を得ることができよう．
（ここで将来の地震発生確率の算出には，地震活動間隔に適用する確率分布の密度関数 f (x|θ)）にパラメータの推定値 θ̂ を代入したもの f (x|θ̂) を用いる．現は一般にパラメータを θ とおく．
時点で最新の地震活動から x 年が経過しているとき，現時点から将来 Δ 年間における地震発生
確率は，次式により求まる．

活断層で繰り返される地震の点過程モデルとその長期確率予測

 x+Δ
(2.5)

Fx (Δ|θ̂) = x ∞x

f (x|θ̂)dx

f (x|θ̂)dx

89.

右辺の分母と分子はそれぞれ，最新の地震活動から次の地震活動までの期間が x 年以上となる
確率，x 年以上 x + Δ 年未満となる確率であり，それらの比を取ることで最新活動から既に x
年経過しているときに x + Δ 年経過する前に地震が起こるという条件付き確率になっている．
（2.5）式の密度関数を（2.2）式
例えば地震活動間隔の確率分布に BPT 分布を適用する場合には，の f (x|μ, α) とおいて，前節の手法で得られるパラメータ μ, α の推定値を代入することで地震発生確率が求まる．
以上の方法によって，2015 年 1 月 14 日現在では，日本内陸の 110 の活断層帯と日本列島周辺の海域における将来の地震発生確率が地震調査委員会により評価され公表されている．最終
的には，地震調査委員会はこの地震発生確率の予測結果を，マグニチュードの予測結果および強震動評価（予測した地震が発生したときの各地の揺れの大きさの評価）と組み合わせて，「全国
地震動予測地図」（http://www.jishin.go.jp/main/chousa/14 yosokuchizu/index.htm）にまとめている．全国地震動予測地図では，各地で将来 30 年以内に一定以上の震度の揺れに見舞われる確率が日本地図に色分けで示されており，自分の住んでいる地域でどの程度地震の危険がある
かをひと目で確認できるようになっている．このように，統計モデルを用いた地震確率の予測
は，あくまで総合的な地震危険評価の中の一部分に過ぎないが，それでもその根幹を成す重要なパートを担っていることは言うまでもない．
3.

繰り返し地震予測における不確定性とその統計的扱い

前章で解説してきた地震予測手法は，利用するデータが十分かつ確実なものである限りは，信頼できる予測を与えるものと考えられる．しかしながら，過去の大地震に関するデータは，その収集の困難さから，決して十分であるとは言えず，さらに相当の不確かさを抱えている．繰
り返し地震の予測を実践するにあたっては，こうした原因から生ずる様々な「不確定性」
を取り扱わなければならない．地震調査委員会による繰り返し地震予測では，一部の不確定性を予測に反映させているものの，決して統計的な取り扱いがなされているわけではない．統計モデルを適用して地震発生確率を評価している以上，不確定性も統計的枠組みの中で取り入れられるべきである．そこで，本章では，実際の繰り返し地震予測において現れる様々な不確定性を説明し，著者自身の研究成果や知見も交えて，それらを統計的に扱った予測手法を解説していきたい．
3.1 地震活動年代の不確定性
統計モデルの適用に必要となる活断層の過去の地震活動年代を得るための手段として，過去の歴史資料などから地震の発生を読み解く方法と，地質調査により活断層に残っている地震の痕跡を特定する方法とがある．歴史資料から特定される前者の地震は歴史地震と呼ばれるが，これまでの研究の成果により近代以降の主要な歴史地震が数多く把握されており，さらに遡って千年以上前の地震までもが一部特定されている．ただし，震源となった活断層を史料のみから特定するのは難しく，活断層と推定震度分布との位置関係や地盤変動の痕跡などによる裏付けを行う必要がある．一方，後者の地質調査による地震の特定は，活断層の地層に表れている地
震ずれの痕跡を専門家の知見により識別し，地震痕跡を挟む上下の地層内の堆積物（木炭など）に残存している炭素同位体の構成比率から地震活動年代を絞り込むという方法が取られる．地質調査では地下深くまで穴や溝を掘って地層内の地震痕跡や堆積物を探るボーリング調査およびトレンチ調査によって，歴史地震よりはるかに昔の数万年前の地震も特定することができる

統計数理

90
の いねだけ第 63 巻

第 1 号 2015

は ね やま表 1．野 稲 岳 - 万年 山 断層帯の地震確率評価．
Table 1. Earthquake forecast on the Noinedake-Haneyama fault zone.

が，その活動年代は 1 点には定まらず，地震痕跡の前後に堆積された試料の推定年代で挟まれた期間内の地震活動として特定される．活断層の地震予測ではこのような過去の地震活動年代の不確定性の取扱いを避けることができず，地震調査委員会では平均地震活動間隔の推定および
将来地震確率の予測にあたり，以降の例で示していくように，この不確かな活動期間内で取りうる活動間隔の平均値および地震確率の予測値の範囲を求めて，その範囲をもって評価結果として提供している．こうした地震活動年代の不確定性の取扱い方は，ある意味で保守的であると言えるが，複数の予測値を提示することで，読む人に多少の混乱を与えることも懸念される．
の いねだけ は ね やま
まず初めに，表 1 に示された野 稲 岳 - 万年 山断層帯の地震予測を例にとって，地震活動年代の不確定性の取扱い方とそれによる地震確率評価の違いを見ていく．地震調査委員会ではこの
野稲岳 - 万年山断層帯の平均的な地震活動間隔について，断層の隆起量から推定される長期的な（間接法）
ずれの速度と，地震活動区間の長さから推定される地震 1 回あたりのずれの量との比によって，約 4000 年と評価している．一方で，過去の地震活動年代については，最新活動のみ，
その年代を評価しており，地震痕跡の上下層に残る堆積物の推定年代より，3900 年前から 6 世紀までの期間内に活動したものと幅をもって評価している．この地震確率の予測には，平均活動間隔と最新活動年代が判明していることから，最新活動からの経過年数に依存する統計モデ
ルである BPT 分布更新過程が適用される．このとき，最新活動年代を上記の期間内で動かし式により求めると，最新活動年代が 3900 年前のとき最ながら将来 100 年間の地震確率を（2.5）
大値の約 9%を，6 世紀末のとき最小値の約 0.001%をとることから，地震調査委員会は将来 100
年間の地震確率評価を 0.001〜9%と幅をもって与えている．こうした幅をもった予測は，今後の調査で最新活動年代がさらに絞り込まれた場合に，改定された地震確率評価がこの幅を逸脱することはないため，保守的な見積り方であるといえる．しかし，0.001〜9%という幅広な評価値からでは地震危険の高さが判断しづらいという弊害もある．また，仮に最新活動年代が全く不明であった場合にはポアソン過程が適用され，平均活動年代の推定値 4000 年から将来 100 年
間の地震確率評価は 2.5%と 1 点で与えられるため，最新活動年代の情報が増えることで確率評
価の不確かさが拡がっているという矛盾もある．
こうした地震活動年代の不確定性に対しては，その不確定性を確率分布として表現すること
により統計的な扱いが可能となり，一貫した予測を行うことができる．Ogata（1999）は，ベイズ統計の枠組みで不確かな地震活動年代に事前分布を与えて，パラメータの推定および将来の地震予測を行っている．一般に，図 3 に示されるように，観測期間 [S, T ] 内に起こった地震活
動年代 t1 , . . . , tn (S ≤ t1 < · · · < tn ≤ T ) が正確に特定できず，確率密度関数 ψ(t1 , . . . , tn ) を持つ事前分布に従うとしたき，将来 Δ 年間の地震確率のベイズ予測を次式により与えることができる．
FT −tn (Δ|θ̂)L(θ̂|t1 , . . . , tn )ψ(t1 , . . . , tn )dt1 · · · dtn(3.1)
Fψ (Δ|θ̂) =.
L(θ̂|t1 , . . . , tn )ψ(t1 , . . . , tn )dt1 · · · dtn

活断層で繰り返される地震の点過程モデルとその長期確率予測91

図 3．繰り返し時間予測における観測データと予測期間．なお，地震発生年代 t1 , . . . , tn は独立ではなく同時事前分布に従う．
Fig. 3. Observed data and time horizon for recurrent earthquake forecast. Occurrence
dates t1 , . . . , tn are not independent and jointly distributed.

ここで，L(θ|t1 , . . . , tn ) はパラメータ θ の尤度関数であり，以下の式で与えられる．
 ∞
n−1

(3.2)
L(θ|t1 , . . . , tn ) = f FRT (t1 − S|θ)f (ti+1 − ti |θ)f (t − tn |θ)dt.
i=1

T

この右辺の f FRT (x|θ) は，地震発生時点ではないある時点から次に地震が起こるまでの経過の確率密度関数で，地震活動間隔の累積分布関数
時間（Forward recurrence time と呼ばれる）
x
F (x|θ) = 0 f (y|θ)dy を用いて次式で表される（Daley and Vere-Jones, 1988）(3.3)

1 − F (x|θ)
f FRT (x|θ) =  ∞.
1 − F (x|θ)dx0

ただし，断層の年代推定は最も古い地震の発生地点で終わることが多く，（3.2）式にて期間の始点 S
FRTが最初の地震活動年代 t1 である場合（S = t 1 ）には f(t1 − S|θ) が不要となり L(θ|t1 , . . . , tn ) =∞n−1（3.2）式の尤度関
i=1 f (ti+1 − ti |θ) T f (t − tn |θ)dt となる．分布パラメータ θ の推定値 θ̂ は，数を最大化する値（最尤推定量）などとして求められる．（3.1）式の L(θ̂|t1 , . . . , tn )ψ(t1 , . . . , tn ) はある地震活動年代の組 (t1 , . . . , tn ) に対する尤度と事前分布の積であり，その地震活動年代の起式は，それぞれの地震活動年代 (t1 , . . . , tn ) から求まる将来の地
こりやすさを表している．（3.1）（最新の地震発生からの経過年数は T − tn である）を，その起こりやすさで
震確率 FT −tn (Δ|θ̂)
加重平均をとったものと解釈される．
地震活動年代の確率分布は，最も単純には活動した可能性のある期間内の一様分布を適用することが考えられる．地質調査によって地震活動年代に関して得られる情報は，地震痕跡の上下層に残る堆積物の推定年代すなわち地震活動年代の上限と下限であるため，地震活動年代の確率分布としてこの上限と下限に挟まれた区間内の一様分布を用いることは自然な考え方といえる．ただし，堆積物の推定年代も一定の誤差を伴っているため，厳密にはこの上限と下限は確率的に変動し，地震活動年代の分布はその変動分だけ一様分布から歪むことになる（Biasi et．多数の堆積物の推定年代を利用する場合，地震活動年代の分布は複雑になるため，
al., 2002）
直接的に分布関数を算出する代わりに，堆積物の年代と地震活動年代の分布標本をマルコフ連鎖．最新活
モンテカルロ法によって抽出する方法なども提案されている（Hilley and Young, 2008）
動年代の分布に一様分布を仮定した場合の地震確率の算出式については，地震調査委員会の報告書（http://www.jishin.go.jp/main/choukihyoka/01b/chouki020326.pdf）でも示されている．表
1 の野稲岳 - 万年山断層帯のケースに戻って考えると，判明している過去の地震活動が最新活動のみであり，最新活動年代の分布が x1 = 1400 年前（6 世紀末）から x2 = 3900 年前の範囲におけ

92

統計数理

第 63 巻

第 1 号 2015

表 2．十日町断層帯西部の地震確率評価．
Table 2. Earthquake forecast on the West Tokamachi fault zone.

る一様分布であると仮定したとき，（3.1）
式の地震確率はより簡潔な次の式へと帰着される．
 x1 +Δ
1 − F (x|θ̂)dx
x +Δ
case1
.
(3.4)
F(x1 ,x2 ) (Δ|θ̂) = 1 − 2x1
1 − F (x|θ̂)dxx2
この式に従って野稲岳 - 万年山断層帯の将来 100 年間の地震確率を求めると 2.5%となり，最新活動年代を使わずにポアソン過程から求めた結果とほぼ同じになる．
もう一つ，特殊なケースとして表 2 の十日町断層帯西部の地震予測を考える．この活断層では，断層のずれ情報（間接法）により平均活動間隔は 3300 年と推定されるものの，最新活動が特定されておらず，少なくとも 3100 年前から現在まで地震活動がなかったことだけが判明している．よって最新活動年代は 3100 年前以前であるとしか分からない状況にある．この断層に対して地震調査委員会は BPT 分布更新過程を当てはめ，将来 100 年間の地震確率を最新活動年代がちょうど 3100 年前であるときの約 10%かそれ以上であるという評価を下している．最新活動
年代の下限がないために，地震確率の上限は与えられないという考え方を採っているが，このようなケースでも更新過程の理論を用いれば地震確率を 1 点に予測することができる．
十日町断層帯西部のケースは，3100 年前からの Forward recurrence time が次の地震発生時刻となる．ただし，現在までの x = 3100 年間は地震が起こっていないため，将来 Δ 年間の地震（2.5）
式に当てはめて次のように求まる．
確率は（3.3）式の密度関数を
 x+Δ
 x+Δ FRTf(x|θ̂)dx
1 − F (x|θ̂)dx= x ∞.
(3.5)
Fxcase2 (Δ|θ̂) = x ∞FRT (x|θ̂)dx
f1− F (x|θ̂)dx
x
x
case1
なお，
（3.4）式の確率算出式において x2 を無限大とする極限 limx2 →∞ F(x,x(Δ|θ̂) をとることで
2)（3.5）
式を用いて十日町断層西部における将来 100 年
も（3.5）式を導出することができる．この間の地震確率を算出すると，表 2 にあるように 13%が得られる．これは地震調査委員会による評価値の下限 10%よりも大きいものの，下限とそれほど差のない数値ともいえる．なお，ポア
ソン過程を適用した場合の地震確率は 3%となるため，最新活動年代を用いた BPT 分布を適用することで確率評価は大きく上昇している．最新活動から既に平均活動間隔 3300 年に近い期間を経過しているため，BPT 分布によると次の地震が十分に起こりうる時期にあるものと考えら
れる．

3.2 地震活動可能性の不確定性
前節では地震の活動年代に関する不確定性について議論してきたが，より根本的に過去の地震の活動可能性が不確かであるケースも多く存在する．地層に残る地震痕跡が曖昧であったり，地層が残っていないなどの理由で，一部の年代で地震活動の有無あるいは回数が特定できない
ことがある．このようなケースにおいて地震調査委員会は，前節と同じ考え方で可能性のある
過去の地震活動パターンの範疇で取り得る平均活動間隔および地震発生確率の範囲を示す方法

活断層で繰り返される地震の点過程モデルとその長期確率予測

93

や，あるいは各地震活動パターンで個別に地震発生確率を評価してそれぞれの結果を提示する方法を採っている．
過去の地震活動パターンが複数考えうる場合の統計的な取扱いは一概には言いづらいが，例えば想定される全パターンの尤度を合算した上で推定，予測を行う方法や，専門家の知見により定めた各パターンの生起確率を重みとして各パターンにおける地震確率評価の加重平均をとる方法などが考えられる．本稿では，こうした考え方の相違の余地が少ない一つの特殊なケースを例に，地震活動可能性の不確定性の統計的な取扱いを示したい．表 3 にあるサロベツ断層
帯では，5100 年前から 4500 年前の間に地震があったことが判明しているものの，それ以降の地
震活動の有無が判断できないという状況にある．平均活動間隔は断層のずれ情報（間接法）より4000 年から 8000 年の間と推定されている．ここでも地震調査委員会は BPT 分布更新過程を当てはめ，最新活動年代が最も古い 5100 年前でかつ平均活動間隔が 4000 年である場合の地震確率約 10%を上限として，下限は最新活動年代を直近にとることでいくらでも小さくなるために示さないという方針をとっている．このようなケースにおいても，点過程の理論から不確定性を取り扱うことができる．x 年前に地震活動があって以降，地震活動の有無が不明である場合の，将来 Δ 年間の地震確率は次式により算出される．
∞  x+Δ(3.6)
Fxcase3 (Δ|θ̂) =f (k) (y|θ̂){1 − F (x + Δ − y|θ̂)}dy.
k=1

x

(k)

ただし，f (y|θ) はある地震から k 回目の地震が起こるまでの時間間隔で，帰納的に f (1) (y|θ) =y（3.6）式右辺のf (y|θ), f (k+1) (y|θ) = 0 f (k) (z|θ)f (k+1) (y − z|θ)dz, k = 1, 2, . . . と定義される．
各項は，x 年前から数えて k 回目の地震が将来 Δ 年間に起こり，かつその次（k+1 回目）の地
震は Δ 年後以降に起こる確率となっており，総和をとることで何回目かの地震が将来 Δ 年間に起こる確率を表している．地震活動間隔に BPT 分布を適用する場合，BPT 分布の再生√性から f (k) (y|μ̂, α̂) = f (y|μ̂k, α̂/ k) となり，また平均 μ に比べて予測期間 Δ が十分短いとき
F (x + Δ − y|μ̂, α̂) はほぼ 0 となるため，（3.6）
式はより簡単に次式で近似計算できる．
∞  x+Δ√(3.7)
Fxcase3 (Δ|μ̂, α̂) ≈f (y|μ̂k, α̂/ k)dy.
k=1x

これを表 3 のサロベツ断層帯に当てはめて計算すると，平均活動間隔を 4000 年としたとき 1.8%，
8000 年としたとき 0.1%となり，いずれも地震調査委員会の評価した地震確率の上限 10%よりはるかに小さく，むしろポアソン過程適用時の確率に近い数値が得られた．これは，直近の比較的長い期間で地震活動が不明であるために，過去の地震活動に依存しないポアソン過程の結果に近づいたものと考えられる．
3.3 過去の地震活動データの不足によるパラメータ推定の不確定性
活断層地震の予測において最も単純かつ深刻な問題はデータの不足であろう．とりわけ日本
表 3．サロベツ断層帯の地震確率評価．
Table 3. Earthquake forecast on the Sarobetsu fault zone.

94

統計数理

第 63 巻

第 1 号 2015の内陸活断層の繰り返し地震はおおよそ千年〜一万年以上の間隔を置いて発生しているので，1
つの地震活動間隔を得るのに必要な過去 2 回分の地震活動年代データを得るのですら容易ではない．地表に地層が表れていない場合，地下深くまで穴を掘るボーリング調査やトレンチ調査が必要となり，調査に必要な予算や人手の問題もある．過去の地震活動データの取得には多大な費用と労力が費やされており，今なおさらなるデータ収集のための調査が行われている．
データが不足することで問題となるのが，統計モデルにおけるパラメータの推定誤差である．
ポアソン過程や BPT 分布更新過程を適用する場合，平均活動間隔のパラメータが必要となるが，地震調査委員会は過去の地震活動間隔あるいは断層のずれの量・速度を用いてそれを推定していることを前章で既に述べた．ここで，地質調査に基づいて特定された過去の地震活動年
代には幅があるため，その場合，地震調査委員会はその幅の範囲内の活動間隔の最大値・最小値を用いて，平均活動間隔も幅をもって推定している．しかし，ほとんどの活断層ではせいぜい 1，2 個しか過去の活動間隔データが得られていないため，そもそも過去の活動間隔と真の平均活動間隔パラメータとの間に大きな乖離がある可能性があり，それを統計的に考慮する必要がある．また，間接法により平均活動間隔を断層のずれ情報から求める場合にもやはり推定誤差が発生する．図 4 では間接法で推定された平均活動間隔と過去の活動間隔およびその平均値を比較しているが，地震調査委員会の推定した変動係数から予測される活動間隔の範囲（図 4 の点線で挟まれた範囲）
をいくつかのデータが大きく逸脱しているのが分かる．このようなデータと真のパラメータとの間の誤差を織り込まなければ，地震調査委員会のような幅をもたせた推定であっても統計的な区間推定とは言えない．
さらに，BPT 分布更新過程や時間予測モデルを適用する場合に必要となる活動間隔のばらつき（変動係数）のパラメータ α についても，地震調査委員会の推定値 α = 0.24 を全国一律に適用することには異論も出ている．前述したとおり，地震調査委員会の α の推定値は当時の 4 つの陸域活断層の活動履歴から最尤推定されたものであるが，その後の地質調査および検討により
当該活断層の活動履歴は更新され当時と大きく異なるものとなった．また，それ以外の活断層の調査も進展するに連れて，当該推定値では当てはまりが悪い活断層がいくつも存在することが明らかになっている．そうした変化を踏まえて，Kumamoto and Hamada（2005）では国内 23の陸域活断層の活動履歴データから BPT 分布の変動係数を α = 0.49 と推定し直しており，これは地震調査委員会の推定値と大きく異なる値となった．活断層によっても変動係数 α の最尤推定値にも大きなばらつきがあり，本来は活断層ごとに異なる変動係数をもつと考えるべきであろう．
以上のようにモデルパラメータの推定誤差が大きい状態では，パラメータの推定値をモデルに代入する予測手法はときに偏った結果をもたらす．なぜなら，代入したパラメータの誤差が大きいと，予測の最終目的である将来の地震確率評価が大きくぶれる場合があるからである．
例えば，BPT 分布更新過程で考えた場合，同じ平均千年周期の活断層でも，変動係数 α を地震調査委員会の 0.24 でおくのと，Kumamoto and Hamada（2005）の推定値 0.49 でおくのとでは，
図 5 のように将来の地震発生確率の高低が全く変わってしまうこともある．図 5 において前回の地震から 400 年経過した時点の将来 100 年間の地震確率は，α = 0.24 とすると 0.2%ほどしかないのに対して，α = 0.49 とすると 7%にまで上昇する．一方，1600 年経過時点では α = 0.24とすると約 50%となるが，α = 0.49 とすると約 20%と大きく下がる．このように変動係数 α の推定値は，地震危険の認識を変えてしまうほどの影響を持っていることが分かる．
このような予測の偏りを無くすには，パラメータの推定誤差を事後分布として予測に織り込んだベイズ予測が有効である．ベイズ予測では，地震活動間隔の分布パラメータに対して，一般的な活断層全体の母集団分布（事前分布）を仮定して，分布パラメータも確率変数として扱う．
ここでパラメータを一般化して θ とおき，事前分布の確率密度関数を π(θ) とおくと，過去の地活断層で繰り返される地震の点過程モデルとその長期確率予測

95

図 4．各断層の間接法（断層ずれ）から推定した平均活動間隔と過去の活動間隔との対比（Nomura．灰色の縦線，横線は地震調査委員会による推定値の幅を示す．点線は間接
et al., 2011）
法の平均活動間隔と変動係数 α = 0.24 に対する上側 2.5%点および下側 2.5%点を示す．
Fig. 4. Recurrence time intervals on a fault against the expected intervals from theestimated ratio of slip per event to a fault’s long-term slip rate. Gray bars aredata given in range by the Earthquake Research Committee of Japan. Dottedlines show 2.5% and 97.5% percentiles of the BPT prediction intervals from theexpected time interval assuming α = 0.24.

震活動年代 t1 , . . . , tn が与えられたときの事後分布の密度関数 π(θ|t1 , . . . , tn ) は（3.2）式の尤度関
数 L(θ|t1 , . . . , tn ) を用いて次のように求まる．
(3.8)π(θ|t1 , . . . , tn ) = L(θ|t1 , . . . , tn )π(θ).
L(θ|t1 , . . . , tn )π(θ)dθ

そして，将来の地震発生確率は，（3.1）
式中の尤度関数をこの事後分布の密度関数に置き換えてパラメータ θ についても積分することで次式のように求められる．

FT −tn (Δ|θ)π(θ|t1 , . . . , tn )ψ(t1 , . . . , tn )dt1 · · · dtn dθ

.
(3.9)
Fπ,ψ (Δ|θ) =π(θ|t1 , . . . , tn )ψ(t1 , . . . , tn )dt1 · · · dtn dθ
もし設定した事前分布が活断層母集団のパラメータ分布を上手く近似できていれば上記の手法
によって，推定値を代入した予測よりも推定誤差に対してロバストな予測を行うことができる．
パラメータの事前分布には，全てのパラメータ値をある意味で平等に取り扱う無情報な事前分布を用いる手もあるが，データが少ない以上は事前分布の段階からパラメータを合理的に絞りでは，多数の活断層地震データを集約して，そこか込める方が望ましい．Nomura et al.（2011）
ら得られるパラメータ情報（尤度）
を用いて活断層全体のパラメータ事前分布を構成している．
さらに，平均活動間隔 μ の事前分布に関しては，断層ずれから導出される間接法の平均活動間

96

統計数理

第 63 巻

第 1 号 2015

図 5．変動係数 α の異なる BPT 分布における経過年数別の将来 100 年間地震確率の比較．平
均 μ は共通で 1000 年としている．将来 100 年間の地震確率 F (x, x + 100|μ, α) の計（2.2）式を参照．
算式は（2.5）式およびFig. 5. Comparison of earthquake forecasts in the next 100 years between the diﬀerent coeﬃcients of variation α for the BPT distribution. Probability forecasts
F (x, x + 100|μ, α) are calculated along the formulas in（2.2）and（2.5）.

隔 T = U/V（U は地震 1 回によるずれの量，V は 1 年あたりのずれ速度）を事前分布に組み込む
ことによって，データの不足を補うことが試みられている．Nomura et al.（2011）による地震発
生確率のベイズ予測では，主に活動間隔の変動係数 α の推定値の差を原因として，いくつかの活断層で地震調査委員会と大きく異なる評価結果が得られている．
（A）に示した
なお，活動間隔の変動係数（ばらつき度合い）のパラメータ α については，図 6
α の推定値の分布に見られるように，ある程度の地域的な傾向が存在することが Nomura et al.
（2011）などにより指摘されている．この地域的な傾向は，海域の海溝型地震と陸域の内陸活断層地震でまず大きく分かれており，海溝型地震の変動係数 α は全般的に低い傾向にある．さらに，陸域の中でも，近畿地方〜中部地方にかけての活断層で変動係数 α は特に高い傾向にある．
（B）にある活断層分布の密集度と関連しており，活断層が密にこの陸域における地域差は，図 6分布しているところでは α が大きい傾向にあることが分かる．これは，前述した BPT 分布の物理的背景であるブラウン緩和振動過程において，（2.4）式に示したようにパラメータ α が周辺
地震による応力かく乱の大きさ σ に依存しており，活断層が密な地域ほどそのかく乱が大きくなるためと説明付けることができる．このような考えを踏まえると，将来的には地域別あるいは断層の空間的特性に基づいて α の事前分布を与えることが望ましいであろう．このように事前分布の情報量を高めることで，特に地震活動間隔データの少ない活断層における予測精度を大いに高めることが期待される．
以上のように，データやパラメータの不確かさを織り込んで予測を立てることにより，デー
タやパラメータに大きな誤差があった場合にもある程度安定した予測性能が発揮できる．なお，
Nomura et al.（2011）では将来の地震確率について 1 点での予測だけでなく，事後分布に基づいた幅のある予測区間も与えている．地震発生確率に不確かさによる変動幅も与えることで，現状の予測値がどの程度の不確かさを抱えていて，さらなる調査により不確かさを解消した場合に地震確率がどの程度変化しうるかを示すことができる．
4.

小規模な繰り返し地震の予測と解析

4.1 地震活動の変化を踏まえた予測
ここまでの話では，社会的に被害をもたらす大規模な繰り返し地震を取り扱ってきたが，より小さなスケールの繰り返し地震も当然ながら存在する．小規模な地震は震源断層が正確に定活断層で繰り返される地震の点過程モデルとその長期確率予測

97

（A）

（B）

図 6．
（A）断層ごとに示された BPT 分布の変動係数 α のベイズ推定量（Nomura et al., 2011）．
（活断層研究会, 1991）．
（B）日本内陸の活断層分布図
Fig. 6.（A）Bayesian estimates of the coeﬃcient of variation α in the BPT distribution
for the respective active faults（Nomura et al., 2011）.（B）Map of active faults
on the inland of Japan.

まらないため，繰り返し地震であることは自明には分からないが，地震波形の相似性などを利．地用して同じ断層で繰り返されている地震系列を特定することができる（Nadeau et al., 1995）
震の発生頻度は，マグニチュードが小さくなるに従い指数関数的に増加するグーテンベルグ・に従うため，小規模な繰り返し地震は大規模なもの
リヒター則（Gutenberg and Richter, 1944）に比べてはるかに多く観測することができ，さらに，大規模なものよりはるかに短い周期で繰

98

統計数理

第 63 巻

第 1 号 2015

り返されている．それゆえ，小繰り返し地震は短期の予測でも高い発生確率が得られ，短期間で予測と実績とを比較して予測手法を検証することが可能である．予測手法の有効性が小繰り返し地震について実証されれば，大地震の予測への活用も期待できるため，近年では小繰り返
し地震の予測実験も行われている（Okada et al., 2011）．しかし一方で，小規模の地震予測には
大地震とは異なる課題がある．それは，広域的な地震活動の変化などに応じて，小繰り返し地震の周期あるいは頻度が変化するという点である．小規模な地震活動は，周辺の大地震による余震活動や，非地震性すべりと呼ばれる一時的なプレート運動の加速などの影響を受けてしば
しば変化している．これまで紹介してきた統計モデルはいずれもプレート運動による応力の蓄積が一定の速度で進んでいくことを前提とした定常なモデルであったが，上記の現象によりプ
レート運動が加速すると応力の蓄積速度も加速されてしまうため，小規模な繰り返し地震デー
タに対する定常なモデルの当てはまりはときに非常に悪くなる．
こうした繰り返し地震活動の変化に対応するため，Nomura et al.（2014）は 2 章で紹介した
BPT 分布の更新過程を時間的に非定常なモデルへと拡張して小繰り返し地震に適用している．
（2.1）式の応力推
応力蓄積率の変化を時間関数 ν(t) で表現すると，BPT 分布更新過程における移は次式の形へと変わる．
(4.1)
dS(t) = λ · ν(t)dt + σ · ν(t)dW (t).
この応力変分の式では，ドリフト係数が λ · ν(t)，拡散係数（拡散項の係数の二乗）が σ 2 · ν(t) と，それぞれ蓄積率関数 ν(t) に比例している．これは，拡散項は周辺地震活動による独立増分を表すため，周辺地震の発生数が蓄積率関数 ν(t) に比例して変化すると，分散の加法性から拡散係数も同様に変化するという考えに基づいている．このようなモデルでは，繰り返し地震の発生間隔はそのまま BPT 分布には従わなくなるが，蓄積率 ν(t) に応じた時間進行になるよう時間
t
変換 t = 0 ν(s)ds を施すことにより，時間変換後の発生間隔が BPT 分布に従うようになる．
このように，時間変換によって通常の更新過程に帰着される非定常な点過程は Trend renewal（Lindqvist et al., 2003）と呼ばれる．
process
将来の応力蓄積率の変化が予測できれば，上記のモデルを用いて小繰り返し地震を予測することができる．平常時は応力蓄積率があまり変化しないが，例えば大地震後であれば余震活動の時間減衰則である大森・宇津の法則（宇津, 1957; Utsu, 1961）から将来の応力蓄積率の予測を立てて，余震の減衰を考慮した繰り返し地震の予測ができる．Nomura et al.（2014）は，図 7（A）のアメリカ西部に位置する San Andreas 断層の小繰り返し地震群を解析するにあたり，2004
年 Parkﬁeld 地震（M6.0）
後の繰り返し地震データに大森・宇津の法則を当てはめ，推定された（B）は 2004 年 Parkﬁeld 地震
時間減衰を広域的な応力蓄積率の時間変化として用いている．図 7（M6.0）後の繰り返し地震について累計発生数と系列別の発生時刻を示しているが，Parkﬁeld 地
震直後に集中的に地震が起こっており，繰り返し地震の発生間隔も本震からの時間経過に従っ（C）は，大森・宇津の法則で推定した余震の時間減衰て大きく変化している．これに対して図 7（C）の時間変式を用いて地震発生頻度が平準化されるよう時間変換を施している．すると，図 7換後の繰り返し地震は系列ごとにほぼ一定の周期で発生している様子が見て取れ，時間変換を加えれば更新過程が上手く当てはまることを示している．
4.2 繰り返し地震によるプレート運動のモニタリング
以上の研究は繰り返し地震自体の予測を目的としてきたが，そうではなく，小規模な繰り返し
地震活動の変化から，大地震の前兆を見出して地震予測に繋げる研究もある．例えば，前述した非地震性すべりはプレート境界のある領域で数日〜数年の間，プレート運動が継続して加速する現象であるが，この現象は大地震の発生前にその震源域付近でよく見られることが知られ

活断層で繰り返される地震の点過程モデルとその長期確率予測

図 7．
（A）アメリカ San Andreas 断層 Parkﬁeld 地区の地図．点が各繰り返し地震系列の震源，星（M6.0）の震源を示す．（B）2004 年 Parkﬁeld 本震後の繰り返し地
が 2004 年 Parkﬁeld 地震
震データ．上図の横線は同系列の繰り返し地震を示し，点の位置と大きさは各繰り返し地震の発生時刻およびマグニチュードを表す．下図は繰り返し地震の累計発生数の推移．点線は累計発生数に大森・宇津の法則を当てはめた曲線を示す．（C）同データの時間スケールを地震発生頻度が平準化されるよう変換した図．
Fig. 7.（A）Map of the Parkﬁeld segment at the San Andreas Fault. Locations of the repeating sequences are shown by black circles. Also shown is the epicenter of the 2004
Parkﬁeld earthquake（M6.0）.（B）Occurrence data of repeating earthquakes after the2004 Parkﬁeld main shock. Upper panel shows their occurrence times and magnitudes.
Earthquakes of each repeating sequence are arranged in each horizontal line. Lower
panel shows the cumulative number of occurrences in all the repeating earthquakes.
Dotted line is its ﬁtted curve by the Omori-Utsu law.（C）Panels show the same dataas（B）in the transformed time, where the occurrence rate is considered to be constant.

99

統計数理

100

第 63 巻

第 1 号 2015ている．このような現象を解析するのに，従来から GPS（Global Positioning System）観測デー
タによる地殻変動の解析が主要な役割を担ってきた．しかし，GPS 観測点は主に地上に分布していることから，海底下深くのプレート境界における地殻変動の推定精度には限界がある．そこで，プレート境界上で起こっている小繰り返し地震活動の変化を解析して，プレート運動の変動を捉える手法が Nadeau and McEvilly（1999）によって提案された．繰り返し地震の周期の変化は，震源近傍のプレート運動の速度変化を反映しているため，GPS 観測網の届かない海域
あるいは地下深い領域のプレート運動を直接的に推定することが可能となる．
Nomura et al.（2014）はプレート境界上に分布する多数の繰り返し地震系列から，プレート運動の時空間的変動を推定する時空間更新過程モデルを提案している．これは，前節の応力蓄積率の時間関数 ν(t) に位置座標 (x, y) を加えた時空間関数 ν(x, y, t) へと拡張させて，各震源位置(xj , yj ), j = 1, . . . , m における繰り返し地震系列に応力蓄積率の時間推移を νj (t) = ν(xj , yj , t)とした前節の非定常更新過程モデルを当てはめるというものである．応力蓄積率の時空間関数
ν(x, y, t) は，3 次元 B-スプライン関数により表現し，平滑化制約を課した罰則付き尤度関数の最大化により推定することで，繰り返し地震間の時間的・空間的な合間を連続的で滑らかな関数で補間して，蓄積率の時空間的推移を可視化することができる．Nomura et al.（2014）は，こ（A）に示された米国 San Andreas 断層 Parkﬁeld 地区における繰り返し地のモデルを用いて図 7
震群を解析した結果，1993 年頃からの非地震性すべりや 2004 年 Parkﬁeld 地震（M6.0）の余震活動について，時空間的に推移する様子を捉えている．
また，応力の蓄積率は基本的にプレート運動のすべり速度と比例するため，地震のマグニ
チュードとすべり量との関係の経験式（Nadeau and Johnson, 1998）を利用して，小繰り返し地
震からプレート運動の速度分布を推定することもできる．2.2 節で紹介した時間予測モデルの考え方に拠れば，各繰り返し地震のすべり量 U と，そこから次の地震までの時間間隔 T から，その期間の平均的なすべり速度が V = U/T と求まる．この時間予測モデルを前述した更新過程と同様にして時空間モデルへと発展させることで，応力蓄積率の代わりにプレートすべり速度の時空間分布を推定することが可能となる．野村 他（2013a, 2013b）はこの時空間モデルを，東北日本太平洋沖のプレート境界で観測された 1993 年から 2011 年までの小繰り返し地震に適用して，当該地域におけるプレートすべり速度の時空間的推移を推定している．この地域では解析期間中にいくつかのマグニチュード 7 以上のプレート境界地震（1994 年 12 月三陸はるか沖地（M8.0），2004 年 11 月釧路沖地震（M7.1），2008 年 5 月茨城震（M7.6），2003 年 9 月十勝沖地震が発生しており，これらの大地震の前後では余震活動等によってすべり
県沖地震（M7.0）など）
速度が広範囲で不連続に変化するため，ここでは各大地震の起こる前と後とで期間を分けて解析することにより大地震前後におけるプレートすべり速度の様々な特徴的変化を捉えている．
（2013b）の結果のうち，2008 年以降のすべり速度の推移を示している．茨城
図 8 では，野村 他（非地
県沖から福島県沖にかけて，およそ 2008 年から 2009 年までの期間，広域的なすべり加速
震性すべり）が起こっている．また，2010 年後半から岩手県沖にて始まったすべり加速が，2011
年東北地方太平洋沖地震の震源域に向けて徐々に南下していく様子も確認できる．これらはいずれも東北地方太平洋沖地震の震源域の周辺で発生したすべり加速であり，東北地方太平洋沖地震の一つの引き金あるいは前駆的現象と捉えることができる．こうした解釈は，大地震が起こってから後付けで考えられたものではあるが，経験を蓄積していくことで将来的にすべり速度の変動から大地震の可能性が議論できるようになると期待される．
5.

終わりに

以上，本稿では主に日本における繰り返し地震の予測に焦点を当てて，予測手法の解説と現

活断層で繰り返される地震の点過程モデルとその長期確率予測

図 8．2008 年から 2011 年東北太平洋沖地震までの 3 ヶ月ごとの時空間すべり速度分布のス．各パネルの黒点は各時点の前後に発生した小繰り返
ナップショット（野村 他, 2013b）
し地震の震央を，最後のパネルの青星は 2011 年東北地方太平洋沖地震の震央を示す．
Fig. 8. Snapshots of space-time loading rate at every three months from 2008 till the2011 Tohoku earthquake. Black dots in each panel indicate the locations of small
repeating earthquakes occurred around each date. Blue star in the last panel
indicates the epicenter of the 2011 Tohoku earthquake.

101

統計数理

102

第 63 巻

第 1 号 2015

状の課題の説明，最近の研究成果や著者の知見の紹介を行ってきた．日本の陸域で起こる大地
震はほぼ全て既存の活断層を震源とする繰り返し地震であるため，繰り返し地震の予測手法は
大地震予測にとって重要なパートを占めている．繰り返し地震予測には主にデータに起因する課題が付きまとうものの，3 章に述べたようにデータやパラメータ推定の不確定性を取り入れたベイズ統計モデルを適用することで，偏りの小さい信頼できる地震確率評価を提供することができる．
一方で，大きな繰り返し地震の予測ではかなり長期のスパンを取らないと大した地震確率が得られず，大地震の発生前に十分に高い確率で短期予測できるようになるにはまだ相当のブレー
クスルーが必要である．保険業や建設業などの業種では主に長期の地震予測が利用されるが，防災や減災の観点ではやはり短期の地震予測が必要となる．その点では，最後に紹介した小繰り返し地震によるプレート運動のモニタリングは大地震の前兆を短中期的に捉えられる可能性をもっており，地震予測の確率的な効率（確率利得）を今後大きく向上させることが期待される．
謝辞

本論文の執筆にあたって，査読者より丁寧でかつ有益な助言をいただき，心より感謝の意を表す．本研究の一部は，文部科学省科学研究費補助金（課題番号 26870193）の助成を受けたもの
である．

参

考

文

献

Biasi, G. P., Weldon, R. J. I., Fumal, T. E. and Seitz, G. G.（2002）. Paleoseismic event dating andthe conditional probability of large earthquakes on the southern San Andreas fault, California,
Bulletin of the Seismological Society of America, 92, 2761–2781.
Daley, D. and Vere-Jones, D.（1988）. An Introduction to the Theory of Point Processes, Springer-Verlag,
New York.
Gutenberg, B. and Richter, C. F.（1944）. Frequency of earthquakes in California, Bulletin of the Seismological Society of America, 34, 185–188.
Hagiwara, Y.（1974）. Probability of earthquake occurrence as obtained from a Weibull distribution
analysis of crustal strain, Tectonophysics, 23, 313–318.
Hilley, G. E. and Young, J. J.（2008）. Deducing paleoearthquake timing and recurrence from paleoseismic data, part 1: Evaluation of new Baysian Markov-chain Monte Carlo simulation methods
applied to excavations with continuous peat growth, Bulletin of the Seismological Society ofAmerica, 98, 383–406.
活断層研究会（1991）.『新編 日本の活断層—分布図と資料』, p.437, 東京大学出版会, 東京.
Kumamoto, T. and Hamada, Y.（2005）. Examination of aperiodicity parameters for the Brownian
Passage Time model using intraplate paleoearthquake data in Japan, 活断層研究, 25, 9–22.
Lindqvist, B. H., Elvebakk, G. and Heggland, K.（2003）. The trend-renewal process for statisticalanalysis of repairable systems, Technometrics, 45, 31–44.
Matthews, M. V., Ellsworth, W. L. and Reasenberg, P. A.（2002）. A Brownian model for recurrent
earthquakes, Bulletin of the Seismological Society of America, 92, 2233–2250.
Mogi, K.（1962）. Study of elastic shocks caused by the fracture of heterogeneous materials and its
relations to earthquake phenomena, Bulletin of Earthquake Research Institute, University ofTokyo, 40, 125–173.
Nadeau, R. M. and Johnson, L. R.（1998）. Seismological studies at Parkﬁeld VI: Moment release rates

活断層で繰り返される地震の点過程モデルとその長期確率予測103

and estimates of source parameters for small repeating earthquakes, Bulletin of the Seismological Society of America, 88, 790–814.
Nadeau, R. M. and McEvilly, T. V.（1999）. Fault slip rates at depth from recurrence intervals of
repeating microearthquakes, Science, 285, 718–721.
Nadeau, R. M., Foxall, W. and McEvilly, T. V.（1995）. Clustering and periodic recurrence of microearthquakes on the San Andreas Fault at Parkﬁeld, California, Science, 267, 503–507.
Nishenko, S. P. and Buland, R.（1987）. A generic recurrence interval distribution for earthquake forecasting, Bulletin of the Seismological Society of America, 77, 1382–1399.
Nomura, S., Ogata Y., Komaki, F. and Toda, S.（2011）. Bayesian forecasting of recurrent earthquakes
and predictive performance for a small sample size, Journal of Geophysical Research, 116,
B04315, doi:10.1029/2010JB007917.
野村俊一, 内田直希, 尾形良彦（2013a）. 時空間非定常 BPT モデルと北日本太平洋沖沈み込み帯の小繰り

の解析例
（1）
Nonstationary space-time BPT model, and analysis
返し地震（1993.7.15–2011.3.10）of repeating earthquakes on the upper surface of the subducting Paciﬁc Plate from July 15,
1993 till March 11, 2011, 地震予知連絡会会報, 90, 49–53.
野村俊一, 内田直希, 尾形良彦（2013b）. 時空間非定常 BPT モデルによる北日本太平洋沖沈み込み帯—特に主要地震の前と後の変化についての解析例（2）の小繰り返し地震（1993.7.15–2011.3.10）—Nonstationary space-time BPT model, and statistical analysis of repeating earthquakes on
the upper surface of the subducting Paciﬁc Plate from July 15, 1993 till March 10, 2011, on
emphasis of before and after major earthquakes, 地震予知連絡会会報, 90, 54–57.
Nomura, S., Ogata. Y. and Nadeau, R. M.（2014）. Space-time model for repeating earthquakes andanalysis of recurrence intervalson the San Andreas Fault near Parkﬁeld, California, Journal ofGeophysical Research, 119, 7092–7122.
Ogata, Y. （1999）. Estimating the hazard of rupture using uncertain occurrence times of paleoearthquakes, Journal of Geophysical Research, 104, 17995–18014.
Okada, M., Uchida, N. and Aoki, S.（2011）. Statistical forecasts and tests for small interplate repeating
earthquakes along the Japan Trench, Earth, Planets, and Space, 64, 703–715.
Shimazaki, K. and Nakata, T.（1980）. Time-predictable recurrence model for large earthquakes, Geophysical Research Letters, 7, 279–282.
宇津徳治（1957）. 地震のマグニチュードと余震の起こりかた, 地震第 2 輯, 10, 35–45.
Utsu, T.（1961）. A statistical study on the occurrence of aftershocks, Geophysical Magazine, 30, 521–
605.
Utsu, T.（1984）. Estimation of parameters for recurrence models of earthquakes, Bulletin of Earthquake
Research Institute, University of Tokyo, 59, 53–66.

104

Proceedings of the Institute of Statistical Mathematics Vol. 63, No. 1, 83–104 (2015)

Point Process Models for Recurrent Earthquakes at Active Faults
and Their Long-term Forecast
Shunichi Nomura
Graduate School of Information Science and Engineering, Tokyo Institute of Technology

There are many active faults in inland Japan that have risks of future catastrophic
earthquakes. However, their occurrence probabilities vary widely because their earthquakes recur quasi-periodically and their cycles are very diﬀerent. This paper introducessome point processes such as renewal processes for evaluating hazards of such ‘recurrent
earthquakes’. Since inland active faults in Japan have very long cycles of activity, forecasting is diﬃcult because of the scarcity and unreliability of historical data. For stabilityof predictive performance, Bayesian prediction may be appropriate to deal with the uncertainty of data and parameters caused by these problems. Some illustrations are shownon forecasting for active faults in Japan.
In recent years, small earthquakes repeating on subduction zones of crustal plates
have also attracted much attention for earthquake prediction. Their recurrence cycles
are relatively short and so they can be used for middle-term forecasting experiments tovalidate prediction methods. Small repeating earthquakes are also analyzed for monitoring interplate slip because their recurrence times reﬂect the quasi-static slip rate aroundtheir fault patches. Here, a space-time model extended from a renewal process is used toestimate the spatio-temporal distribution of slip rate on plate boundaries.

Key words: Recurrent earthquakes, long-term forecast, Brownian Passage Time distribution, renewal
processes, Bayesian prediction.