統計数理（2021）
第 69 巻 第 2 号 145–163©2021 統計数理研究所

特集「Hawkes 過程の新展開と応用」［総合報告］

統計地震学における ETAS モデル—その進展とホークス型モデル
庄 建倉1 ・尾形 良彦2（受付 2021 年 2 月 15 日；改訂 7 月 26 日；採択 7 月 26 日）

要旨

統計地震学における標準的な地震活動のモデルとして，Epidemic-Type Aftershock Sequence（ETAS）モデルについて，その歴史，理論，方法，応用における新しい進展，および一般的な
ホークス過程に与える影響について紹介する．
キーワード：ETAS モデル，ホークス過程，統計地震学，点過程．
1.

はじめに

統計的地震学は，地震の発生様式を記述する統計的モデルや統計的推論手法を開発し，地震データから物理的メカニズムを定量的に理解し，高い確率利得で地震の確率予測をめざす学問である．他にも，地球物理学的インバージョンなどのベイズ法など，伝統的な地震学の研究に必要な統計的手法の開発も含まれている．統計的地震学は地震物理学，地震予測，地震工学，
防災科学など地震学の学際的な分野に広く活用されている．
統計的地震学の歴史は 150 年に亘り，物理モデルと統計モデルの間のギャップを埋め，地震の予測および被害の軽減を目指した学問であるが，過去数十年の間に急速に発展した．中でも
1970 年代に，確率論
（点過程）
，情報理論
（エントロピー），統計科学（リスク解析）の境界分野で育まれた共通の概念「条件付き強度関数」
は今日の研究の礎となる重要概念である．これは一定の時間区間で地震が生起する確率の微分関数（予測発生率）で，時間・空間的に変動する地震発
生率を記述するために自然な扱いが可能である．そのうえで予測性能の評価が「確率利得」として定量化できる．条件付き強度関数の具体的な一例として，尾形（Ogata, 1985, 1988）は ETAS
（Epidemic-Type Aftershock Sequence）
モデルを提案した．これは，今日では，伝統的な統計地
震学のポアソン過程モデルに代わる標準的な地震活動の帰無仮説モデルとして，よく引用されている．他の幾多の地震活動モデルとの比較対象となっているだけでなく，地震活動の異常現象の検出や仮説研究に採用されている．
ETAS モデルの先駆けであるホークスの点過程（Hawkes, 1971a, 1971b）
は，早くから「条件付き強度関数」で特徴付けられ，当時として画期的な確率的点過程へのアプローチであったが，最近，金融・保険の分野で注目されるまで論文引用回数がそれほど多くなかった．ETAS モデル
は数学的にはホークス過程の地震系列への具体化の一例であるが，既に数十年の歴史がある．
本稿では，ETAS モデルが統計地震学の分野でどのように発展してきたのか，また，他分野で
1

統計数理研究所：〒 190–8562 東京都立川市緑町 10–32

統計数理研究所 名誉教授：〒 190–8562 東京都立川市緑町 10–3

統計数理

146

第 69 巻 第 2 号 2021のホークス過程の統計的推論の理論や手法の発展にどのような影響を与えてきたのか，その歴史を解説する．
2.

初期の歴史

地震統計学は 19 世紀後半に世界に先駆けて日本地震学会が創設されるとともに始まった．
近代的な地震計の設置により，世界中の地震の発生を検出できるようになったため，その発生時刻や震源位置を計算し，比較的正確な地震カタログを作成することができるようになった．
この時期から長い間，地震研究における統計学の応用は，図示による記述統計，線形回帰や相関解析などの一般的な統計手法であり，個別の研究に散在していた（安芸, 1956）．
2.1 グーテンベルク・リヒターの法則と大森・宇津式
Gutenberg と Richter は，地震の規模（マグニチュード）の頻度に関する式
log10 (N ≥ m) = a − bmを発見した（Gutenberg and Richter 1944）．ここで，N (≥ m) は，与えられた地域と期間においてマグニチュード m 以上の地震の数であり，定数 b は b 値と呼ばれ，今日に至るまで地震学の論文で最も頻繁に引用されている．マグニチュードの下限が mc のとき，これは確率分布として次のような指数分布になる．
(2.1)

P {magnitude ≥ m | magnitude ≥ mc } ≈N (≥ m)= 10−b(m−mc ) = e−β(m−mc ) ,N (≥ mc )ここで β = b ln 10 である．本稿で上式は頻出するので，以後 G-R 則と略記する．
は 1891 年の濃尾地震の余震発生頻度の減衰について，先ず，物理
他方，大森（Omori, 1894）
学に良く見られる減衰率として指数関数をデータに当てはめたが，満足のいく結果が得られなかった．結局，彼は 1 日当たりに発生する余震の数がn(t) = K(t + c)−1で良く記述できることを発見した．ここで，t は本震発生からの時間であり ，K と c は定数である．
宇津（1957）は，余震頻度と経過時間の両対数のグラフを使って，余震頻度の減衰が一般的にn(t) = K(t + c)−pとなることを示した．宇津はこの式を「修正大森式」と名づけたが，今日では「大森・宇津式」と呼ばれている．
2.2 統計的地震学における点過程の導入
Vere-Jones and Davies（1966）
および Vere-Jones（1970）の点過程による統計解析は，統計地
震学のマイルストーンである．すなわち，地震発生時刻の系列を記述するために点過程を用いることを提案し，関数解析やスペクトル解析を生成するツールを開発した．その例として，
地震群の記述に Neyman-Scott モデル（Neyman and Scott, 1953, 1958）や Lewis-Barlett モデルの 2 種類の点過程を採用している．Neyman-Scott 過程は銀河宇宙（Bartlett, 1963; Lewis, 1964）の位置分布の記述に，Lewis-Barlett モデル過程は降雨過程のモデル化に用いられていた．
しかし，これらを地震活動のモデルとして実装した「トリガーモデル」
は余震を誘発できる地とできない地震（第 2 種地震）を予め分けなければならない．それらは暗黙的に震（第 1 種地震）

統計地震学における ETAS モデル — その進展とホークス型モデル

147

本震と余震を区別したものであるが，Vere-Jones（1970）の英国王立統計学会での招待討論で，宇津は後述するように，二次余震の存在を指摘している．したがって地震系列データにおける
第 1 種地震と第 2 種地震の種分けには組み合わせ論的複雑性（e.g., Baudin, 1981）があり，これを乗り越えなければならない．なので，明示的な形での条件付き強度を書くことが出来ない
ため，尤度法や AIC で適合度を議論できない．その代わりにモーメント法に基づいた統計量のグラフ化で推論や適合度解析が適用された．特記すべきことに Hawkes の学生 Adamopoulos（1976）は，トリガーモデルと後述のホークスモデルをスペクトル尤度（Whittle, 1962）で推定した．尤度関数の代わりにスペクトル尤度を計算すれば，AIC でモデル比較も可能だが，精度は点過程の最尤法に比べて良くない．
条件付き強度関数はλ(t)dt = P{t と t + dt の時間に事象が発生する | 時刻 t 迄の観測データ}.

によって定義され，過去の地震過程の歴史や外部観測データを条件に，直近未来の地震発生の期待確率を与えたものであり，その重要性に鑑み Vere-Jones（1973, 1975, 1978）は地震学コ
ミュニティに点過程の最尤法の採用を推奨している．条件付き強度関数を使用する利点は，推定やシミュレーションを含む予測の目的のために自然な概念であるということである（Ogata,．
1978, 1981）
地震発生のデータにマグニチュードがあると，マーク付き点過程になる．例えば，余震の規模も含めて予測する単純なモデルは
λ(t, m) =Ks(m)(t + c)pになる．ここで t は本震発生からの経過時間で，s(m) は G-R 則に対応するマグニチュードの確率密度関数である．このモデルは G-R のマグニチュード頻度関数と大森・宇津式を組み合わせたもので，1990 年代から，余震活動を予測するために Reasenberg and Jones（1989, 1994）によって実装された．
2.3 ホークス過程の誕生
ほぼ同時期に，Hawkes は一連の論文で自己および相互誘発モデルを条件付き強度関数で与え，点過程の理論スペクトル関数を計算した（Hawkes, 1971a, 1971b）．要するにホークス過程
は，背景の定常 Poisson 過程のもと，過去に発生した全ての事象の誘発効果で生成された事象からなる点過程である．それぞれの事象は，それが背景事象であろうと誘発された事象であろ
うと，或る確率規則に従って，事象の発生を順次誘発（励起）する．このモデルは，次のような形の条件付き強度を持っている．
λ(t) = μ +g(t − ti ).

