平成21年度年次報告

（1）実施機関名

名古屋大学

（2）研究課題（または観測項目）名

プレート境界地震のための地殻活動予測シミュレーション・データ同化システムの構築

（3）最も関連の深い建議の項目

• • 1.地震・火山現象予測のための観測研究の推進
  • （2）地震・火山現象に関する予測システムの構築
  • （2-1）地震発生予測システム
  • ア．地殻活動予測シミュレーションとデータ同化

（4）その他関連する建議の項目

• 1.地震・火山現象予測のための観測研究の推進
• （2）地震・火山現象に関する予測システムの構築
• （2-1）地震発生予測システム
• イ．地殻活動予測シミュレーションの高度化

（5）本課題の５か年の到達目標

本課題では，第２次計画までに開発された地殻活動予測シミュレーションモデルをさらに改良し発展させるとともに、地震活動、地殻変動等の観測データの情報をシミュレーションに取り入れるためのデータ同化システムのプロトタイプを開発し、プレート境界地震の発生履歴の再現やプレート境界の摩擦特性の推定等を行う。さらに、これらのシステムおよびモニタリングシステムを構成要素とする地震発生予測システムの全体像を設計し、プロトタイプシステムを実際のデータに適用してプレート境界におけるすべりの時空間発展や地震発生を予測するデータ同化・予測実験を行う。こうした研究を通して本格的な地震発生予測システムの構築に向けて準備を整える。

（6）本課題の５か年計画の概要

(a) 地殻活動予測シミュレーションモデルの開発研究および予測実験
日本列島全域を対象とした地震発生シミュレーションの為の数値計算コードの改良と最適化を行なう．また，プレート境界面の摩擦特性の推定に向けた地殻変動解析モデルの構築を進める。
また，沈み込むプレートの3次元形状を考慮し、規模依存の破壊エネルギー摩擦特性分布を用いた、半無限均質弾性媒質中における南海トラフおよび千島海溝巨大地震発生サイクルシミュレーションを行い、歴史地震に見られる、大きな発生間隔・規模・東西セグメントの破壊時間差を再現するモデル作成を行うとともに，シミュレーションコードの高速化を進める。さらに，すべり応答関数計算を高度化するため，GeoFEMコードを改良して均質粘弾性媒質モデル作成の効率化を図る。
一方，これまでのシミュレーションであまり考慮されていなかった．応力や摩擦パラメータ等の短波長不均一が地震サイクルや予測に及ぼす影響をシミュレーションにより調べる。
さらに，プレート境界と内陸地震の相互作用の理解と活動度予測のための、計算効率の高い不連続体セルモデルシミュレーションモデル構築において、成層粘弾性媒質中でのすべり応答関数を用いた粘弾性地震発生サイクルシミュレーションに取り掛かる。

平成21年度は，南海トラフ、日本海溝、千島海溝等のプレート境界を対象に、シミュレーションモデルを作成する。フォワード計算により、応力相互作用の評価を行う。
平成22年度は，不均質粘弾性構造、動的破壊過程、プレート境界地震と内陸地震の相互作用等を考慮して、モデルの高度化を進める。相互作用評価を継続する。
平成23年度は，モデル高度化を継続する。過去の大地震発生系列を再現するシミュレーションにより、摩擦特性等を推定する。近年の大地震の破壊過程、余効すべり過程のシミュレーションにより、摩擦特性等を推定する。
平成24年度は，シミュレーションによる摩擦特性等の推定の継続。プレート境界と内陸の相互作用のシミュレーション。摩擦パラメータ推定の不確定性が予測精度に及ぼす影響の評価を行う。
平成25年度は，これまでに推定された摩擦パラメータ等に基づき、予測シミュレーションを試行する。予測シミュレーションの誤差を評価する。シミュレーションに基づき、西南日本内陸の地震活動の変化について議論する。

