論文要約：Phase Neural Operator for Multi-Station Picking of Seismic Arrivals

https://agupubs.onlinelibrary.wiley.com/doi/full/10.1029/2023GL106434

ライセンス　CC BY-NC 4.0 Deed | Attribution-NonCommercial 4.0 International | Creative Commons

Abstract

目的: 従来の地震波到達時間測定に基づくアプローチを改善し、地震ネットワーク全体を同時に分析することで、地震波の位相をより正確に識別する新しいアルゴリズムを開発する。
手法: 「Neural Operator」という最近開発された機械学習パラダイムに基づいた「Phase Neural Operator」と呼ばれるモデルを導入。このモデルは、空間的および時間的な文脈情報を活用して、任意の地震ネットワークのジオメトリに対して位相を同時にピックアップする。
結果: このアプローチは、多くの地震を検出し、より多くの位相到達をピックアップすることができ、測定精度も大幅に向上するなど、主要なベースラインアルゴリズムよりも優れた性能を示した。
結論: 世界中で収集されている大量の地震データセットをフル活用することの潜在的な利点を示す一例として、人工知能の領域全体で見られる類似のトレンドに従っている。

Introduction

地震学において、地震波の到達時間を特定し、地震を識別する基本的な作業が地震相の検出とピッキングです。
歴史的には、人間の分析者が手動で地震信号と地震相の到達時間をラベリングしていましたが、この作業は主観的で時間がかかり、誤りやすいものでした。
精度の高い自動地震検出方法が開発されてきましたが、ほとんどの機械学習検出アルゴリズムは単一ステーション検出戦略を使用しており、弱い振幅のイベントを検出できない、または局所的なノイズ信号を誤って検出することがあります。
最近のディープラーニングの進歩は自動フェーズピッキングアルゴリズムの精度と効率を大幅に向上させましたが、単一ステーションのニューラルフェーズピッカーの性能向上は急速に頭打ちになっています。
本論文では、Phase Neural Operator（PhaseNO）と呼ばれるアルゴリズムを紹介し、ネットワーク全体での地震検出とフェーズピッキングを目的とした一般的なアプローチを提案します。PhaseNOは、任意のジオメトリを持つ複数のステーションで異なるフェーズの到達時間を正確に測定できるように、地震波場の無限次元関数表現を学習します。

Method: Phase Neural Operator

PhaseNOは、ネットワーク全体でのフェーズピッキングのためのオペレーター学習モデルであり、有界な物理領域上の無限次元関数空間間のオペレーターを学習するよう設計されています。入力関数は空間と時間の任意の点の集合で観測された地震波場$f(x,y,t)$であり、出力関数は各点$x,y,t$でのP波とS波の到着可能性を示す確率マスク$g(x,y,t)$です。
Neural Operatorは離散化不変性を持つため、入力と出力関数を毎回異なる（任意の）メッシュで離散化しても、モデルを再学習する必要がなく、入力物理領域内の任意の点で評価することができます。
PhaseNOは、時間情報にはFourier Neural Operator（FNO）層を、空間情報にはGraph Neural Operator（GNO）を組み合わせて、地震ネットワークデータの数学的構造を自然に扱います。モデルはFNOとGNOを順番に接続し、複数回繰り返すオペレーターレイヤーのブロックで構成され、地震ネットワーク内のすべてのステーション間で十分な時空間情報の通信と交換を可能にします。

スキップ接続を使用してブロックを接続し、U字型アーキテクチャを形成します。スキップ接続は、深い層を通過せずにモデルの左部分のFNO結果を右部分のGNO結果と直接結合し、収束を改善し、より深く、過パラメータ化されたモデルを可能にします。

Results

Performance Evaluation