i : ti <t

ただし μ は背景
（バックグラウンド）率を表し ，g(t) は自己励起（self-exciting）効果を表す．こ∞の点過程が安定（定常）であるために条件 ρ = 0 g(u)du < 1 を必要とする．ここで，ρ は臨界パラメータと呼ばれる．つまり，任意の事象から直接的に誘発された事象の平均数が 1 以下であることを示している．具体的にホークス過程が伝染病の感染過程（分枝確率過程）（Kendall,1949）に対応していることが示されている（Hawkes and Oakes, 1974）．
では 2 本の点過程の相互誘発項を持つモデルを提案してい同時に， Hawkes（1971a, 1971b）
る．その一成分の過程の条件付き強度は

148

統計数理

λ(t) = μ(t) +第 69 巻 第 2 号 2021

g(t − ti ) +i : ti <th(t − sj )j : sj <tと書ける．ここで，h(t) は，外部過程 sj : i = 1, 2, . . . , Ns からの，励起の応答関数を表す．
3.

ETAS モデル以前の統計的地震学

ETAS モデルが提案される前の研究として（1）
大森・宇津式の提案と余震解析への利用，（2）
相互誘発ホークスモデルによる地震活動解析への適用が挙げられる．本節では，これらを解説する．
3.1 大森・宇津式による余震データの解析
1950 年代以降，大森・宇津式は余震活動の解析に広く用いられてきた．宇津は余震頻度対経過時間の両対数表示で国内外の 200 以上の余震系列について作成し，p 値が 0.6 から 2.5 の範囲で，中央値は 1.1 であることを示した．また p 値の推定値と本震のマグニチュードや下限マ
グニチュードとの間には相関はないことを発見し，余震列ごとに b 値が下限マグニチュードに関わらず不変な定数であることを示している（Utsu, 1962, 1971; Utsu et al., 1995）．さらに宇津は，余震の統計に関連したいくつかの経験則を提案しており，その中には余震の面積を本震
マグニチュードでスケーリングする宇津・関の法則（宇津・関, 1955）も含まれている．
宇津は，本震だけでなく，大きな余震が更なる余震（2 次余震）を誘発する場合があることを観測した．このような現象を次のような「多重の大森・宇津モデル」
で示した．

 Ki H(t − ti )K+.
p(t − t0 )(t − ti + ci )piNT

λ(t) =

i=1

ここで t0 は本震の発生時刻，ti , i = 1, 2, . . . , NT は有意に大きな余震の発生時刻を示し，H はHeaviside 関数である．この観測は，地震の誘発に関して，本震と余震は違う性質を持つとした
当時までの地震学的常識を覆す事実を定式化したもので，画期的といえる．
一般の地震列に多重の大森・宇津式を適用する上での一つの難点は，どの地震が他の事象を誘発するかを判断することである．最大の余震は二次的な余震を伴うことが多いが，常は赤池情報量規準（AIC, Akaike, 1974 参照）
を用にそうとは限らない．Ogata（1983a, 1983b）
いて，大森・宇津減衰の二次的余震の関数の重ね合わせが何項必要なのかをデータ解析で示は，1965 年のアリューシャン列島地震の余震列している．また Ogata and Shimazaki（1984）{ti : i = 1, 2, . . . , N } が λ(t) = K/(t + c)p の単純な大森・宇津式に従わないことを，時間変換
ti → τ i = ti

λ(u)du.
0を使って示している．すなわち，変換時刻系列 {τi } は最大余震前後で標準的なポアソン過程
から大きく逸脱しているが，顕著な 2 次余震を含む多重大森・宇津モデルの累積関数は 2200 日まで直線的に推移した後上方に逸脱し，この余震活動が背景の常時地震活動に埋もれたことをは，同じ期間のデータに ETAS モデルを当てはめれば，そ示している．後ほど，Ogata（1989）の後も逸脱せず直線的に推移していることを示した．
が大森・宇津の式を基準にして変換された時刻系列の静穏化
なお，Matsu’ura（1986, 1991）
が，大きな余震が発生する前に見られる場合があることを明示した．すなわち，Ogata andで使用した変換時間の発生率が標準的なポアソン過程の発生率から下方に逸Shimazaki（1984）
脱する場合を多数の余震例で示した．このような現象は「相対的」静穏化と呼ばれている．この

統計地震学における ETAS モデル — その進展とホークス型モデル

149

概念は，従来から直感的・経験的に議論されている地震活動の静穏化を明確に定義したものの検で，これは ETAS モデルに基づく地震活動の相対的静穏化（Ogata, 1985, 1988, 1992 など）
出に広く用いられている．
3.2 点過程による地震活動のトレンド，季節性，および相関解析
尾形は一連の論文（尾形, 1981; Ogata and Akaike, 1982; Ogata et al., 1982）で，Hawkes の相互励起過程を

λ(t) = a0 + PJ (t) + CK (t) +g(t − ti ) +i : ti <th(t − sj )j : sj <tの形に拡張して地震活動変化の解析に応用した．ここで ，a0 は定数項 ，PJ (t) は背景地震活動
K+ b2k sin 2kπ} は季の長期的トレンドを表す J 次の多項式，そして Ck (t) = k=1 {b2k−1 cos 2kπ
T0
T0
節変化の項である．g(t) と h(t) は地域内の地震による自己励起項と地域外地震からのそれぞれの励起項を表し，Laguerre 型の多項式g(t) = e−αt

K1

k=0

pk tk ,

h(t) = e−βtK2
q k tk ,k=0を採用している．ここで，K1 と K2 は非負の整数であり，p1 , p2 , . . . , pK1 , α, q1 , q2 , . . . , qK2 , β は推定するパラメータである．対数尤度計算の詳細については，Vere-Jones and Ozaki（1982），
Ogata and Akaike（1982），計算マニュアル（Utsu and Ogata, 1997; Ogata, 2006）を参照された
い．トレンド項の個数 J ，季節変化項の個数 K ，および励起項の個数 K1 と K2 は赤池情報量規
準（AIC）で決定できる．初期の研究では，これらのモデルは主に，様々な地域での地震活動同士の誘発効果を調査するために使用された．たとえば，尾形（1981），Ogata et al.（1982），Ogataは，沈み込むプレートに沿った深発地震と浅発地震の発生の一方向への移and Katsura（1986）の統計的因果関係を探査し，様々な地域での地震活動度
動性（Mogi, 1968, 1973; 宇津, 1975）の季節変化が年降雨量の変化パターンと対応することを示した．同様に，Ma and Vere-Jones（1997）は，ニュージーランドの浅い地震活動の季節性を発見した．Matsumura（1986）は，長期間の不均質データに対して，上記の長期トレンドを検出率の増加としてモデル化し，世界中の地震活動の季節性の有無を分析した．年降水量変化が顕著な中緯度内陸部では地震活動の季節性が有り，反対に低緯度地域や海域地帯では季節性が見られないことを示した．これらは，後述するような，流体貫入による誘発地震の可能性の傍証となっている．
また，このような統計的因果関係のモデルは，各種の地球物理学的観測変数の先駆的異常現象の検出や統計的有意性を議論するためにも使用された（Nishizawa and Nagao, 1994; Zhuang．
et al., 2005b, 2013a）
4.

ETAS モデル

4.1 ETAS モデルの誕生
多重の大森・宇津式で説明されているように，余震活動は多くの場合，二次余震を明確に含んでいる．しかし，一次余震との区別は一般的に明瞭でなく，地震系列データから二次余震を分離するのは難しい（特殊なパラメータ化によるトリガーモデルで，ある程度分離可能（Ogata,であるが予測には向かない)．
2001）
Ogata（1988）は，地震発生をモデル化するにあたって，誘発する地震と誘発される地震を区別しないで，各地震は多かれ少なかれ誘発すると考え，大森・宇津式をマグニチュードの効

150

統計数理

第 69 巻 第 2 号 2021