(b) 地殻変動データを用いたデータ同化手法の開発と同化実験
プレート境界地震の余効すべりなどゆっくりとした断層運動を対象として，GPSデータなど地殻変動データから, 速度・状態依存摩擦構成則で用いられる摩擦パラメータや, すべり速度等の初期値を推定する手法を開発する．単純なバネ・ブロックモデルや複数のセルモデルを用いた研究を進めるとともに，２次元の連続体モデル（１次元断層）による解析手法について検討する．

平成21年度は，地殻活動予測システムの全体設計、パラメータ最適化手法、逐次データ同化手法の検討を行う。
平成22年度は，地殻活動予測システムの全体設計を完了し、パラメータ最適化手法、逐次データ同化手法の開発、試験を行う。
平成23年度は，データ同化の解析事例を増やし、手法の高度化を行う。
平成24年度は，予測シミュレーションと連携したデータ同化実験のためのシステム開発を行う。
平成25年度は，予測シミュレーションとデータ同化システムを連携させた予測システムのプロトタイプを稼働させる。

(c)地震活動データの活用手法の開発
繰り返し地震をもれなく抽出し、すべり量分布の推定の精度を上げるため、幅広い規模の地震に適用できる新たな繰り返し地震抽出基準の策定を行う。具体的には、釜石沖等のすでに知られていて、現在の波形相似性の基準では抽出されない中規模繰り返し地震を用い、波形の相似性からこれらの繰り返し地震を同定できる基準を調べる（解析周波数帯、ウインドウ等の工夫）。次に、この基準を用い、過去の波形記録を系統的に探索する。
関東地方及び東海・東南海地域の微小地震活動度の長期的変化を高信頼度で推定する。そのためには、時間的・空間的に均質な地震カタログを作成する必要がある。東京大学地震研究所の微小地震観測網で1994年から2007年までに観測されたM3以上の地震の震源とMの再決定を実施する。その際、使用する観測点は1994年当時の配置に固定し、新しい観測点を導入したことによる影響を除去する。（北大・地震研）
測地データに加えて小繰り返し地震のデータを用いて、余効すべりの時間発展を推定する手法を確立する。特に、異種データの重み付けや滑りの空間分布に対する penalty 項の与え方について種々検討を行い、最適なものを見いだす。

平成21年度は，波形の相関、詳細な震源決定や、発生間隔等をもとに、小繰り返し地震を抽出する最適な手法の開発を開発する。また、地震活動から応力変化を推定する手法を開発する。
平成22年度は，小繰り返し地震とGPSデータを用いてプレート間すべりの推定を行う手法を開発する。 また、地震活動データから応力変化の推定を試み、その有効性を検証する。
平成23年度は，開発した手法を用いて、プレート境界のすべりの時空間発展を推定する。また、推定された応力変化がプレート境界のすべりによるものと仮定して、すべりの時空間発展について知見を得る。
平成24年度は，地震活動データから推定されたすべり変化と数値シミュレーションを比較することにより、プレート境界の摩擦特性、応力について知見を得る。
平成25年度は，地震活動データを用いたデータ同化手法について検討する。

（7）平成21年度成果の概要

(a) 地殻活動予測シミュレーションモデルの開発研究および予測実験