この研究では、PhaseNOの性能をEQTransformer（Mousavi et al., 2020）、PhaseNet（W. Zhu & Beroza, 2018）、およびEdgePhase（Feng et al., 2022）という3つの主要なベースラインモデルとの比較によって評価されている。PhaseNOは、1984年から2019年にかけての北カリフォルニア地震データセンター（NCEDC）からの地震データセットでトレーニングされ、PhaseNetと同じトレーニングデータセットが使用されている。2020年の期間にわたるサンプルテストデータセットで、PhaseNOと各ベースラインモデルが評価された。このテストデータセットには、5,769件のイベントに対する43,700のP/Sピックが含まれている。
各サンプルの時間窓は、それぞれの事前トレーニングされたモデルに基づいて選択されており、PhaseNOとPhaseNetには30秒、EQTransformerとEdgePhaseには60秒が設定されている。ピックの位置は、時間窓の中央30秒にランダムに配置されている。全てのモデルにおいて、P波とS波のピックは、予測された確率分布のピークから、事前に定められた閾値を設定することによって決定されている。各モデルは、最適なF1スコアを最大化する閾値を使用し、最良の条件下で比較されている。
PhaseNOは、P波とS波の両方において最も高いF1スコア（それぞれ0.99および0.98）を達成している。これは、テストされた全てのモデルの中で最も高い最適閾値（P波に対して0.70、S波に対して0.65）を持つことに加えて、PhaseNOが他の方法よりも地震到達を検出し、ピックする際の信頼度が高いことを示している。真のピックが利用できない場合、PhaseNOは他の方法よりも偽検出を最小限に抑え、同じ事前に定められた閾値でより多くのピックを提供することができる。

PhaseNOは、ほぼ全ての信号対雑音比（SNR）レベルで、他のディープラーニングモデルと比較してより多くの真の陽性を検出し、偽の陰性と偽の陽性のピックを少なくしている。さらにピックを生成するにもかかわらず、PhaseNOはP波とS波の両方で最小の平均絶対誤差を達成している。全体として、PhaseNOは6つの指標すべてにおいて最高の性能を発揮しており、わずかな例外が1つある。PhaseNOのP波残差の標準偏差は、PhaseNetよりも0.01秒（1タイムステップ）大きい。PhaseNOによって新たに検出されたフェーズは、信号対雑音レベルが低いため、より挑戦的なケースである可能性が高く、その結果、標準偏差が若干増加する可能性がある。
複数の代表的なイベントに対する各ニューラルフェーズピッカーの予測確率分布を比較しました。PhaseNOは、異なるイベントの大きさ、計測器の種類、波形の形状に対して非常によく機能します。

PhaseNetはいくつかの偽陽性ピックを生成しますが、これらは複数の受振器を使う方法（PhaseNOおよびEdgePhase）によって除去されます。しかし、EdgePhaseは多くの偽陰性も生成します。
時間と空間の情報を複数回入れ替えることで、PhaseNOは偽のピックを効果的に防ぎながら、真のピックの検出能力を向上させます。
PhaseNOは、他の受振器からの情報を利用して、低SNRの波形上でピックを成功させます。
S波は通常、P波のコーダに存在し、検出がより困難です。そのため、S波ではP波よりも人間のアナリストによるラベリングエラーが予想されます。

Application to the 2019 Ridgecrest Earthquake Sequence

2019年のリッジクレスト地震シーケンスにおけるPhaseNOの検出性能と汎化能力がテストされました。
南カリフォルニア地震ネットワーク（SCSN）の20ステーションにおけるEH、HH、HNセンサーの連続波形データが、2019年7月4日15:00:00から7月10日00:00:00までの期間でダウンロードされました。これは合計36の異なるセンサーに相当します。
波形データは100Hzのサンプリングレートで1時間ごとのストリームに分割されています。このデータセットは、多数のイベントが重なっているため、より難しいデータセットであるとされています。
計算量を小さくするた、、地震ネットワーク全体が2つに分割され、データの各時間ごとに2つのグラフが構築されました。36のノードが18のノードを含む2つのグラフにランダムに分割されています。
PhaseNOは一度に18ノード上で地震相の確率を予測しています。連続データは30秒の時間窓に切り出され、10秒のオーバーラップを持ち、18ノードにおいて1時間のデータに対して180の予測がなされました。
PhaseNOはPhaseNetよりも多くの意味のある移動パターンを持つピックを検出しています。これは複数のステーション間で学習された波形の一貫性によるものです。