果の重み付けとして誘発効果を重ねた．このアイデアは Hawkes の自己励起過程と同様に，伝に由来している．それ故これは感染型の余震系列
染病の感染過程（分枝過程）（Kendall, 1949）（Epidemic Type Aftershock Sequence; ETAS）
モデルと名付けられた．
ETAS モデルを条件付き強度関数で表現すると，考慮する地震系列のマグニチュードの下限
を m0 としたとき，
λ(t) = μ +κ(mi )g(t − ti )i : ti <tである．ここで g(u) は大森・宇津式の正規化形式 p−1(1 + ct )−p ，つまり誘発する地震を起点
c
として，その余震が誘発される非定常ポアソン過程の強度関数を正規化した確率密度関数である．これらに対する重み付け関数 κ(m) = A exp[α(m − m0 )] は地震マグニチュード m の大小に拠って直接的に誘発される地震の期待数である．余震の大きさ（マグニチュード）は誘発する地震よりも小さくなっている必要はない．その後たまたま大きな地震が誘発された場合，前者は「前震」と呼ばれる．上記の条件付き強度は，Ogata（1988）では履歴のマグニチュードを既に与えられたものとして説明変数に組み込んでいるが，それらのマグニチュード列は外的データ入力変数と考えられ，それら自身の分布はモデル化されていない．しかし多くの研究者は，地震系列のシミュレーション実験などで理論的議論をするために，独立な指数分布つまり G-R 則の発生頻度を仮定している．すなわちマーク付き点過程の条件付き強度関数として下記の形のマーク付きモデルλ(t, m) = s(m) μ +κ(mi )g(t − ti )i : ti <tを考えている．ここで s(m) = βe−β(m−m0 ) , (m ≥ m0 ) は，G-R 則の確率密度関数形式である．
このマーク付きモデルのもとで，ETAS モデルは地震活動を記述するための標準モデルとして．
受け入れられている（Huang et al., 2016 の特集号参照）
このマーク付き標準 ETAS モデルで計算すると，臨界パラメータは，任意のイベントによって直接的に誘発された地震群の平均数，

 ∞

=

s(m)κ(m)dm =
m0

Aβ
β−α

α
が必要十分条β件である． ≥ 1 の場合 ，t が増加すると，単位時間間隔内のイベントの発生率が発散して無限大になる．さらなる ETAS モデルの漸近特性の詳細については，Helmstetter and Sornette（2002），Zhuang and Ogata（2006），Zhuang et al.（2013b）および Saichev and Sornette（2007）を参照されたい．
但し，地震マグニチュード時系列の独立同分布の G-R 則は実証されているわけでなく，近似的なものと理解すべきである．非独立や非同一の確率法則のモデル化は大地震の確率予測の観点から議論されている（Utsu, 1970a, 1970b; Ogata and Katsura, 2014; Ogata et al., 2018; 野．
村・尾形, 2020; Zhuang et al., 2004）となる． < 1 の場合 ETAS モデルは安定する．このため β ≥ α と A < 1 −4.2 時間モデルから時空間モデルへ
ETAS モデルが時空間 ETAS モデルに一般化される前に，Musmeci and Vere-Jones（1992）
は，時空間拡散点過程モデルを使用してイタリアの地震活動を分析した．それらのモデルの条
件付き強度関数は

統計地震学における ETAS モデル — その進展とホークス型モデル

λ(t, x, y) = μ(x, y) + A eαmi e−c(t−ti )ti <tとλ(t, x, y) = μ(x, y) + A

2πσx σy (t − ti )


ti <t

exp

−

1
2(t − ti )151

(x − xi )2(y − yi )2
+
2
σx
σy2

Aeαmi e−c(t−ti ) (t − ti )2 Cx Cy
2π [(x − xi )2 + (t − ti )2 Cx2 ][(y − yi )2 + (t − ti )2 Cy2 ]の 2 種類である．ここで，A, α, σx , σy , Cx ，および Cy は定数である．これらによると，固定位
置 (x, y) で t → ∞ のとき，それぞれのモデルの余震発生率は漸近的に t−1 e−ct または t−2 e−ct
に従って時間とともに減衰する．Kagan（1991）と Rathbun（1993）も独自の時空間モデルを与えた．
これらに対して，現在多用されている時空間 ETAS モデル（Ogata, 1998）の条件付き強度λ(t, x, y) = μ(x, y) +κ(mi )g(t − ti )f (x − xi , y − yi ; mi ),i : ti <tは地震活動に調和的である．ここで，
κ(m) = Aeα(m−m0 )はマグニチュード m の地震から誘発される地震の期待数である．そして時間関数g(t) =(p − 1)cp−1 (t + c)−p ,0,

t > 0;
t≥0は誘発に至る時間分布の確率密度関数（大森・宇津式）
で，空間的広がりに関して
f (x, y; m) =

1
f0
πσ(m)

x2 + y 2σ(m)はマグニチュード m の地震から誘発される地震の相対的位置に関する密度関数である．
f0 (x, y) と σ(m) は様々な関数が考慮されているが，通常 f0 は二次元正規分布の確率密度関数
−q
− ω2
1
f0 (ω) = 2D
またはスケール D を含む逆ベキの関数 f0 (ω) = q−1
1 + Dω2が考えら
2 e 2DD2γα
れている．そしてスケーリング関数は σ(m) ∝ κ(m) または σ(m) ∝ [κ(m)] である．以上でD, q ，および γ は定数パラメータである．Zhuang et al.（2004），Zhuang（2006）
および Ogataγに示されているように，逆ベキの関数 f0 と σ(m) ∝ [κ(m)] α は経験的に地and Zhuang（2006）

震データによく適合する．空間的応答関数 x2 + y 2 を √ 1 2 σσ21 x2 − 2ρxy + σσ12 y 2 で置き換
1−ρ

えれば余震域の非等方形状を近似する（Ogata, 1998; Ogata et al., 2003b; Ogata and Zhuang,
2006）．
時空間 ETAS モデルは地震活動解析に広く使用されている（Ogata, 1998; Ogata et al., 2003a;
Zhuang et al., 2002, 2004; Console et al., 2003; Helmstetter et al., 2003; Lombardi et al., 2010;．これらは，主要地震国の研究機関や政府機関によって，主要な地震予Guo et al., 2015b など）．とくにアメリカ合衆国地質調査所
測モデルとして採用されている（Schorlemmer et al., 2018）（USGS）は，ETAS モデルを第 3 次全カリフォルニア地震確率予測モデル（UCERF3）の短期予測に UCERF3-ETAS モデルとして採用した（Field et al., 2017）．
4.3 ETAS モデルに関連する方法論の開発
4.3.1 確率的除群法
Zhuang et al.（2002）
は以下のような「確率的除群法」を開発した．時空間 ETAS モデルによ

統計数理

152

第 69 巻 第 2 号 2021

ると j 番目の地震の発生時空間座標 (tj , xi , yj ) での背景（常時）地震発生率への寄与の割合は，μ(x ,y )ϕj = λ(tj ,xj j ,yj j ) である．これは，j 番目の地震が背景地震である確率である．同様に，j 番目の地震が i 番目の地震によって誘発される確率は，
ρij =κ(mi )g(tj −ti )f (xj −xi ,yj −yi ;mi ),λ(tj ,xj ,yj )0,i<jOtherwiseとなる．確率 ϕj と ρij の各地震 j を選択すると，背景地震活動過程や i 番目の地震によって誘発された地震活動の過程をそれぞれ分離し，全体を確率分枝過程と考えることができる．
地震活動静穏化などを定常ポアソン仮説で検定するために，地震の群れを取り除く必要があった．そのような常時地震データを得るための除群法は，従来から多数提案されたが，地震群の複雑な混合の解釈の違いで，合理的な除群法の合意に関して困難に面している．たとえば，二次余震を含む余震系列の除群法は極端に難しい．これに対して上記の ETAS モデルを用いた分離アルゴリズム（確率的除群法）
は，時空間の地震カタログから背景の地震活動を決定論的ではなく，確率論的に分類するため，多数の除群地震カタログをシミュレーションで作成できる好ましい方法になった（van Stiphout et al., 2012）．
確率 ϕj と ρij を直接扱うことで乱数を変えて確率的除群を繰り返すことにより，様々なバー
ジョンの除群地震パターンを取得できる．このような多数パターンの実現は，地震の群を除群する際の不確定性を示すものであり，地震クラスターの統計的特性や静穏化の有意性を評価す
るのに役立つ．
例えば，ある地域の累積常時地震活動度をS(t) =