・準静的地震発生サイクルシミュレーションコードの高速化（京大）
　地震発生サイクルシミュレーションでは、プレート境界をN個のセルに分割し、各セルでの運動方程式と摩擦則をカップルさせて、すべりの発展を追いかける。この際応力計算において、すべり応答関数行列とすべり（速度）ベクトル積を計算する必要があるが、その計算量および計算に必要なメモリー量はO(N*N)である。南海トラフなどの大領域における地震発生サイクルシミュレーションではNが大きく、計算時間および使用メモリーの大きな大規模計算になる。また摩擦パラメータの推定（データ同化）等には、摩擦構成則の強い非線形性の故、多くの繰り返し計算を必要とする。こういったことから、地震発生サイクルシミュレーション計算の高速化・省メモリー化が必須である。これまでの高速化の手法としてFFTを用いたものが挙げられるが、空間的対称性を仮定した周期的境界条件を用いていて、適応範囲が限られている。そこで、天体物理学におけるN体長距離相互作用の高速計算用に開発された高速多重極法（FMM）を用いて、無限均質弾性媒質中でのすべり応答関数行列とすべり(速度)ベクトル積の計算の高速化の検討を行なった。その結果、O(N)の計算時間およびメモリー量で計算でき、高速化が可能であると分かった(図１)。今後は、実際の地震発生サイクルシミュレーションに組み込み、計算を実行する予定である。さらに、半無限均質弾性媒質中での高速多重極展開に関する検討を行った。この場合、多重極に展開できない項が含まれており、さらなる近似を用いる必要が判明し、現状ではまだうまくいっていない。他の方法として、自由表面である地表変位も求めて、無限媒質における多重極展開を用いる方法を検討中である。

・不均質粘弾性媒質中でのすべり応答関数の高度化（東大地震研・京大）
　粘弾性不均質媒質中でのすべり応答関数計算において、不均質粘弾性媒質モデル有限要素（FEM）モデル作成を作成する必要がある。この効率化を図るため、これまで開発してきた有限要素法（GeoFEM）コードの改良を行った。GeoFEMでは、これまで、３次元６面体体要素しか取り扱えず、複雑なプレート形状のメッシュ作成には有効なメッシュジェネレータソフトが少なく、非常に時間がかかった。そこで、メッシュ作成が容易な３次元４面体要素を扱えるように改良を行った。また、２次元問題も扱えるように、２次元３角形要素も扱えるような改良も行った。

・プレート境界地震と内陸地震の相互作用を含む不連続セルモデルシミュレーション（京大）
プレート境界地震と内陸地震の相互作用の理解と活動度予測のための、計算効率の高い不連続体セルモデルシミュレーションモデル構築する必要がある。そのため、成層粘弾性媒質中での点震源に対するすべり応答関数を用いて、矩形断層セルに対するすべり応答関数を効率よく計算する手法を開発した。

(b) 地殻変動データを用いたデータ同化手法の開発と同化実験

・GPSデータに基づくすべり欠損分布の推定（名大）
　測地データインバージョン解析手法を西南日本のGPSデータに適用して，ユーラシア−フィリピン海プレート境界の固着−すべり状態を推定し，南海・東南海・東海地震の震源域のすべり遅れレートの詳細な分布を明らかにした（図２）。その結果、東海から紀伊半島、四国にかけてすべり欠損量が次第に大きくなり、また、すべり欠損域の西端が豊後水道にも及ぶことから、宝永地震のような連動型地震では、従来よりも震源域の想定を広げる必要性のあることが示唆される。

・ 測地データに基づくすべり欠損とアスペリティ分布に関する定量的検討（名大）
測地データの逆解析から推定されるすべり欠損分布はアスペリティ分布を反映したものであるが、両者の間にはプレート境界面の摩擦特性が介在しており、必ずしも一対一の対応関係ではない。そこで、これら２つの分布の関係を定量的に検証するため、摩擦特性（アスペリティ）の分布とプレート沈み込みを仮定して、プレート境界面のすべり欠損分布と、それにより生じる地表の地殻変動分布を理論計算し、この計算値に誤差を加えて合成したデータの逆解析から得られるすべり欠損分布の推定結果と、当初仮定したアスペリティ分布との比較を行った。その結果、すべり欠損分布から想定される地震モーメントの蓄積速度は、正しいすべり欠損分布を反映するものの、当初仮定したアスペリティの面積と比較すると最大３倍程度過大評価となる可能性があることが分かった（図３）。また、観測網とプレート境界面の位置関係によらず、推定されたすべり欠損分布の重心は、アスペリティ分布の重心よりも浅めに推定されることも分かった。逆に、測地データからすべり欠損が推定できない場所においては、実際にアスペリティが存在しない可能性が高い。以上の結果から、すべり欠損の結果に基づいて数値シミュレーションなどを行う際に注意すべき点が明らかとなった。