予測後、P波とS波の両方について0.3の閾値を使用してフェーズピックを決定しました。PhaseNOは693,266のP波と686,629のS波の到着時間を検出しましたが、PhaseNetは同じ閾値と同じステーションで542,793のP波と572,991のS波の到着時間を見つけました。
検出されたピックの精度は、SCSNから手動でレビューされたピックとの到着時間を比較することによって評価されました。SCSNとPhaseNOの間のピック残差の標準偏差は、P波で0.10秒、S波で0.14秒でした。
これらの標準偏差は、PhaseNetのもの（P波で0.08秒、S波で0.13秒）よりも若干高かったです。新たに検出されたピックは信頼性が低い困難なケースであるため、PhaseNOがより大きな旅行時間差を示すのは合理的です。
候補となるフェーズ検出をGaMMAを使用してフェーズ関連付けによりイベントに変換しました。イベントごとに最低17のピックを設定して、低品質の関連付けを除外しました。これにより、PhaseNetはイベントごとに37.54のピックを持つ21,748イベントを検出し、PhaseNOはイベントごとに39.37のピックを持つ26,176イベントを検出しました。
GaMMAによって計算されたPhaseNOとPhaseNetによって検出されたイベントのマグニチュードは、両者ともに線形のガッテンベルグ・リヒター分布を示しました。結果として、Ross et al.（2019）によるテンプレートマッチングカタログよりもマイクロ地震が少なくなっています。

イベントの位置はGaMMAによって決定され、PhaseNOとPhaseNetのカタログ間で一般的に一致しており、PhaseNOによる追加イベントが実際の地震の合理的な検出であることを確認しています。
手動でレビューされたSCSNカタログをベースラインとして扱い、どれだけの地震が正確に回収されたかを評価しました。2つのイベントが3秒以内に発生した場合、それらを一致しているとみなします。
この基準により、Shelly、Ross et al.、PhaseNetはそれぞれ約81%、86%、88%のイベントに一致しました。比較して、PhaseNOのカタログは26,176イベントの合計でSCSNカタログの約94%（11,389のうち10,673）に一致し、追加イベントを含め、PhaseNOの最も高いリコールスコアを示しました。
PhaseNOは、PhaseNet、SCSN、Shellyのテンプレートマッチングカタログよりも一貫して多くのイベントを検出し、もう一つの詳細なテンプレートマッチングカタログ（Ross et al., 2019）に報告された地震の数に近づきました。
PhaseNOは、Mw 7.1の本震が発生したときに最も多くのイベントを検出し、その後の地震活動率が徐々に減少することで、複雑な地震シーケンスを照らす方法の力を示しました。
私たちのカタログは、ステーションの数と関連付けアルゴリズムの点でSCSNおよびテンプレートマッチングカタログと異なります。しかし、PhaseNOとPhaseNetからのピックは正確に同じステーションで検出され、GaMMAで関連付けられ、最も公平な比較を提供します。
2つのポストプロセッシングハイパーパラメーター、フェーズピッキングの閾値とイベントに関連付けられるピックの最小数が、カタログ内の地震の総数を制御します。より低い閾値とより小さい関連付け最小限は、より多くのイベントを提供しますが、偽陽性イベントが増える可能性があります。
PhaseNOは、同じハイパーパラメーターを使用してPhaseNetよりも一貫して多くのイベントを検出し、フェーズピッキングのために時間的情報に加えて空間情報を活用することの重要性を指摘しています。

Discussion and Conclusions

固定されたモデルアーキテクチャを持つPhaseNOは、任意のジオメトリを持つ地震ネットワークを処理できます。これは、北カリフォルニア地震ネットワークでのトレーニングと、再トレーニングなしで南カリフォルニア地震ネットワークでのモデル評価によって示されました。
PhaseNOは、一つのステーションの地震計測値ではなく、全地震ネットワークからの文脈を使用することで、人間の学習と意思決定を模倣します。これにより、特にP波のコーダに隠れがちなS波を含む低SNRデータのフェーズピッキングが大幅に改善されます。
PhaseNOは、時間情報からフェーズを識別する独特の能力を持ちます。よく知られているトランスフォーマーアーキテクチャは、Neural Operatorsの特別なケースと見なすことができます。
EdgePhaseと比較して、PhaseNOはGNOレイヤーを複数使用し、ネットワーク全体でのカーネル統合を可能にするNeural Operatorを使用して、豊かな空間特徴を抽出します。
PhaseNOの主な制限点は、一つの予測におけるステーションの数に対するメモリ使用量の依存です。空間情報はグラフ内のすべてのノードペア間で交換されるため、計算コストはノードの数に対して二次的にスケールし、複雑さはO(n^2)です。
広範囲にわたる地震ネットワークをカバーする場合、k-meansクラスタリングに基づいてステーションを選択することができます。この方法により、予測手順を大幅に加速し、特に多くのステーションが存在し、計算リソースが限られている場合にメモリ使用量を節約できます。