ϕii : ti <tと定義している．この関数を用いることで，ある地域の常時地震活動が何時何処で静穏化して．
いるかどうかを調査できる（例えば，Wang et al., 2010; Guo et al., 2017; Zhuang et al., 2017）
や Nishikawa et al.（2019）は，確率的除群法を用いて，最近では Nishikawa and Ide（2017, 2018）
沈み込み帯の群発地震や「ゆっくりすべり」の研究を展開している．また Ueda and Kato（2019）
は確率的除群法を用いて地震活動の季節性の時空間的詳細を議論している．
4.3.2 時空間 ETAS モデルの推定
Zhuang et al.（2004）
は，各クラスタリング成分の「残渣分析」（residual analysis）を実行する
ため，確率的除群法で分離した分枝過程をシミュレーションして再構成するアルゴリズム「確率的再構成法」を提案した．このアルゴリズムは，一般的なホークスモデル（Zhuang, 2006;，または複雑な背景強度を持つホークスモデル
Marsan and Lengliné, 2008; Mohler et al., 2011）のノンパラメトリックな推定手法として発展した．
（Zhuang and Mateu, 2019）
Zhuang et al.（2002）
は，時空間 ETAS モデルの背景率関数 μ(x, y) を求めるため，確率的除群法で得られた地震をノンパラメトリックに推定し，それを背景率関数にした ETAS モデルを再推定し，これで元データを再除群をする．空間的に一様な発生率から始めて，この操作を反復してモデル推定を収束させ，背景率関数を求める方法を提案した．同様に，確率的除群法にと Li et al.（2020）はパラメータ推定のための期待値最大
基づき，Veen and Schoenberg（2008）
化アルゴリズムを開発した．
4.4 ETAS モデルの応用
ETAS モデルを適用して，以下のように多くの興味深い課題が研究されている．

統計地震学における ETAS モデル — その進展とホークス型モデル

153

4.4.1 余震列の Båth の法則
先ず余震統計でよく知られている経験則である Båth の法則とは，本震と最大余震のマグニ
チュード差 M0 − M1 の中央値は 1.2（Båth, 1965）であると主張するものである．より一般的ながもたらしたもので，次のようになる．
結果は，Utsu（1961, 1970a）
M0 − M1 = c1 M0 + c2

ここで ，c1 と c2 は或る定数である．この法則は余震列によって引き起こされる可能性のある最大被害を予測するのに役立つ．Helmstetter and Sornette（2003a)，Saichev and Sornette（2005a)，Zhuang and Ogata（2006），Vere-Jones and Zhuang（2008）および Luo and Zhuang（2016）は，余震に ETAS モデルと G-R マグニチュード分布を当てはめ，Båth 法則が任意のクラスター内の最大イベントの極値分布（二重指数分布）となることを示した．
4.4.2 前震仮説の検証
G-R 独立マグニチュード列を仮定した時空間 ETAS モデルでシミュレーションしたとき，地震群の中で最大になった地震の直前の地震群が，伝統的な前震の時・空間・マグニチュード列の統計的性質の殆どを再現したという Helmstetter and Sornette（2003b）の結果は，意外性を伴って，多くの地震学者の興味を引いた．余震の性質で構成された ETAS モデルで前震現象が説明できるかどうか，つまり，伝統的な前震概念と無関係なのかどうかは興味深い問題である．
元来 ETAS モデルにおける，そのような前震確率の理論は，Saichev and Sornette（2005a）とZhuang and Ogata（2006）
によって示されていた．これに続いて，Zhuang et al.（2008）は日本，
ニュージーランド，南カリフォルニア地域の地震の前震確率を理論計算やシミュレーションで得られた値と比較し，Helmstetter and Sornette（2003b）と同様の結論を出している．その後，多くのシミュレーション研究（例：Marzocchi and Zhuang, 2011; Lippiello et al., 2012; Bouchonが行われたが，多様な結論が示された．
et al., 2013; Ogata and Katsura, 2014）
4.4.3 人工的な誘発地震
廃液処理や貯水池の貯留などの人為的な活動に起因する誘発地震が注目されている．流体貫入による誘発地震とストレス蓄積による自然地震の違いを定量化し，地震ハザードを評価する
ために ETAS モデルが適用されている（Lei et al., 2008, 2020; Llenos and Michael, 2013, 2016;．Llenos and Michael（2013）
は，Ogata（1992）が提案した変化点検出技術を用
Jia et al., 2020）
いて，米国オクラホマ州とアーカンソー州の両方で常時地震発生率と誘発パラメータの変化を発見し，流体注入後に発生した群と自然現象によって引き起こされた群とは異なると結論づけた．Bachmann et al.（2011）は，スイスのバーゼルの EGS 地点の地震活動データに ETAS モ
デルを適用することで，常時地震発生率が各時間窓における揚水履歴と良好な相関関係を示しは，同様の手法を用いて，福島県柳津西山熱水鉱区ていることを発見した．Eto el al.（2013）における誘導地震活動についても同様の結論を導き出した．ETAS モデルやその拡張形を使っや火山性の誘発群発地震（Kumazawa et al., 2016）の研究でた，自然地震（Ogata et al., 2003a）
は，既知・未知の原因による地震活動パターンの変化を検出するのに有効であることが示されている．
4.4.4 ETAS モデルと速度・状態依存の摩擦則
ETAS モデルは，時間変換された時系列で地震活動率の変化点を検出すること（Ogata,，または残渣解析を用いて背景地震発生率の変化点を検出することができる
1988, 1989, 1992）．このような変化を説明するために，「すべり速度と状態に依存す（Zhuang et al., 2005a, 2019）
る摩擦則」が用いられてきた（例えば，Ogata, 2010; Jia et al., 2014, 2018）．

154

統計数理

第 69 巻 第 2 号 2021

4.5 時空間 ETAS モデルのさらなる拡張と最近の進歩
4.5.1 位置依存の ETAS パラメータ
Ogata et al.（2003a）と Ogata（2004）
は様々な場所での地震活動クラスタリング構造の違いを説明するために，階層的時空間 ETAS（HIST-ETAS）モデルを開発した．モデル内の各パラメー
タが場所の関数（震央を頂点とするデローネ三角形網上の部分的線形関数）となる．これらに滑らかさの事前分布を仮定して，パラメータの空間的変動の形でクラスター特性の地域性を推定する．そのために開発した，ペナルティ付き尤度を備えたベイズ的推定ソフトウェアを公開し．特に HIST-ETAS モデルは常時地震活動度の解析から長期予測に応用た（Ogata et al., 2021）．非ベイズのアプローチ
されている（e.g. Ogata et al., 2018; 尾形, 2020; Ogata and Omi, 2020）
は時空間 ETAS モデルの空間変動を推定するための重み付き尤度法を開として Zhuang（2015）
発した．
4.5.2 自己相似 ETAS モデル
Vere-Jones（2005）は，ETAS モデルの下限マグニチュードの選定に起因する問題を回避するために，自己相似性を満たす ETAS モデルを開発し，完全な自己相似的な特徴を実現した．
しかし統計学的には，本震直後の余震や微小地震のデータの欠測問題のため困難に直面していにとどまっており，実際る．このモデルは理論的な議論（Saichev and Sornette, 2005b も参照）の地震活動データにはまだ適用されていない．自己相似性を導入するための別の展開として，モデル（Turcotte et al., 2007; Holliday et al., 2008）があるが， こ
分岐型余震シーケンス（BASS）の過程で発生する全ての地震の G-R 則が破壊されてしまうため，成功しているとは言えない（Zhuang, 2013）．
4.5.3 地震群の非等方性と震源深度の影響
大きな地震の余震クラスターは通常，本震の震央の周りでは等方性ではなく，本震の断層に沿って分布している．Hainzl et al.（2008）は，時空間 ETAS モデルのシミュレーションで，この性質を無視すると，偏ったパラメータ推定値が得られ，特に α パラメータは過小推定量になる可能性があることを示した．各余震クラスターの非等方性の効果に対処するために，元来Ogata（1998, 2004）と Ogata et al.（2003b）
は楕円形等高線の空間分布関数を使っている．さらに Ogata et al.（2018）は深さを考慮した 3 次元 HIST-ETAS モデルでは楕円体等高面の空間分
布関数を採用し，関東地方直下のデータにあてはめている．
は，大地震の破壊域形状による余震配置への影響を考慮
最近 Guo et al.（2015b, 2017, 2019）
し，断層に沿った二次元有限震央域に基づく ETAS モデルを開発した．このモデルでは，大地震の発生源は，もはや空間点震源ではなく，空間的に広がる震源域である．このモデルの定式化では，地表に投影面積 S の有限の震央域を持つ地震の場合，空間の応答関数が