・アジョイント法によるデータ同化手法の開発（京大）
　京都大学では、地殻変動データから断層面上の断層パラメータを推定する手法として、アジョイント法を適用する方向で検討している。21年度は本格的なデータ同化を行う準備として、自由度の小さいセルモデルを用いて余効すべりを模したsynthetic testを行った。
　まず、2003年十勝沖地震をモデルケースとし、本震および余効すべり域を、それぞれ１つおよび２つの断層セルで近似し、それぞれに十勝沖地震に類似したイベントが発生するように摩擦パラメータの組（これが真の値となる）を与え、速度・状態依存摩擦構成則（状態の発展方程式を含む）と準静的な場合の運動方程式とを連立させて、各断層セルのすべり速度・応力・状態変数を、地震の数サイクルにわたって計算した。本実験では、断層セルのすべり速度が直接観測できるとし、上記によって得られたすべり速度に白色雑音を加えたものを人工観測データとし、これを基に余効すべり領域の断層セルの摩擦パラメータの推定を行った。
　アジョイント法では、観測データ（o）と、ある一組の摩擦パラメータの組（m）から数値モデルを通して計算されるデータ（c）の残差二乗和と、何らかの方法で得られた摩擦パラメータに関する先験情報（b）と上記 m との残差二乗和の線形和を評価関数（J(m)）とし、mをいろいろ変化させて評価関数が最小となるmの最適値をさがす。mとしては、摩擦パラメータ(a,b,L)だけでなく、シミュレーションの初期値（すべり速度、応力、状態変数から２つ）も推定でき、これは先行研究ではなされていない。
　アジョイント法を適用するにあたり、自由度の低い系での人工データを用いて、各摩擦パラメータに対する評価関数の振る舞いをグリッドサーチで調べた。その結果、例えばGPS日座標値のようなデータが長期間あってもLの値の拘束は困難だが、地震直後における高サンプリングのデータがあれば拘束が可能であることがわかった。このような考察は系の自由度が高い場合でも有効である。人工データに対して実際にアジョイント法を適用すると、摩擦パラメータの推定のみをアジョイント法で行う場合、シミュレーション変数の初期値が正しければ摩擦パラメータもほぼ正しい値が得られるが、間違った初期値を与えれば、摩擦パラメータの値も誤った結果が得られる。また、摩擦パラメータと初期値を同時に推定すれば、ともにほぼ正しい値が推定された。シミュレーション変数の初期値に関しては、すべり速度はインバージョン解析から推定できるが、もう一つの変数を先験的に知ることは難しい。したがって、本研究で示されたように、パラメータだけでなく初期値も同時に推定することが重要である。

・余効変動に基づく摩擦パラメータの推定（東大地震研）
摩擦法則を用いた地震サイクルの数値シミュレーションを行うためには，断層面上の摩擦特性を表わす摩擦パラメータを与える必要がある．速度・状態依存摩擦法則を用いたシミュレーションでは，L, as, bsの3つのパラメータを与える必要がある．これらのパラメータの値は室内実験で測定されているが，現実の断層に対する値が推定された例はほとんどない．例外として，いくつかの研究で，数カ月から数年間の余効変動を記録した1日毎のGPS時系列から(a-b)sの値のみが推定されてきた．ここでは，地震後数時間の余効変動が3つの摩擦パラメータの推定に有効であることを示し，これらのパラメータを推定するためのインバージョン手法を開発した．さらに，この手法を2003年十勝沖地震後5時間のGPS時系列に適用し，3つのパラメータを独立に推定した．フォワードモデルとしては，速度・状態依存摩擦法則に従う1自由度のばね・ブロックモデルを用いた．摩擦パラメータ推定のための逆問題はベイズ推定の枠組みで定式化し，逆問題の解である摩擦パラメータの事後確率分布をマルコフ連鎖モンテカルロ(MCMC)法により推定した．推定されたLは室内実験で測定された値より10-1000倍大きくなった(図4)．また，推定されたasと(a-b)sの値は，a及びa-bとして実験室での測定値，間隙水圧として静水圧を仮定して得られた値より1桁から2桁小さくなった．これは現実のプレート境界に対するa及びa-bの値が実験値よりも小さいか，間隙水圧が静水圧よりも高いことを示している．