Supplemental

Neural Operators

ニューラルオペレーターは、無限次元関数空間間のマッピングを行うニューラルネットワークの一般化です。この新しいモデルクラスは、オペレーターに対する普遍近似定理を満たすことが証明されています。
ここでは、波動場からフェーズピックへのマップを学習するための新しいアーキテクチャを提案します。ニューラルオペレーターは一般に、入力関数をより大きな共域にマッピングするリフティングオペレーター(P)から始まります。
これらの関数は非線形カーネル統合オペレーターで反復的に操作され、最終的には隠れた表現を出力関数にマッピングする射影オペレーター(Q)を通過します。PとQは完全に接続されたニューラルネットワークでパラメータ化され、物理領域上で点ごとに作用します。
反復カーネル統合の基本式は、線形積分オペレーターと非線形活性化関数の組み合わせです。各積分オペレーターは以下の形式を持ちます：

ここで、vとuはそれぞれ中間入力と出力関数であり、κはカーネル関数です。オペレーターへの入力としてを定義します。カーネルをパラメータ化する方法は複数あります。
地震ネットワークによって記録された地震波形を入力関数fとして扱い、時間領域では規則的なメッシュ、空間領域では不規則なメッシュで離散化します。私たちのアーキテクチャでは、空間と時間でカーネル関数を別々に計算します。

Fourier Neural Operators