f (x, y; S) =S

w(x − u, y − v)τ (u, v)du dv

Sτ (u, v)du dvになる．ここで ，(u, v) は震央域 S 上の位置であり，w(x, y; mi ) =

q−1
πD 2

1+

x2 + y 2
D 2

−q

,

ここで D は新しく導入されたパラメータである．さらに τ (u, v) は震央域に沿った本震の非一様な余震生産性の密度関数であり，分母の積分は分子を正規化するためのものである．さらは，震源の深さのデータも組み込んだ 3 次元有限震源域の時空間に Guo et al.（2015a, 2018）
ETAS モデルも提案した．そのモデルでは空間応答関数が水平面と深さの直積で

統計地震学における ETAS モデル — その進展とホークス型モデル

155

f (x, y, z; mi , zi ) = f (x, y; mi )h(z, zi )となる．ここで f (x, y; mi ) は前述の二次元空間応答関数と同じで，深さ方向の応答関数を以下のように与えている．すなわち，h(z, zi ) = HB η

1Zi
H − Zi
+ 1,
+1
H
H



z
H

η ZHi 

1−

z
H

i
η H−ZH,

ここに B(p, q) = 0 tp−1 (1 − t)q−1 dt はベータ関数，そして定数 H は考慮する深さの最大値である．Zhuang et al.（2019）は，3 次元有限震源域時空間 ETAS モデルをイタリアのカタログに適合させた場合，二次元点震央の時空間 ETAS モデルよりも優れていることを示した．
5.

他分野でのホークス過程の応用

地震学における時空間 ETAS モデルの成功は，他の分野でホークス型モデルの適用を促進した．F. R. Schoenberg と彼の UCLA 研究グループは，2000 年代以降，山火事（e.g., Peng et al.,2005），犯罪データ（Molher et al., 2011），ソーシャルネットワーク（Fox et al., 2016; Zipkinなど，多くの分野にホークス型モデルを適用し，その結果はアメリカの主
et al., 2015 など）
要メディアでも注目された．その他に近年のデータ分析として，テロリストの行動パターン，経済における価格変動（e.g., Bacry et al., 2014），およびゲノムまた（e.g., Tench et al., 2016）はニューロン活動（e.g., Truccolo et al., 2005）などがある．
いずれの分野においても，その理論や方法論の大部分は，元来統計数理研究所における点過程解析法の開発や ETAS モデルの応用の成果として，地震データの研究で開発されたものである．いまのところ，他分野でのホークス過程の応用は，主にパラメータ推定と結果説明のためのものに留まっていると思われる．
6.

まとめと将来の展望

ホークス型誘発モデルは，個々の事象や粒子間のクラスタリング効果（正の相互作用）を検出し，それらの間の潜在的な因果関係を決定するのに役立つことから，自然科学と社会科学の両方の点過程データ分析で最も人気のあるモデルの 1 つとなっている．現在では，観測技術やデータ蓄積の急速な発展により，点過程データ解析においてもビッグデータ問題は対象になっている．多重の事象系列や，膨大な数の発生点を含む長時間の事象系列では，事象間の群れや誘発効果を定量化して予測するための迅速な解析・予測ツール開発や一般的な新しい枠組みが望まれている．その為の原型であるホークス型モデルはこの目的に適合するはずである．
ホークス型点過程の実装と応用の中で，ETAS モデルは統計的推論の発展において重要な役割を果たしてきた．これらの豊富な統計的技術により，ETAS モデルは統計地震学において重要な役割を果たすだけでなく，ホークス型点過程を他の多くの分野に広げ，深化することがで
きている．
謝

辞

筆者らは David Vere-Jones 氏から長期間受けた励ましと手厚い支援に感謝している．本研究は，学術振興会の科学研究費補助金（JSPS Kakenhi 17H00727 と 19H04073）により支援された．
また，レビューアからの有益なコメントにも感謝する．

156

統計数理

参

第 69 巻 第 2 号 2021

考

文

献

Adamopoulos, L. (1976). Cluster models for earthquakes: Regional comparisons, Mathematical Geology, 8, 463–475.
Akaike, H. (1974). A new look at the statistical model identication, IEEE Transactions on AutomaticControl, 19(6), 716–723.
安芸敬一 (1956). 統計地震学の現状, 地震 第 2 輯, 8(4), 205–228.
Bachmann, C., Wiemer, S., Woessner, J. and Hainzl, S. (2011). Statistical analysis of the induced Basel2006 earthquake sequence: Introducing a probability-based monitoring approach for enhanced
geothermal systems, Geophysical Journal International, 186, 793–807, doi:10.1111/j.1365246X.2011.05068.x.
Bacry, E. and Muzy, J.-F. (2014). Hawkes model for price and trades high-frequency dynamics, Quantitative Finance, 14(7), 1147–1166.
Bartlett, M. S. (1963). The spectral analysis of point processes, Journal of the Royal Statistical Society.
Series B (Methodological), 25(2), 264–296.
Båth, M. (1965). Lateral inhomogeneities in the upper mantle, Tectonophysics, 2, 483–514.
Baudin, M. (1981). Likelihood and nearest neighbor distance properties of multidimensional Poissoncluster processes, Journal of Applied Probability, 18, 879–888.
Bouchon, M., Durand, V., Marsan, D., Karabulut, H. and Schmittbuhl, J. (2013). The long precursory
phase of most large interplate earthquakes, Nature Geosciences, 6, 299–302.
Console, R., Murru, M. and Lombardi, A. M. (2003). Rening earthquake clustering models, Journal ofGeophysical Research, 108(B10), 2468.
Eto, T., Asanuma, H., Adachi, M., Saeki, K., Aoyama, K., Ozeki, H., Mukuhira, Y. and Markus Häring,M. (2013). Application of the ETAS seismostatistical model to microseismicity from geothermalﬁelds, GRC Transactions, 37, 149–153.
Field, E., Porter, K. and Milner, K. (2017). A prototype operational earthquake loss model for California based on UCERF3-ETAS, a rst look at valuation, Earthquake Spectra, 33(4), 1279–1299.
Fox, E. W., Schoenberg, F. P. and Gordon, J. S. (2016). Spatially inhomogeneous background rate estimators and uncertainty quantication for nonparametric Hawkes pointprocess models of earthquake occurrences, The Annals of Applied Statistics, 10(3), 1725–1756.
Guo, Y., Zhuang, J. and Zhou, S. (2015a). A hypocentral version of the space–time ETAS model,
Geophysical Journal International, 203(1), 366–372.
Guo, Y., Zhuang, J. and Zhou, S. (2015b). An improved space-time ETAS model for inverting the
rupture geometry from seismicity triggering, Journal of Geophysical Research: Solid Earth,120(5), 3309–3323. 2015JB011979.
Guo, Y., Zhuang, J., Hirata, N. and Zhou, S. (2017). Heterogeneity of direct aftershock productivityof the main shock rupture, Journal of Geophysical Research: Solid Earth, 122(7), 5288–5305.
Guo, Y., Zhuang, J. and Hirata, N. (2018). Modelling and forecasting 3d-hypocentre seismicity in the
Kanto region, Geophysical Journal International, 214(1), 520–530.
Guo, Y., Zhuang, J. and Ogata, Y. (2019). Modeling and forecasting aftershocks can be improvedby incorporating rupture geometry in the etas model, Geophysical Research Letters, 46(22),12881–12889.
Gutenberg, B. and Richter, C. F. (1944). Frequency of earthquakes in California, Bulletin of the
Seismological Society of America, 34, 184–188.
Hainzl, S., Christophersen, A. and Enescu, B. (2008). Impact of earthquake rupture extensions on parameter estimations of point-process models, Bulletin of the Seismological Society of America,
98(4), 2066–2072.
Hawkes, A. G. (1971a). Point spectra of some mutually exciting point processes, Journal of the Royal