(c)地震活動データの活用手法の開発

・シミュレーションによる応力状態の評価（東大地震研）
プレート内地震の応力降下量はプレート境界地震よりも大きいことが知られている．これについて定性的に説明を試みた研究はいくつかあるが，定量的な説明に成功した研究はない．ここでは，プレート内の断層とプレート境界とではすべり特性の不均一性の違いにより応力集中過程が異なると考えられることに着目して数値シミュレーションを行い，応力降下量の違いの定量的説明を試みた（Kato, 2009）．プレート内断層のモデルでは，地震性すべり域の周囲に強度無限大のバリアを仮定した．一方，プレート境界地震のモデルでは，地震発生領域（アスペリティ）の摩擦特性を速度弱化とし，その周囲では速度強化の摩擦特性であると仮定した．これらモデルに一定速度のせん断をかけると，プレート内地震のモデルではほぼ一様にせん断応力が増大するが，プレート境界地震のモデルでは周囲の速度強化域で発生する非地震性すべりのために，地震発生域の縁で応力集中が生じる（図５）．プレート境界地震モデルでは，応力集中のために早めに地震が発生するが，地震発生時にアスペリティ内部の応力は十分に大きくなっていないため，平均的な応力降下量は小さくなる．様々なパラメータで数値シミュレーションを行った結果，プレート内地震モデルの応力降下量は，プレート境界地震モデルよりも数十%程度大きくなった．

・相似地震解析によるプレート固着状態の検出（東大地震研）
日本列島全域に展開されているテレメータ地震観測点で観測された地震の波形データの蓄積を行い，プレート沈み込み境界域を中心に相似地震・小繰り返し地震の抽出を行った．その結果，根室半島沖から房総半島沖までの太平洋プレート沈み込み帯及び，日向灘沖から八重山諸島沖のフィリピン海プレート沈み込み帯において，比較的定常的な間隔で発生する小繰り返し地震を多数発見した．特に，奄美大島より南のフィリピン海プレートでは，海溝軸に近い深さ10-20kmと，プレート境界型地震発生域の深部端と考えられる深さ30-50km付近にそれぞれ，海溝軸に平行に分布する特徴が見られた．これらから推定されるすべり速度はおよそ6cm/yrから12cm/yrであり，南西部に行くほど速くなっていた．これはGPSによって推定されたプレート間の相対運動速度ともほぼ一致しており，この地域での固着が弱い傾向を示していると考えられる．

・小繰り返し地震抽出のための手法検討（東北大）
岩手県釜石沖のM4.8の繰り返し地震を含む地震クラスターには，大小様々な規模の地震が含まれ，複数の繰り返し地震グループも存在する。これまでに，このクラスターにおいて波形の相関を用いた数十mの精度での震源決定と，コーナー周波数に基づく10%程度の誤差でのすべり域（断層半径）の推定を行った。ここで得られた震源間距離，断層サイズデータを用い，小繰り返し地震を抽出するための最適な手法の検討を行った。
地震間の滑り域の重なり具合や規模の違いと波形のコヒーレンスの関係を調べた結果，2つの地震のうち小さい方の地震のコーナー周波数付近のコヒーレンスを用いると，繰り返し地震とそうでないものの区別が明瞭にできること分かった。この地震の規模に合った周波数帯域で解析を行う新しい手法により，これまでよりも正確に幅広い規模にわたって繰り返し地震を抽出することができると考えられる。