時間領域の規則的なメッシュに対して、カーネルをフーリエ空間でパラメータ化し、高速フーリエ変換を使用してカーネル積分オペレーターを計算し、ほぼ線形の複雑さで効率的な計算を実現します。
畳み込み定理から、次の式が得られます
ここでFとF^{−1}はフーリエ変換とその逆変換を示します。Rϕは、ϕによってパラメータ化されたκのフーリエ変換です。(時間領域での畳み込み操作は周波数領域での乗算）
各層で局所的に作用する活性化関数σを用いて、単一のFNO層の更新は次のようになります:、ここでWは局所的な線形オペレーターです。
実際には、フーリエ級数を最大モード数で切り捨て、Rをいくつかの低いモードでパラメータ化します。入力vから始まり、1つのFNO層には2つの並列分岐が含まれます：一方の分岐はフーリエ空間でカーネルを計算し、グローバル統合を実行し、もう一方は入力に対して点ごとの線形変換Wを適用します。2つの分岐からの結果は、σを適用する前に加算されます。
PhaseNOは、U-NOアーキテクチャに類似した7つのFNO層を利用します。各FNO層のモード数は24、12、8、8、12、24、24です。
各ノードで離散化されたuの幅（チャネル数）は、Rϕの次元に応じて変化します。各FNO層で、離散化されたuは48×3000、96×750、192×200、96×750、48×3000、48×3000、48×3000の次元を持ちます。
すべてのFNO層は、最後の層を除いて、ガウス誤差線形単位による非線形性を含みます。最後の層では活性化関数は適用されません。簡単のため、最後のFNO層は図1には描かれていません。

Graph Neural Operators and the message passing framework

カーネル統合は、グラフニューラルネットワークの一般化されたメッセージパッシングによるメッセージの集約として見ることができます。空間領域の$f(x,y,t)$は地震ネットワークのジオメトリに基づいて離散化されるため、カーネルをGNOでパラメータ化し、メッセージパッシングフレームワークで実装します。
任意のジオメトリを持つ地震ネットワークをグラフとみなします。地震ネットワーク内の各ステーションはグラフのノードです。ノードの特徴v(x)が与えられたとき、ノードxiの値v(xi)を次の平均集約によって値u(xi)に更新します

ここでτとϕは多層パーセプトロン(MLP)のような異なる微分可能関数を示します。N(x_i)はxiの近傍の数を表します。
グローバルな依存関係を捉えるために、グラフ内の各ノードをグラフ内のすべてのノードに接続するグラフを構築し、総エッジ数を$n^2$にします。つまり、N(xi)は地震ネットワーク内のすべてのステーションを含み、自己ループも含むことを意味します。
エッジ特徴はe_ij=ϕ(v(x_i),v(x_j))によって計算されます。
は学習可能なパラメータのセットを持つ非線形関数です。4cニューロンを含む1つの隠れ層を持つMLPをϕとして選択します。
すべてのエッジ特徴（メッセージ）が利用可能になると、ターゲットノードx_iはすべてのメッセージを収集し、平均操作でそれらを集約します。最後に、v(x_i)と集約されたメッセージの連結を入力として使用して、別のMLPでx_iのノード特徴を更新します。
メッセージパッシングにより、隣接ノード間で情報の交換が可能になり、隣接ノードによって共有される関連信号が強化されます。

Training and test datasets

高度なディープラーニングモデルのアーキテクチャを完全に活用するためには、十分な品質と量を持つトレーニングデータセットが必要である。地震ネットワークの波動場特性を考慮に入れ、グラフ型サンプルに特化した効果的なデータ拡張戦略が考案されている。初めに、グラフレベルでイベントを積み重ね、異なるステーションでの到着パターンを各イベントに対して保持することを目指している。この方法により、グラフ内の異なるステーションの波形は異なる数のイベントで構成される場合がある。次に、地震はべき乗則に従い、ほとんどのイベントが小さく、全地震ネットワークではなく数ステーションでのみ観測される可能性がある。そのため、グラフ内のランダムな位置に仮想ステーションを追加し、ノイズ波形を用いてPhaseNOを正則化することが重要である。
トレーニングデータセットは、1984年から2019年までのNCEDCからの三成分地震波形と、STanford EArthquake Dataset（STEAD）からの三成分ノイズ波形で構成されている。地震データはイベントごとにダウンロードされ、人間のアナリストによってP波とS波の到着時間がピックされたステーションが含まれている。一部のセグメントが欠落している場合はゼロで埋められる。各成分は、平均を除去し、標準偏差で割ることによって正規化される。フェーズ到着をラベル付けするために、三角形の形状をした確率関数が使用される。手動でピックされたP/S到着での確率は1とラベル付けされ、手動ピックの前後で0に線形に減少する。各ピックに対して、0より大きい確率の持続時間は0.4秒で、最高確率は時間窓の中央に中心される。単一ステーションの地震計測値を1つのサンプルとして扱う代わりに、地震ネットワーク内のすべてのステーションでグラフを構築し、そのグラフを1つのサンプルとして使用する。
トレーニング中にデータ拡張が実行される。異なるステーションでの移動パターンを保持するために、個別にダウンロードされたイベントは以下のステップで積み重ねられる：
- イベントAを記録するすべてのステーションからランダムにステーションAが選択される。
- ステーションAによって記録されたすべてのイベントからランダムにイベントBが選択される。
- 2つのイベントの振幅にαとβ（0.1 < α < 0.9, 0.1 < β < 0.9, α + β = 1）の重みがランダムに割り当てられる。
- 2つのイベントを積み重ねるために、ランダムに時間シフトが選択される。
- 両方のイベントが同じステーションで記録されている場合、イベントAとイベントBが積み重ねられる。
- 1つのイベントのみを記録する他のステーションでの波形が保持される。
地震ネットワーク内では通常、複数のステーションが両方のイベントを記録している。イベントBはイベントAを記録した一つのステーションに基づいて選択されるが、より多くのイベントを積み重ねることが可能である。
トレーニングデータセットの約66%のサンプルに2つまたは3つのイベントが含まれている。最大16の仮想ステーションをランダムな位置に生成し、これらにノイズ波形を割り当てる。
ノイズデータはSTEADの235Kのノイズサンプルからランダムに選択される。地震波形のみを含むサンプルが6.25%であり、残りのサンプルには地震ネットワーク内の異なるステーションで記録された地震波形とノイズ波形が含まれている。
一つのグラフ型サンプル内の各ステーションのイベント数は0から3の範囲である。地震ネットワークには三成分と一成分の地震計が含まれている場合があり、ランダムに選択された複数のステーションを一成分ステーションとみなす。
各サンプルには異なる数のステーションが含まれ、計算コストを節約するために一つのグラフ内のステーション数は5以上32以下に保持される。全ステーションの30秒の波形をランダムな開始時間で切り出し、ウィンドウ内のフェーズの位置が変化する。
サンプリングレートは100Hzであり、入力波形と出力確率は各ステーションの各成分について3000のデータポイントを持つ。トレーニング用には合計57Kのグラフがある。
テストデータセットには、2020年のNCEDC地震データセットで構築された5769のサンプルが含まれており、データ拡張なしでトレーニングデータセットと同様の方法で波形が前処理される。テストデータセットの各サンプルには、複数のステーションで記録された一つのイベントのみが含まれている。