統計地震学における ETAS モデル — その進展とホークス型モデル

157

Statistical Society: Series B (Statistical Methodology), 33(3), 438–443.
Hawkes, A. G. (1971b). Spectra of some self-exciting and mutually exciting point processes, Biometrika,
58(1), 83–90.
Hawkes, A. G. and Oakes, D. (1974). A cluster process representation of a self-exciting process, Journalof Applied Probability, 11(3), 493–503.
Helmstetter, A. and Sornette, D. (2002). Subcritical and supercritical regimes in epidemic models ofearthquake aftershocks, Journal of Geophysical Research, 107(B10), 2237.
Helmstetter, A. and Sornette, D. (2003a). Bath’s law derived from the Gutenberg-Richter law and fromaftershock properties, Geophysical Research Letters, 30(20), 2069.
Helmstetter, A. and Sornette, D. (2003b). Foreshocks explained by cascades of triggered seismicity,
Journal of Geophysical Research, 108(B10), 2457.
Helmstetter, A., Ouillon, G. and Sornette, D. (2003). Are aftershocks of large Californian earthquakes
diﬀusing?, Journal of Geophysical Research, 108(B10), 2483.
Holliday, J. R., Turcotte, D. L. and Rundle, J. B. (2008). Self-similar branching of aftershock sequences,
Physica A: Statistical Mechanics and its Applications, 387(4), 933–943.
Huang, Q., Gerstenberger, M. and Zhuang, J. (2016). Current challenges in statistical seismology, Pureand Applied Geophysics, 173(1), 1–3.
Jia, K., Zhou, S., Zhuang, J. and Jiang, C. (2014). Possibility of the independence between the 2013
Lushan earthquake and the 2008 Wenchuan earthquake on Longmen Shan Fault, Sichuan,
China, Seismological Research Letters, 85(1), 60–67.
Jia, K., Zhou, S., Zhuang, J., Jiang, C., Guo, Y., Gao, Z. and Gao, S. (2018). Did the 2008 Mw
7.9 Wenchuan earthquake trigger the occurrence of the 2017 Mw 6.5 Jiuzhaigou earthquake in
Sichuan, China?, Journal of Geophysical Research: Solid Earth, 123(4), 2965–2983.
Jia, K., Zhou, S., Zhuang, J., Jiang, C., Guo, Y., Gao, Z., Gao, S., Ogata, Y. and Song, X. (2020).
Nonstationary background seismicity rate and evolution of stress changes in the Changning salt
mining and shale gas hydraulic fracturing region, Sichuan basin, China, Seismological Research
Letters, 91(4), 2170–2181.
Kagan, Y. (1991). Likelihood analysis of earthquake catalogues, Journal of Geophysical Research, 106,
135–148.
Kendall, D. G. (1949). Stochastic processes and population growth, Journal of the Royal StatisticalSociety, Series B, 11, 230–264.
Kumazawa, T., Ogata, Y., Kimura, K., Maeda, K. and Kobayashi, A. (2016). Background rates ofswarm earthquakes that are synchronized with volumetric strain changes, Earth and Planetary
Science Letters, 442, 51–60.
Lei, X., Yu, G., Ma, S., Wen, X. and Wang, Q. (2008). Earthquakes induced by water injection at 3km depth within the Rongchang gas eld, Chongqing, China, Journal of Geophysical Research:
Solid Earth, 113, B10310.
Lei, X., Su, J. and Wang, Z. (2020). Growing seismicity in the Sichuan basin and its association withindustrial activities, Science China Earth Sciences, 63, 1633–1660.
Lewis, P. A. W. (1964). A branching poisson process model for the analysis of computer failure patterns,
Journal of the Royal Statistical Society. Series B (Methodological), 26(3), 398–456.
Li, C., Song, Z. and Wang, W. (2020). Space–time inhomogeneous background intensity estimatorsfor semi-parametric space–time self-exciting point process models, Annals of the Institute of
Statistical Mathematics, 72, 945–967.
Lippiello, E., Godano, C. and de Arcangelis, L. (2012). The earthquake magnitude is inuenced byprevious seismicity, Geophysical Research Letters, 39(5), 051083.
Llenos, A. L. and Michael, A. J. (2013). Modeling earthquake rate changes in Oklahoma and Arkansas:

158

統計数理

第 69 巻 第 2 号 2021

Possible signatures of induced seismicity, Bulletin of the Seismological Society of America,
103(5), 2850–2861.
Llenos, A. L. and Michael, A. J. (2016). Characterizing potentially induced earthquake rate changes in
the Brawley seismic zone, Southern California, Bulletin of the Seismological Society of America,
106(5), 2045–2062.
Lombardi, A. M., Cocco, M. and Marzocchi, W. (2010). On the increase of background seismicity rateduring the 1997–1998 Umbria-Marche, central Italy, sequence: Apparent variation or uid-driven
triggering?, Bulletin of the Seismological Society of America, 100(3), 1138–1152.
Luo, J. and Zhuang, J. (2016). Three regimes of the distribution of the largest event in the critical etas
model, Bulletin of the Seismological Society of America, 106, 1364–1369.
Ma, L. and Vere-Jones, D. (1997). Application of M8 and Lin-lin algorithms to New Zealand earthquake
data, New Zealand Journal of Geology and Geophysics, 40, 77–89.
Marsan, D. and Lengliné, O. (2008). Extending earthquakes’ reach through cascading, Science,
319(5866), 1076–1079.
Marzocchi, W. and Zhuang, J. (2011). Statistics between mainshocks and foreshocks in Italy andSouthern California, Geophysical Research Letters, 38, L09310.
Matsumura, K. (1986). On regional characteristics of seasonal variation of shallow earthquake activities
in the world, Bulletin of Disaster Prevention Institute, Kyoto University, 36, 43–98.
Matsu’ura, R. S. (1986). Precursory quiescence and recovery of aftershock activity before some large
aftershocks, Bulletin of the Earthquake Research Institute, The University of Tokyo, 61, 1–65.
Matsu’ura, R. S. (1991). Case 1: Precursory quiescence and recovery of aftershock activity before
some large aftershocks, Evaluation of Proposed Earthquake Precursors (ed. M. Wyss), 8–11,
American Geophysical Union, Washington DC, http://dx.doi.org/10.1029/SP032p0008.
Mogi, K. (1968). Migration of seismic activity, Bulletin of the Earthquake Research Institute, The
University of Tokyo, 46, 53–74.
Mogi, K. (1973). Relationship between shallow and deep seismicity in the western Paciﬁc region,
Tectonophysics, 17, 1–22.
Mohler, G. O., Short, M. B., Brantingham, P. J., Schoenberg, F. P. and Tita, G. E. (2011). Self-exciting
point process modeling of crime, Journal of the American Statistical Association, 106(493),
100–108.
Musmeci, F. and Vere-Jones, D. (1992). A space-time clustering model for historical earthquakes,
Annals of the Institute of Statistical Mathematics, 44, 1–11, doi:10.1007/BF00048666.
Neyman, J. E. and Scott, E. L. (1953). Frequency of separation and interlocking of clusters of galaxies,
Proceedings of the National Academy of Sciences of the United States of America, 39, 737–743.
Neyman, J. E. and Scott, E. L. (1958). A statistical approach to problems of cosmology, Journal of theRoyal Statistical Society, Series B (Methodological), 20, 1–43.
Nishikawa, T. and Ide, S. (2017). Detection of earthquake swarms at subduction zones globally: Insights
into tectonic controls on swarm activity, Journal of Geophysical Research, 122, 5325–5343,
doi:10.1002/2017JB014188.
Nishikawa, T. and Ide, S. (2018). Recurring slow slip events and earthquake nucleation in the source
region of the M7 Ibaraki-Oki earthquakes revealed by earthquake swarm and foreshock activity,
Journal of Geophysical Research, 123, 7950–7968, doi:10.1029/2018JB015642.
Nishikawa, T., Matsuzawa, T., Ohta, K., Uchida, N., Nishimura, T. and Ide, S. (2019). The slow earthquake spectrum in the Japan Trench illuminated by the S-net seaﬂoor observatories, Science,
doi:10.1126/science.aax5618.
Nishizawa, O., Lei, X. and Nagao, T. (1994). Hazard Function Analysis of Seismo-electric Signals in
Greece, Terra Scientic Publishing Company, Tokyo.