・均質な地震カタログに基づく地震活動評価（北大）
地震活動の長期変化から広域応力場の時間変化を検出するために、1994年1月1日から2007年12月31日までの間に関東地方および東海・東南海地域で発生したM3.5以上の地震の震源およびMを再決定した。再決定には、東京大学地震研究所が関東甲信越に展開している微小地震観測点、名古屋大学と東北大学が展開している観測点、および気象庁の東海沖と房総沖の海底地震計など、合計23か所を使用した。記録された約3500個のイベント波形ファイルのP波・S波到着時と最大振幅を、全て手動により注意深く再検測し、それらのデータを使用して震源とMを再決定した。選定した23か所の観測点は期間中に観測条件が変化していないこと、再検測は検測経験豊富な一人の作業員が行ったことなどから、作成された地震カタログは時間的・空間的に極めて均質であると考えられる。
　上記の地震カタログを解析した結果以下のことが分かった。
１．東海地方で発生した深さ20km以浅の地震70個の積算度数分布を見ると、2000年7月以降発生レートが増加し、その後徐々に低下、2005年頃に元のレートに戻った。この変化は東海スロースリップと同期しているように見える。
２．2005年1月から2007年1月にかけての2年間、茨城県南部から東京湾にかけて顕著な地震活動の活発化が、銚子沖では顕著な静穏化が観測された。活発化域と静穏化域の間の深さ約40kmの太平洋プレート上面付近で長期的スロースリップが発生し、それに伴う応力変化によって地震発生レートが変化した可能性がある。

（8）平成21年度の成果に関連の深いもので、平成21年度に公表された主な成果物（論文・報告書等）

• Fukuda, J., K. M. Johnson, K. M. Larson, and S. Miyazaki (2009), Fault friction parameters inferred from the early stages of afterslip following the 2003 Tokachi-oki earthquake, J. Geophys. Res., 114, B04412, doi:10.1029/2008JB006166. Fukuyama, E., R. Ando, C. Hashimoto, S. Aoi, and M. Matsu’ura　(2009), A physics-based simulation of the 2003 Tokachi-oki, Japan, earthquake to predict strong ground motions, Bull. Seismol. Soc. Am., 99, 3150-3171. Hashimoto, C., A. Noda, T. Sagiya, and M. Matsu’ura (2009), Interplate seismogenic zones along the Kuril-Japan trench inferred from GPS data inversion, Nature Geoscience, 141-144. Kato, N. (2009), A possible explanation for difference in stress drop between intraplate and interplate earthquakes, Geophys. Res. Lett., 36, L23311, doi:10.1029/2009GL040985. Mitsui, N., T. Hori, S. Miyazaki, amd K. Nakamura (2010), Constraining interplate frictional parameters by using limited terms of synthetic observation data for afterslip: a preliminary test of data assimilation, Theoretical and Applied Mechanics Japan, 58, in press.

（9）平成22年度実施計画の概要

(a) 地殻活動予測シミュレーションモデルの開発研究および予測実験
多重極展開を用いた計算高速化を進めるとともに、不均質粘弾性構造、動的破壊過程、プレート境界地震と内陸地震の相互作用等を考慮して、モデルの高度化を進める。相互作用評価を継続する。
(b) 地殻変動データを用いたデータ同化手法の開発と同化実験
地殻活動予測システムの全体設計を完了し、アジョイント法やグリッドサーチを用いたパラメータ初期値の最適化手法、逐次データ同化手法の開発、試験を行う。
(c)地震活動データの活用手法の開発
小繰り返し地震とGPSデータを用いてプレート間すべりの推定を行う手法を開発する。 また、地震活動データから応力変化の推定を試み、その有効性を検証する。

（10）実施機関の参加者氏名または部署等名

名古屋大学大学院環境学研究科：鷺谷威、橋本千尋，伊藤武男

（11）問い合わせ先

• 部署名等
  名古屋大学大学院環境学研究科
• 電話
  052-789-3046
• e-mail
  web-master@seis.nagoya-u.ac.jp
• URL