統計地震学における ETAS モデル — その進展とホークス型モデル

159

野村俊一, 尾形良彦 (2020). 多様な予測方式に対識別モデルする前震とその予測性能評価, 地震予知連絡会会報, 103 (12-8), 361–366, https://cais.gsi.go.jp/YOCHIREN/report/kaihou103/12_08.pdf.
Ogata, Y. (1978). The asymptotic behavior of maximum likelihood estimators for stationary point
processes, Annals of the Institute of Statistical Mathematics, 30, 243–261.
尾形良彦 (1981). 事象発生の因果解析 — 地震の地域的関連性を測る, 統計モデル：モデル構成の新しい, 数理科学, 3 月号, 30–36.
波（赤池弘次 編）
Ogata, Y. (1981). On Lewis’ simulation method for point processes, IEEE Transaction on Information
Theory, IT-27, 23–31.
Ogata, Y. (1983a). Estimation of the parameters in the modiﬁed Omori formula for aftershock frequencies by the maximum likelihood procedure, Journal of Physics of the Earth, 31, 115–124.
Ogata, Y. (1983b). Likelihood analysis of point processes and its application to seismological data,
Bulletin of the International Statistical Institute, 50, 943–961.
Ogata. Y. (1985). Statistical models for earthquake occurrences and residual analysis for point processes, Research Memorandum, No. 288, The Institute of Statistical Mathematics, Tokyo.
Ogata, Y. (1988). Statistical models for earthquake occurrences and residual analysis for point processes, Journal of the American Statistical Association, 83(401), 9–27.
Ogata, Y. (1989). Statistical model for standard seismicity and detection of anomalies by residual
analysis, Tectonophysics, 169(1-3), 159–174.
Ogata, Y. (1992). Detection of precursory seismic quiescence before major earthquakes through astatistical model, Journal of Geophysical Research, 97, 19845–19871.
Ogata, Y. (1998). Space-time point-process models for earthquake occurrences, Annals of the Instituteof Statistical Mathematics, 50(2), 379–402.
Ogata, Y. (2001). Exploratory analysis of earthquake clusters by likelihood-based trigger models,Festscrift Volume for Professor Vere-Jones, Journal of Applied Probability, 38A, 202–212.
Ogata, Y. (2004). Space-time model for regional seismicity and detection of crustal stress changes,
Journal of Geophysical Research, 109(B3), B03308, doi:10.1029/2003JB002621.
Ogata, Y. (2006). Statistical analysis of seismicity - updated version (SASeis2006), Computer Science
Monographs, No.33, The Institute of Statistical Mathematics, Tokyo, https://www.ism.ac.jp/
editsec/csm/.
Ogata, Y. (2010). Space-time heterogeneity in aftershock activity, Geophysical Journal International,
181(3), 1575–1592, doi:10.1111/j.1365-246X.2010.04542.x.
尾形良彦 (2020). 階層的時空間 ETAS モデルによる短期・中期予測, 地震予知連絡会会報, 103(12-13),
385–387, https://cais.gsi.go.jp/YOCHIREN/report/kaihou103/12_13.pdf.
Ogata, Y. and Akaike, H. (1982). On linear intensity models for mixed doubly stochastic Poisson andself-exciting point process, Journal of the Royal Statistical Society, 44, 102–107.
Ogata, Y. and Katsura, K. (1986). Point-process models with linearly parameterized intensity for
application to earthquake data, Journal of Applied Probability, 23A, 291–310.
Ogata, Y. and Katsura, K. (2014). Comparing foreshock characteristics and foreshock forecasting inobserved and simulated earthquake catalogs, Journal of Geophysical Research: Solid Earth,
119, 8457–8477.
Ogata, Y. and Omi, T. (2020). Statistical monitoring and early forecasting of the earthquake sequence: Case studies after the 2019 M6.4 Searles Valley earthquake, California, Bulletin of the
Seismological Society of America, 110, 1781–1798.
Ogata, Y. and Shimazaki, K. (1984). Transition from aftershock to normal activity, Bulletin of the
Seismological Society of America, 74(5), 1757–1765.
Ogata, Y. and Zhuang, J. (2006). Space-time ETAS models and an improved extension, Tectonophysics,
413(1-2), 13–23.

160

統計数理

第 69 巻 第 2 号 2021

Ogata, Y., Akaike, H. and Katsura, K. (1982). The application of linear intensity models to the investigation of causal relations between a point process and another stochastic process, Annals of
the Institute of Statistical Mathematics, 34(2), 373–387.
Ogata, Y., Jones, L. M. and Toda, S. (2003a). When and where the aftershock activity was depressed:
Contrasting decay patterns of the proximate large earthquakes in southern California, Journalof Geophysical Research, 108, B62318.
Ogata, Y., Katsura, K. and Tanemura, M. (2003b). Modelling heterogeneous space-time occurrencesof earthquakes and its residual analysis, Journal of the Royal Statistical Society: Series C(Applied Statistics), 52, 499–509.
Ogata, Y., Katsura, K., Tsuruoka, H. and Hirata, N. (2018). Exploring magnitude forecasting of thenext earthquake, Seismological Research Letters, 89, 1298–1304.
Ogata, Y., Katsura, K., Tanemura, M., Harte, D. and Zhuang, J. (2021). Hierarchical space-time
point-process models (HIST-PPM): Software documentation, Computer Science Monographs,
No. 35, https://www.ism.ac.jp/editsec/csm/.
Omori, F. (1894). On the aftershocks of earthquakes, Journal of the College of Science, Imperial
University of Tokyo, 7, 111–200.
Peng, R. D., Schoenberg, F. P. and Woods, J. A. (2005). A space-time conditional intensity model forevaluating a wildre hazard index, Journal of the American Statistical Association, 100(469),
26–35.
Rathbun, S. L. (1993). Modeling marked spatio-temporal point patterns, Bulletin of the International
Statistical Institute, 55, Book 2, 379–396.
Reasenberg, P. A. and Jones, L. M. (1989). Earthquake hazard after a mainshock in California, Science,
243, 1173–1176.
Reasenberg, P. A. and Jones, L. M. (1994). Earthquake aftershocks: Update, Science, 265, 1251–1252.
Saichev, A. and Sornette, D. (2005a). Distribution of the largest aftershocks in branching models of
triggered seismicity: Theory of the universal Båth law, Physical Review E, 71(5), 056127.
Saichev, A. and Sornette, D. (2005b). Vere-Jones’ self-similar branching model, Physical Review E,
72(5), 056122.
Saichev, A. and Sornette, D. (2007). Theory of earthquake recurrence times, Journal of the Geophysics
Research, 112, B04313.
Schorlemmer, D., Werner, M., Marzocchi, W., Jordan, T., Ogata, Y., Jackson, D., Mak, S., Rhoades, D.,
Gerstenberger, M., Hirata, N., Liukis, M., Maechling, P., Strader, A., Taroni, M., Wiemer, S.,
Zechar, J. and Zhuang, J. (2018). The collaboratory for the study of earthquake predictability:
Achievements and priorities, Seismological Research Letters, 89(4), 1305–1313.
Tench, S., Fry, H. and Gill, P. (2016). Spatio-temporal patterns of IED usage by the Provisional Irish
Republican Army, European Journal of Applied Mathematics, 27(3), 377–402.
Truccolo, W., Eden, U. T., Fellows, M. R., Donoghue, J. P. and Brown, E. N. (2005). A point process
framework for relating neural spiking activity to spiking history, neural ensemble and extrinsic
covariate eects, Journal of Neurophysiology, 93, 1074–1089.
Turcotte, D., Holliday, J. and Rundle, J. (2007). BASS, an alternative to ETAS, Geophysical Research
Letters, 34(12), L12303, doi:10.1029/2007GL029696.
Ueda, T. and Kato, A. (2019). Seasonal variations in crustal seismicity in San-in District, southwest
Japan, Geophysical Research Letters, 46(6), 3172–3179, https://doi.org/10.1029/2018GL081789.

宇津徳治 (1957). 地震のマグニチュードと余震の起りかた, 地震 第 2 輯, 10, 35–45, doi:10.4294/
zisin1948.10.1_35.
Utsu, T. (1961). A statistical study on the occurrence of aftershocks, The Geophysical Magazine, 30,
521–605.

統計地震学における ETAS モデル — その進展とホークス型モデル

161

Utsu, T. (1962). On the nature of three Alaskan aftershock sequences of 1957 and 1958, Bulletin of
the Seismological Society of America, 52(2), 279–297.
Utsu, T. (1970a). Aftershocks and earthquake statistics (I): Some parameters which characterize anaftershock sequence and their interrelations, Journal of the Faculty of Science, Hokkaido University, Series 7, 3, 129–195.
Utsu, T. (1970b). Aftershocks and earthquake statistics (II): Further investigation of aftershocks andother earthquakes sequence based on a new classiﬁcation of earthquake sequences, Journal of
the Faculty of Science, Hokkaido University, Series 7, 3, 379–441.
Utsu, T. (1971). Aftershock and earthquake statistic (III): Analyses of the distribution of earthquakes
in magnitude, time and space with special consideration to clustering characteristics of earthquake occurrence (1), Journal of the Faculty of Science, Hokkaido University, Series 7, 3,
379–441.

宇津徳治 (1975). 関東地方の地震と飛騨地方のやや深発地震の相関について, 地震 第 2 輯, 28, 303–311.
Utsu, T. and Ogata, Y. (1997). Statistical analysis of seismicity, Algorithms for Earthquake Statistics
and Prediction (eds. J. Healy, V. Keilis-Borok and W. Lee), International Association of Seismology and Physics of the Earth’s Interior (IASPEI) Library 6, 13–94, IASPEI, Menlo Park,
California.
宇津徳治, 関彰 (1955). 余震区域の面積と本震のエネルギーとの関係, 地震 第 2 輯, 7, 233–240.
Utsu, T., Ogata, Y. and Matsu’ura, R. S. (1995). The centenary of the Omori formula for a decay lawof aftershock activity, Journal of Physics of the Earth, 43(1), 1–33.
van Stiphout, T., Zhuang, J. and Marsan, D. (2012). Seismicity declustering, Community Online Resource for Statistical Seismicity Analysis, doi:10.5078/corssa-52382934, available at
http://www.corssa.org.
Veen, A. and Schoenberg, F. P. (2008). Estimation of space-time branching process models in seismology using an EM-type algorithm, Journal of the American Statistical Association, 103(482),
614–624.
Vere-Jones, D. (1970). Stochastic models for earthquake occurrence, Journal of the Royal StatisticalSociety, Series B (Methodological), 32(1), 1–62 (with discussion).
Vere-Jones, D. (1973). The statistical estimation of earthquake risk, New Zealand Statistician, 8, 7–16.
Vere-Jones, D. (1975). Stochastic models for earthquake sequences, Geophysical Journal of the Royal
Astronomy Society, 42(2), 811–826.
Vere-Jones, D. (1978). Earthquake prediction-a statistician’s view, Journal of Physics of Earth, 26,
129–146.
Vere-Jones, D. (2005). A class of self-similar random measure, Advances in Applied Probability, 37(4),908–914.
Vere-Jones, D. and Davies, R. B. (1966). A statistical survey of earthquakes in the main seismic regionof New Zealand, Part 2, Time series analyses, New Zealand Journal of Geology and Geophysics,
9, 251–284.
Vere-Jones, D. and Ozaki, T. (1982). Some examples of statistical estimation applied to earthquake
data, Annals of the Institute of Statistical Mathematics, 34, 189–207.
Vere-Jones, D. and Zhuang, J. (2008). Distribution of the largest event in the critical epidemic-typeaftershock-sequence model, Physical Review E, 78(4), 047102.
Wang, Q., Jackson, D. D. and Zhuang, J. (2010). Are spontaneous earthquakes stationary in California?, Journal of Geophysical Research: Solid Earth, 115, B08310.
Whittle, P. (1962). Gaussian estimation in stationary time series, Bulletin of the International Statistical Institute, 39, 105–129.
Zhuang, J. (2006). Second-order residual analysis of spatiotemporal point processes and applications in

162

統計数理

第 69 巻 第 2 号 2021

model evaluation, Journal of the Royal Statistical Society: Series B (Statistical Methodology),
68(4), 635–653.
Zhuang, J. (2015). Weighted likelihood estimators for point processes, Spatial Statistics, 14, 166–178.
Zhuang, J. and Mateu, J. (2019). A semiparametric spatiotemporal Hawkes-type point process modelwith periodic background for crime data, Journal of the Royal Statistical Society: Series A(Statistics in Society), 189, 919–942.
Zhuang, J. and Ogata, Y. (2006). Properties of the probability distribution associated with the largest
event in an earthquake cluster and their implications to foreshocks, Physical Review E, 73,046134.
Zhuang, J., Ogata, Y. and Vere-Jones, D. (2002). Stochastic declustering of space-time earthquake
occurrences, Journal of the American Statistical Association, 97(3), 369–380.
Zhuang, J., Ogata, Y. and Vere-Jones, D. (2004). Analyzing earthquake clustering features by using
stochastic reconstruction, Journal of Geophysical Research, 109(3), B05301.
Zhuang, J., Chang, C.-P., Ogata, Y. and Chen, Y.-I. (2005a). A study on the background and clustering seismicity in the Taiwan region by using a point process model, Journal of GeophysicalResearch, 110, B05S13.
Zhuang, J., Vere-Jones, D., Guan, H., Ogata, Y. and Ma, L. (2005b). Preliminary analysis of observations on the ultra-low frequency electric eld in a region around Beijing, Pure and AppliedGeophysics, 162, 1367–1396.
Zhuang, J., Christophersen, A., Savage, M. K., Vere-Jones, D., Ogata, Y. and Jackson, D. D. (2008).
Dierences between spontaneous and triggered earthquakes: Their inuences on foreshock probabilities, Journal Geophysical Research, 113(B12), B11302, doi:10.1029/2008JB005579.
Zhuang, J., Ogata, Y., Vere-Jones, D., Ma, L. and Guan, H. (2013a). Statistical modeling of earthquake
occurrences based on external geophysical observations: with an illustrative application to theultra-low frequency ground electric signals observed in the Beijing region, Imaging, Modelingand Assimilation in Seismology (ed. Y. Li), Volume II, 351–377, De Gruyter, Berlin.
Zhuang, J., Werner, M. J. and Harte, D. S. (2013b). Stability of earthquake clustering models: Criticality and branching ratios, Physical Review E, 88, 062109.
Zhuang, J., Ogata, Y. and Wang, T. (2017). Data completeness of the Kumamoto earthquake sequencein the JMA catalog and its inﬂuence on the estimation of the ETAS parameters, Earth, Planets
and Space, 69, 36, doi:10.1186/s40623-017-0614-6 [open access].
Zhuang, J., Murru, M., Falcone, G. and Guo, Y. (2019). An extensive study of clustering features ofseismicity in Italy from 2005 to 2016, Geophysical Journal International, 216(1), 302–318.
Zipkin, J. R., Schoenberg, F. P., Coronges, K. and Bertozzi, A. L. (2015). Point-process models of social
network interactions: Parameter estimation and missing data recovery, European Journal ofApplied Mathematics, 27(3), 502–529, doi:10.1017/S0956792515000492.

Proceedings of the Institute of Statistical Mathematics Vol. 69, No. 2, 145–163 (2021)

163

The ETAS Model in Statistical Seismology: Its History, Recent
Developments, and Influences on General Hawkes Processes
Jiancang Zhuang1 and Yosihiko Ogata2
1 The Institute of Statistical Mathematics

2 The Institute of Statistical Mathematics, Professor Emeritus

We review the history of the ETAS model, a model for standard seismicity, in statistical seismology. In addition, we review new developments in theory, methods, andapplications of the model, as well as its influences on the developments of general Hawkes
processes.

Key words: ETAS model, Hawkes process, statistical seismology, point processes.