層構造の推定方法、及び層構造の解析を行う解析装置
【課題】逆解析による解の精度を高め、層構造を高い精度で推定する。
【解決手段】本発明の層構造の推定方法では、まず、地震観測点における水平成分及び上下成分の観測データに基づいて、水平／上下スペクトル比を算出する。次に、水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定する。そして、地震観測点における層構造を未知数として、前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定する。
【解決手段】本発明の層構造の推定方法では、まず、地震観測点における水平成分及び上下成分の観測データに基づいて、水平／上下スペクトル比を算出する。次に、水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定する。そして、地震観測点における層構造を未知数として、前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定する。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、層構造の推定方法、及び層構造の解析を行う解析装置に関する。
【背景技術】
【０００２】
従来、任意の測定地点においてボーリングが行われ、このボーリング資料に基づいて測定地点における層構造を測定することが行われている。しかし、ボーリングによる層構造の測定には、ボーリング作業によるコストがかかってしまう。
【０００３】
そこで、地表面における地震観測記録を活用することが考えられている。例えば、地震観測記録の水平／上下スペクトル振幅比に基づいて各層の速度構造の逆解析を行うことが考えられている（非特許文献１）。
【非特許文献１】小林喜久二、植竹富一、真下貢、小林啓美、「地震動のＰ波部分における水平／上下スペクトル比の逆解析による深部地下速度構造の推定」(ESTIMATION OF DEEP UNDERGROUND VELOCITY STRUCTURES BY INVERSION OF SPECTRAL RATIO OF HORIZONTAL TO VERTICAL COMPONENT IN P-WAVE PART OF EARTHQUAKE GROUND MOTION)、Proceedings of 12th World Conference on Earthquake Engineering（ＣＤ−ＲＯＭ）論文番号2658、New Zealand Society for Earthquake Engineering、ニュージーランド、(2000)
【発明の開示】
【発明が解決しようとする課題】
【０００４】
非特許文献１では、水平／上下スペクトル比の絶対値である水平／上下スペクトル振幅比に基づいて、逆解析を行っている。しかし、単に水平／上下スペクトル振幅比を用いて逆解析を行っただけでは、解の一義性を確保することができないおそれがあり、解の精度も低いものとなる。
そこで、本発明では、水平／上下スペクトル比の位相も利用して、逆解析による解の精度を高め、層構造を高い精度で推定することを目的とする。
【課題を解決するための手段】
【０００５】
上記目的を達成するための主たる発明は、地震観測点における水平成分及び上下成分の観測データに基づいて水平／上下スペクトル比を算出し、水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、地震観測点における層構造を未知数として、前記目標データに基づいて逆解析により前記未知数を求めて地震観測点における層構造を推定することを特徴とする。
【０００６】
本発明の他の特徴については、本明細書及び添付図面の記載により明らかにする。
【発明を実施するための最良の形態】
【０００７】
＝＝＝開示の概要＝＝＝
本明細書及び添付図面の記載により、少なくとも、以下の事項が明らかとなる。
【０００８】
地震観測点における水平成分及び上下成分の観測データに基づいて、水平／上下スペクトル比を算出し、
水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、
地震観測点における層構造を未知数として、前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定する
ことを特徴とする層構造の推定方法。
このような層構造の推定方法によれば、逆解析の解の精度を高めることができ、層構造を高い精度で推定することができる。
【０００９】
かかる層構造の推定方法であって、前記目標データは、前記水平／上下スペクトル比の絶対値と、前記水平／上下スペクトル比の位相に基づいて算出されるレシーバーファンクションとを含むことが望ましい。これにより、解の一義性を確保でき、逆解析の解の精度を高めることができる。
【００１０】
かかる層構造の推定方法であって、前記逆解析では、前記目標データを用いて個体を評価する遺伝的アルゴリズムが用いられることが望ましい。これにより、少ない計算量で、目標データに比較的近い解を安定的に算出することができる。また、前記逆解析では、遺伝的アルゴリズムによる解の探索が複数回行われ、得られた複数の解に基づいて探索範囲が設定され、設定された探索範囲に基づいて、再度遺伝的アルゴリズムによる解の探索が行われ、前記未知数の値が求められることが望ましい。これにより、より最適な解を安定的に算出することができる。
【００１１】
層構造の解析を行う解析装置であって、
地震観測点における水平成分及び上下成分の観測データを記憶するメモリと、
前記観測データに基づいて水平／上下スペクトル比を算出し、水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、地震観測点における層構造を未知数として前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定するコントローラと、
を有することを特徴とする解析装置。
このような解析装置によれば、精度の高い層構造解析装置を提供できる。
【００１２】
＝＝＝本実施形態の層構造推定方法の概要＝＝＝
本実施形態では、地震観測点における地震観測記録を活用して、地盤の層構造を推定する。地震観測点の地震観測記録は例えばＫ−ＮＥＴ等において公開されており、本実施形態における層構造の推定の際には、このような公開された地震観測記録を利用することができる。
【００１３】
地震観測点は通常地表に設けられているので、地震観測記録に基づいて層構造の推定ができれば、ボーリング等を必要とせずに深部の層構造を探査できることになる。そして、地震観測点が層構造の探査点にもなるので、全国各地の地震観測点の地震観測記録があれば、各地震観測点において層構造を探査することができることになる。本実施形態では、ある地震観測点における層構造の探査について説明を行なう。
【００１４】
図１は、ある地震観測点における観測データである。観測データには、地震観測点における水平成分の時間変化データと、上下成分の時間変化データが含まれている。図中の点線で囲まれた部分は、初期微動部である。この初期微動部は、地震基盤を伝播するＰ波に起因して形成されるものであり、地震基盤を伝播するＳ波の影響は受けていないと考えられている。但し、観測データの水平成分の波形を見て分かる通り、水平成分にも初期微動部が現れている。
【００１５】
図２は、初期微動部における水平成分の生成理由の説明図である。以下、初期微動部において水平成分が生成される理由について説明する。震源からのＰ波が地震基盤中を伝播し、このＰ波が堆積層と地震基盤との境界面に到達すると、堆積層にＰ波が透過するだけでなく、Ｐ波のエネルギーの一部がＳ波に変換されて、堆積層にＳ波が伝播する。つまり、震源からのＳ波が境界面に到達する前に、堆積層においてＳ波が伝播する状態になる。この堆積層のＳ波が、初期微動部の水平成分として現れるのである。（なお、震源からのＰ波のエネルギーの一部は、境界面で反射もしている。）
地震観測点では堆積層を伝播した波を観測しているので、地震観測点で観測された観測データには、地震観測点における堆積層の波動伝播速度（Ｐ波及びＳ波の波動伝播速度）が反映されている。そこで、観測データに基づいて堆積層の波動伝播速度を算出し、地震観測点における堆積層の構造（速度構造）を推定することが考えられる。
【００１６】
但し、観測データの初期微動部には、地震基盤を伝播するＳ波の影響はないものの、地震基盤を伝播したＰ波の影響が当然含まれている。このため、観測データの初期微動部に基づいて層構造を推定するためには、地震基盤を伝播したＰ波の特性を除去する必要がある。しかし、地震基盤を伝播するＰ波を直接的に測定することは困難である。そこで、本実施形態では、観測データから振動源の特性の影響を除去するため、「水平／上下スペクトル振幅比」を用いている。
【００１７】
図３は、水平／上下スペクトル振幅比の特性の説明図である。同図において、Ｓｐ（ω）は、境界面に入射したＰ波のスペクトルである。Ｈ（ω）は、地震観測点における水平成分のスペクトルである。Ｖ（ω）は、地震観測点における上下成分のスペクトルである。Ｓｐ（ω）は未知であるが、水平成分スペクトルＨ（ω）と上下成分スペクトルＶ（ω）は、初期微動部に基づいて算出できることが知られている。
【００１８】
ここで、水平成分スペクトルＨ（ω）や上下成分スペクトルＶ（ω）は、Ｓｐ（ω）を用いて以下のように表現できる。
Ｈ（ω）＝Ｔｈ（ω）×Ｓｐ（ω）
Ｖ（ω）＝Ｔｖ（ω）×Ｓｐ（ω）
なお、上式において、Ｔｈ（ω）は水平成分の増幅率（水平成分の伝達関数）であり、Ｔｖ（ω）は上下成分の増幅率（上下成分の伝達関数）である。そして、これらＴｈ（ω）及びＴｖ（ω）は、堆積層の構造に依存している。
【００１９】
ここで、水平成分スペクトルＨ（ω）と上下成分スペクトルＶ（ω）の比である水平／上下スペクトル比は、以下のようになる。
Ｈ（ω）／Ｖ（ω）＝Ｔｈ（ω）／Ｔｖ（ω）
つまり、水平／上下スペクトル比は、Ｓｐ（ω）の影響を受けずに、堆積層の構造に依存するものとなる。そこで、本実施形態では、堆積層の波動伝播速度等を未知変数として、地震観測点の観測データから算出された水平／上下スペクトル比Ｈ（ω）／Ｖ（ω）の逆解析により、堆積層の構造の推定を行っている。
【００２０】
加えて、本実施形態では、逆解析の精度を高めるため、水平／上下スペクトル比の位相情報も利用している。具体的には、水平／上下スペクトル比の絶対値（振幅比）と位相とに基づいて、レシーバーファンクションＲＦ（後述）を算出でき、このレシーバーファンクションＲＦの逆解析も行っている。
【００２１】
図４は、本実施形態の層構造推定方法のフロー図である。まず、本実施形態では、逆解析の目標データを作成する（Ｓ００１）。なお、この目標データは、水平／上下スペクトル振幅比とレシーバーファンクションＲＦの２つである。次に、逆解析の対象となる未知数の探索範囲（未知数が取り得る値の範囲）を設定する（Ｓ００２）。そして、この探索範囲において未知数の逆解析を行う（Ｓ００３）。本実施形態では、逆解析に遺伝的アルゴリズムが用いられる。以下、本実施形態の層構造推定方法の各工程を詳しく説明する。
【００２２】
＝＝＝目標データの作成（Ｓ００１）＝＝＝
図５は、目標データ作成（Ｓ００１）の処理のフロー図である。解析装置であるコンピュータは、解析プログラムに従って各処理を実行する。
【００２３】
まず、コンピュータは、所定の地震観測所（ここでは成田）の複数の観測データ（観測データ群）を取得する（Ｓ１０１）。例えば、コンピュータは、Ｋ−ＮＥＴ等において公開されている観測データ群をダウンロードし、コンピュータのハードディスクに格納する。この観測データには、前述の水平成分及び上下成分の時間変化データだけでなく、震源の深さや震央距離（震源から地震観測点までの間の水平方向距離）に関するデータも含まれている。
【００２４】
図６は、成田にある地震観測所（地震観測点）の観測データ群の説明図である。図中の縦軸は震源の深さを示し、横軸は震央距離を示している。また、図中には、観測データ群の各観測データの震源の深さと震央距離に応じた座標に、点が記されている。
【００２５】
次に、コンピュータは、各観測データの見かけ入射角を算出する（Ｓ１０２）。図７は、見かけ入射角の説明図である。図に示す通り、見かけ入射角は、地震観測点から見た震源の角度であり、震源−地震観測点間の方向と地震観測点の鉛直方向との角度である。このため、コンピュータは、観測データに含まれる震源の深さのデータと震央距離のデータとに基づいて、見かけ入射角を算出することができる。図６では、各観測データの見かけ入射角に応じて、点の形状を変えて示している。図中に示す通り、同じ地震観測点の観測データであっても、見かけ入射角は様々であることが分かる。
【００２６】
次に、コンピュータは、見かけ入射角と震央距離が類似する複数の観測データを抽出する（Ｓ１０３）。ここでは、図６に示されるグループＡに属する複数の観測データを抽出する。なお、類似する複数の観測データのみを抽出する理由については、後述する。
【００２７】
次に、コンピュータは、抽出された複数の観測データに基づいて、水平／上下スペクトル振幅比の平均値を算出する（Ｓ１０４）。具体的には、まず、コンピュータは、観測データ毎に、水平成分の時間変化データに基づいて水平成分スペクトルＨ（ω）を算出し、上下成分の時間変化データに基づいて上下成分スペクトルＶ（ω）を算出する。そして、コンピュータは、観測データ毎に、この水平成分スペクトルＨ（ω）と上下成分スペクトルＶ（ω）との比を算出することにより、水平／上下スペクトル比を算出する。そして、コンピュータは、観測データ毎に、複素数を用いて示される水平／上下スペクトル比の絶対値を算出し、これを水平／上下スペクトル振幅比として算出する（図８の左上図参照）。そして、コンピュータは、各観測データの水平／上下スペクトル振幅比の平均値を算出する。水平／上下スペクトル振幅比の平均値は、図８の左上図に示すように、振動数に対する振幅比として表される。
【００２８】
次に、コンピュータは、レシーバーファンクションの平均値を算出する（Ｓ１０５）。具体的には、コンピュータは、前述のＳ１０４において算出された水平／上下スペクトル比をフーリエ逆変換することにより、観測データ毎のレシーバーファンクションＲＦ（図８の左下図参照）を算出する。又は、コンピュータは、前述のＳ１０４において算出された水平／上下スペクトル振幅比（図８の左上図参照）とその位相情報（図８の右上図参照）とに基づいて振幅の時間変化を算出し、観測データ毎のレシーバーファンクションＲＦ（図８の左下図参照）を求めることもできる。したがって、レシーバーファンクションは、水平／上下スペクトル比の位相を反映した情報である。そして、コンピュータは、各観測データのレシーバーファンクションの平均値を算出する。レシーバーファンクションの平均値は、図８の左下図に示すように、振幅の時間変化として表される。なお、このレシーバーファンクションは、境界面に単位パルスのＰ波が見かけ入射角で入射したときの地震観測点における応答と考えることができる。
【００２９】
そして、コンピュータは、上記のＳ１０４及びＳ１０５において算出された「水平／上下スペクトル振幅比の平均値」と「レシーバーファンクションの平均値」を、後述する逆解析の際の目標データとして設定する（Ｓ１０６）。
【００３０】
ところで、目標データを作成する際に、図６に示される全ての観測データを用いずに、一部の観測データのみを抽出する理由について説明する。
【００３１】
図８は、Ａグループの観測データに基づく水平／上下スペクトル比（振幅比と位相）及びレシーバーファンクションの説明図である。図９は、Ｂグループの観測データに基づく水平／上下スペクトル比（振幅比と位相）及びレシーバーファンクションの説明図である。Ａグループの観測データは、見かけ入射角が２０〜３０度であり、震央距離が２５〜３０ｋｍの範囲の観測データである。一方、Ｂグループの観測データは、見かけの入射角が４０〜５０度であり、震央距離が４０〜５０ｋｍの範囲の観測データである。図８と図９とを比較して分かる通り、見かけ入射角や震央距離が異なると、水平／上下スペクトル比（振幅比と位相）やレシーバーファンクションの特性が異なっている。
【００３２】
仮に、目標データの作成の際に、一部の観測データのみを抽出するのではなく、全ての観測データを用いるとすると、特性の大きく異なる複数の水平／上下スペクトル振幅比を平均化したり、特性の大きく異なる複数のレシーバーファンクションを平均化したりすることになる。そうすると、算出された目標データに大きな誤差が含まれるおそれがあり、逆解析の目標データには不向きである。
【００３３】
一方、本実施形態のように、見かけ入射角と震央距離が類似する複数の観測データを抽出するようにすれば、この観測データにより求められる水平／上下スペクトル振幅比やレシーバーファンクションは同様な特性を示すので、これらの平均値には大きな誤差が含まれることはなく、逆解析の目標データに適している。また、観測データにノイズが含まれていても、平均化された水平／上下スペクトル振幅比やレシーバーファンクションにはノイズの影響が低減されている。
【００３４】
＝＝＝逆解析の探索範囲を設定（Ｓ００２）＝＝＝
図１０は、逆解析の探索範囲の設定（Ｓ００２）の処理のフロー図である。解析装置であるコンピュータは、解析プログラムに従って各処理を実行する。
【００３５】
まず、コンピュータは、目標データのレシーバーファンクション（Ａグループの各観測データに基づいて算出されたレシーバーファンクションの平均）に基づいて、ＰＳ−Ｐ時間を算出する（Ｓ２０１）。具体的には、コンピュータは、目標データのレシーバーファンクションの最大ピーク（図８の左下図の三角印）の時間を検出し、この時間をＰＳ−Ｐ時間とする。ここでは、ＰＳ−Ｐ時間は１．１８（ｓ）と算出されたものとする。
【００３６】
なお、ＰＳ−Ｐ時間は、Ｐ波と境界面でＰ波から変換されたＳ波との到達時間差を示す。レシーバーファンクションの特性上、レシーバーファンクションの最大値発生時刻がＰＳ−Ｐ時間として現れる。このＰＳ−Ｐ時間は、堆積層の構造に関係する値であり、経験式（例えばＤｂ＝１．４（ｋｍ／ｓ）×［ＰＳ−Ｐ時間］）に基づいて、地震基盤までの概略深さＤｂを推定することができる。
【００３７】
次に、コンピュータは、堆積層の層モデルの設定を行い（Ｓ２０２）、各層の探索範囲の設定を行う（Ｓ２０３）。図１１は、堆積層の層モデル及び探索範囲の説明図である。
【００３８】
通常、地表面から深い層ほどＳ波速度Ｖｓが速くなるので、０．４ｍ／ｓ〜２．２ｍ／ｓまでの範囲を所定の速度Ｖｓで区切って、層モデルを設定している。このため、層Ｎｏ．４〜Ｎｏ．１１の各層のＳ波速度Ｖｓの探索範囲は、重複することなく連続的な範囲になっている。一方、地表面から比較的浅いＶｓ＝０．４ｋｍ／ｓ以下の層では、既存のボーリング資料などを利用できるので、このボーリング資料が利用できるように層Ｎｏ．１〜Ｎｏ．３が設定され、これらの層のＳ波速度Ｖｓの探索範囲がボーリグ資料に基づいて設定される。このため、層Ｎｏ．１〜Ｎｏ．３のＳ波速度Ｖｓの探索範囲は、一部重複している。
【００３９】
各層のＳ波速度Ｖｓの探索範囲の設定後、各層のＰ波速度Ｖｐの探索範囲が設定される。Ｐ波速度ＶｐはＳ波速度Ｖｓに拘束され、経験式（例えばＶｐ＝１．５５Ｖｓ＋０．８１５）から計算される。そこで、各層のＰ波速度Ｖｐの探索範囲は、Ｓ波速度に基づいて経験式から計算されるＰ波速度Ｖｐに対して、経験式の誤差を考慮した付加範囲を追加して、設定される。この付加範囲は、層の深さに応じて異なっても良い。なお、層Ｎｏ．１〜Ｎｏ．３のＰ波速度Ｖｐの探索範囲の設定の際には、ボーリング資料が利用される。
【００４０】
各層のＳ波速度Ｖｓの探索範囲の設定後、各層の層厚の探索範囲が設定される。一番深い層である層Ｎｏ．１１では、層厚の探索範囲は無限大までに設定される。層Ｎｏ．８〜１０の探索範囲は、レシーバーファンクションから推定される概略深さＤｂに基づいて、０．０１〜Ｄｂ／１．５（ｋｍ）として設定される。同様に、層Ｎｏ．７の探索範囲は０．０１〜Ｄｂ／２．０（ｋｍ）、層Ｎｏ．６の探索範囲は０．０１〜Ｄｂ／３．０（ｋｍ）として、概略深さＤｂに基づいて設定される。これに対し、層Ｎｏ．１〜５では、概略深さＤｂに基づかずに、経験的な範囲又はボーリング資料等を参考にして、探索範囲が設定される。
【００４１】
また、各層のＳ波速度Ｖｓの探索範囲の設定後、各層の減衰特性を示すＱ値の探索範囲が設定される。Ｑ値には周波数ｆ（Ｈｚ）に対する依存性があり、以下のようにモデル化される。
１／Ｑｓ＝１／（Ｑ０・ｆ）＋１／Ｑｉ
そこで、本実施形態では、Ｑ０及びＱｉの探索範囲を各層に設定している。
【００４２】
＝＝＝逆解析（Ｓ００３）＝＝＝
＜逆解析の手法について＞
本実施形態では、逆解析を行う際に、探索範囲の中から、目標データに比較的近い状態になるパラメータ（Ｖｓ、Ｖｐ、Ｈ、Ｑ０、Ｑｉ）の値を探索することになる。このようなパラメータの探索方法として、例えば「全探索法」、「ランダムサーチ」、「最急勾配法（山登り法）」等がある。
【００４３】
全探索法は、各パラメータの探索範囲の全ての可能性について網羅的に調べるものであり、その中で最も目標データに近づいたパラメータを解とするものである。しかし、本実施形態のように多数のパラメータが存在する場合、探索時間が長くなり、現実的ではない。
ランダムサーチでは、探索範囲の中からランダムに解を選び出し、これを所定回数繰り返し、その中で最も目標データに近づいたパラメータを解とするものである。但し、ランダムサーチでは、ランダムに解が選び出されるため、解の精度が運次第である。
最急勾配法では、ある状態から目標データに近づく方向へ変化するようにパラメータの値を変化させ、どのパラメータをどの方向へ変化させても目標データから遠ざかるとき、そのパラメータの値を解とするものである。山登りにたとえると、現在の地点から隣り合う地点の中で最も標高の高い地点へ移動を繰り返し、標高の高い地点に移動できなくなった地点を頂上（解）とするものである。但し、最急勾配法では、初期値の設定によって、求められる解が変わってくる。山登りにたとえると、登り始めの地点によって、辿り着く頂上（解）が変わってくる。このため、初期値の設定によっては、低い頂上（目標データに遠い解）にしか辿り着けないことがある（局所安定）。
そこで、本実施形態では、「遺伝的アルゴリズム」と呼ばれる手法を利用して、探索範囲の中から、目標データに比較的近い状態になるパラメータ（Ｖｓ、Ｖｐ、Ｈ、Ｑ０、Ｑｉ）を探索する。
【００４４】
＜遺伝的アルゴリズムについて＞
図１２は、遺伝的アルゴリズムの説明図である。遺伝的アルゴリズムは、環境に適合するように生物が進化するように、目標データに徐々に近づくようにパラメータの値を変化させて、最適解を得ようとするものである。以下に説明するように、遺伝的アルゴリズムでは、確率論的なデータ処理と適応度による選択などの操作で最適解を得ようとしている。
【００４５】
まず、コンピュータは、初期集団を生成する（Ｓ４０１）。初期集団には、多数の個体がある（本実施形態では５００個）。各個体は２進数で表された遺伝子を持っており、この遺伝子は各パラメータの値を示している。本実施形態では、ある個体の遺伝子は、各層の各パラメータ（Ｖｓ、Ｖｐ、Ｈ、Ｑ０、Ｑｉ）の特定の値を示している。初期集団には多様性のある個体を用意するため、各パラメータの探索範囲に一様に分布するような遺伝子になるように多数の個体を用意する。
【００４６】
次に、コンピュータは、各個体に対して評価を行う（Ｓ４０２）。環境に適合する個体ほど高い評価が与えられる。本実施形態では、ある個体について評価を行うとき、まず、コンピュータは、その個体の遺伝子の示すパラメータ（各層のＶｓ、Ｖｐ、Ｈ、Ｑ０、Ｑｉ）の値に基づいて、水平／上下スペクトル振幅比と、レシーバーファンクションを算出する。そして、コンピュータは、算出された水平／上下スペクトル振幅比及びレシーバーファンクションと、目標データの水平／上下スペクトル振幅比及びレシーバーファンクションとを比較し、目標データに近くなる個体ほど高い評価を与える。本実施形態では、遺伝子を基にして算出された水平／上下スペクトル振幅比と目標データの水平／上下スペクトル振幅比との差分を１〜１０Ｈｚの範囲で積分し、また、遺伝子を基にして算出されたレシーバーファンクションと目標データのレシーバーファンクションとの差分を０〜３秒の範囲で積分し、２つの積分値の合計に基づいて評価が行われる。なお、水平／上下スペクトル振幅比に対しては、１〜５Ｈｚの範囲には５．０の重み付けがされ、５〜１０Ｈｚの範囲には１．０の重み付けがされる。また、レシーバーファンクションに対しては、０．５の重み付けがされる。
【００４７】
次に、コンピュータは、各個体に対して選択・交叉・突然変異を行う。選択により、評価の低い個体は淘汰され、評価の高い個体がより多くの子孫を残すようになる。交叉では、複数の個体（親）から遺伝子を受け継ぐ新しい個体（子）を生成する。交叉では、個体の評価に応じて親となる確率を定めても良い。突然変異では、低い確率で遺伝子の一部を変化させる（例えば、２進数の遺伝子の「０」を「１」に、又は「１」を「０」に変更する）。突然変異により、交叉だけでは生成できない子を生成できるので、集団内の個体の多様性を維持することができる。但し、最も評価の高い個体に突然変異を起こすと、その時点での最良解が消滅することになるので、評価の高い個体をそのまま複製したり、突然変異の対象から外したりしても良い。なお、交叉や突然変異により生成された個体であっても、その遺伝子の示すパラメータの値が探索範囲（図１１参照）を外れないようにする必要がある。
【００４８】
コンピュータは、Ｓ４０２〜Ｓ４０４の処理を、設定世代数まで繰り返し行う。本実施形態では、第１５０世代まで繰り返す。そして、集団の中の個体の中で最も評価の高い個体が抽出され、その個体の遺伝子の示す各パラメータの値が解となる。
【００４９】
遺伝的アルゴリズムによれば、少ない計算量で、目標データに比較的近い解を安定的に算出することができる。
【００５０】
＜逆解析について＞
図１３は、逆解析（Ｓ００３）の処理のフロー図である。解析装置であるコンピュータは、解析プログラムに従って各処理を実行する。
【００５１】
まず、コンピュータは、遺伝的アルゴリズムを用いて、第１次逆解析を行う（Ｓ３０１）。なお、コンピュータは、初期集団を変えて、遺伝的アルゴリズムによる解の探索を１０回行う。これにより、コンピュータは、１０個の解を取得することになる。但し、１０個の解の中から最良な解を抽出するわけではない。
【００５２】
図１４は、Ｓ３０１の第１次逆解析によって得られた１０個の解の説明図である。同図は、それぞれの解の示す各層のパラメータ（ここではＶｓ、Ｖｐ、Ｈ）をグラフにしたものである。図から分かる通り、初期集団のランダムな設定や、ランダムな突然変異の影響により、遺伝的アルゴリズムにより探索された解は、それぞれ異なる値になっている。そこで、次に、コンピュータは、第１次逆解析の結果に基づいて、各パラメータの値の平均値及び標準偏差を算出する（Ｓ３０２）。
【００５３】
図１５は、第１次逆解析結果から算出された平均値及び標準偏差の説明図である。
次に、コンピュータは、平均値及び標準偏差に基づいて、各層の探索範囲を設定する（Ｓ３０３）。図１６は、第２次逆解析における各層の探索範囲の説明図である。各パラメータの探索範囲の下限値は、平均値から標準偏差を引いた値である。また、各パラメータの探索範囲の上限値は、平均値から標準偏差を足した値である。なお、本実施形態において２段階で逆解析を行うのは、速度構造の絞込みを目的としている。このため、第２次逆解析のＱ値の探索範囲は、第１次逆解析のＱ値の探索範囲（図１１参照）と同じである。
【００５４】
次に、コンピュータは、遺伝的アルゴリズムを用いて、第２次逆解析を行う（Ｓ３０１）。遺伝的アルゴリズムの手法は、第１次逆解析の場合と同様である。なお、第２次逆解析でも遺伝的アルゴリズムによる解の探索を１０回行う。これにより、コンピュータは、１０個の解を取得することになる。但し、第２次逆解析では、この１０個の解を平均化し、その値を逆解析結果とする。
【００５５】
図１７は、目標データと最適データの比較の説明図である。同図では、本実施形態の逆解析結果による水平／上下スペクトル振幅比とレシーバーファンクションが実線で示されている。また、同図では、観測データに基づく水平／上下スペクトル振幅比とレシーバーファンクションが点線で示されている（つまり、点線は、目標データである）。図に示される通り、逆解析結果は観測結果の特性に近似している。
【００５６】
図１８は、第２次解析の結果の説明図である。同図では、本実施形態の解析結果（実線）とともに、検層データ（点線）も示している。解析結果は、地震基盤までの深さや、堆積層における波動伝播速度（Ｖｓ、Ｖｐ）等を精度良く検出できることを示している。
【００５７】
以上説明した本実施形態によれば、コンピュータは、ハードディスクにダウンロードした水平成分及び上下成分の観測データに基づいて、地震基盤を伝播したＰ波の特性が除去された水平／上下スペクトルを算出し、この水平／上下スペクトルに基づいて逆解析のための目標データが設定される。目標データには水平／上下スペクトル比の絶対値と位相を反映したデータが含まれており、具体的には、水平／上下スペクトル比の絶対値と、水平／上下スペクトル比の絶対値及び位相に基づいて算出されるレシーバーファンクションとが目標データとなる。このように、本実施形態では、位相を反映したデータが目標データに含まれることにより、単に水平／上下スペクトル比の絶対値だけを目標データとして逆解析した場合よりも、解の一義性を確保でき、逆解析結果の精度を高めることができる。そして、逆解析結果の精度を高めることにより、堆積層の層構造（波動伝播速度、層厚、減衰特性等）を精度良く推定できる。
【００５８】
また、本実施形態では、目標データを用いて個体を評価する遺伝的アルゴリズムが逆解析に用いられている。つまり、本実施形態では、目標データである環境に適合する個体が生存するように、遺伝的アルゴリズムによる逆解析が行われる。本実施形態では、目標データに水平／上下スペクトル比の位相の影響が含まれているので、単に水平／上下スペクトル比の絶対値だけを目標データとして遺伝的アルゴリズムによる逆解析を行った場合よりも、逆解析結果の精度を高めることができる。また、本実施形態では、パラメータの数が多いため、他の探索方法と比べて、少ない計算量で、目標データに比較的近い解を安定的に算出することができる。
【００５９】
遺伝的アルゴリズムでは、その性質上、最終的に得られる解が、ランダムに生成される初期集団の影響を受ける。このため、前述の最急勾配法と同様に、得られた解が局所安定の解の可能性がある。そこで、本実施形態では、まず第１次逆解析として遺伝的アルゴリズムを複数実施して複数の解を取得し、この解の平均値と標準偏差に基づいて第２次逆解析用の探索範囲が設定され、この探索範囲に基づいて第２次逆解析として遺伝的アルゴリズムが実施される。これにより、より最適な解を安定的に算出することができる。
【００６０】
なお、特に図示していないが、本実施形態は、ハードディスクを備えたコンピュータにより実現されている。メモリ（記憶装置）として機能するハードディスクには、解析プログラムと、ダウンロードした観測データが記憶される。コンピュータ（詳しくは、コンピュータのＣＰＵ）は、解析プログラムに従って上記の各処理を行い、観測データに基づいて堆積層の層構造を解析する。このようなコンピュータによれば、精度の高い層構造解析装置を提供できる。
【図面の簡単な説明】
【００６１】
【図１】ある地震観測点における観測データである。
【図２】初期微動部における水平成分の生成理由の説明図である。
【図３】水平／上下スペクトル振幅比の特性の説明図である。
【図４】本実施形態の層構造推定方法のフロー図である。
【図５】目標データ作成の処理のフロー図である。
【図６】成田にある地震観測所（地震観測点）の観測データ群の説明図である。
【図７】見かけ入射角の説明図である。
【図８】Ａグループの観測データに基づく水平／上下スペクトル比（振幅比と位相）及びレシーバーファンクションの説明図である。
【図９】Ｂグループの観測データに基づく水平／上下スペクトル比（振幅比と位相）及びレシーバーファンクションの説明図である。
【図１０】探索範囲の設定の処理のフロー図である。
【図１１】堆積層の層モデル及び探索範囲の説明図である。
【図１２】遺伝的アルゴリズムの説明図である。
【図１３】逆解析の処理のフロー図である。
【図１４】第１次逆解析によって得られた１０個の解の説明図である。
【図１５】第１次逆解析結果から算出された平均値及び標準偏差の説明図である。
【図１６】各層の探索範囲の説明図である。
【図１７】目標データと最適データの比較の説明図である。
【図１８】第２次解析の結果の説明図である。
【技術分野】
【０００１】
本発明は、層構造の推定方法、及び層構造の解析を行う解析装置に関する。
【背景技術】
【０００２】
従来、任意の測定地点においてボーリングが行われ、このボーリング資料に基づいて測定地点における層構造を測定することが行われている。しかし、ボーリングによる層構造の測定には、ボーリング作業によるコストがかかってしまう。
【０００３】
そこで、地表面における地震観測記録を活用することが考えられている。例えば、地震観測記録の水平／上下スペクトル振幅比に基づいて各層の速度構造の逆解析を行うことが考えられている（非特許文献１）。
【非特許文献１】小林喜久二、植竹富一、真下貢、小林啓美、「地震動のＰ波部分における水平／上下スペクトル比の逆解析による深部地下速度構造の推定」(ESTIMATION OF DEEP UNDERGROUND VELOCITY STRUCTURES BY INVERSION OF SPECTRAL RATIO OF HORIZONTAL TO VERTICAL COMPONENT IN P-WAVE PART OF EARTHQUAKE GROUND MOTION)、Proceedings of 12th World Conference on Earthquake Engineering（ＣＤ−ＲＯＭ）論文番号2658、New Zealand Society for Earthquake Engineering、ニュージーランド、(2000)
【発明の開示】
【発明が解決しようとする課題】
【０００４】
非特許文献１では、水平／上下スペクトル比の絶対値である水平／上下スペクトル振幅比に基づいて、逆解析を行っている。しかし、単に水平／上下スペクトル振幅比を用いて逆解析を行っただけでは、解の一義性を確保することができないおそれがあり、解の精度も低いものとなる。
そこで、本発明では、水平／上下スペクトル比の位相も利用して、逆解析による解の精度を高め、層構造を高い精度で推定することを目的とする。
【課題を解決するための手段】
【０００５】
上記目的を達成するための主たる発明は、地震観測点における水平成分及び上下成分の観測データに基づいて水平／上下スペクトル比を算出し、水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、地震観測点における層構造を未知数として、前記目標データに基づいて逆解析により前記未知数を求めて地震観測点における層構造を推定することを特徴とする。
【０００６】
本発明の他の特徴については、本明細書及び添付図面の記載により明らかにする。
【発明を実施するための最良の形態】
【０００７】
＝＝＝開示の概要＝＝＝
本明細書及び添付図面の記載により、少なくとも、以下の事項が明らかとなる。
【０００８】
地震観測点における水平成分及び上下成分の観測データに基づいて、水平／上下スペクトル比を算出し、
水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、
地震観測点における層構造を未知数として、前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定する
ことを特徴とする層構造の推定方法。
このような層構造の推定方法によれば、逆解析の解の精度を高めることができ、層構造を高い精度で推定することができる。
【０００９】
かかる層構造の推定方法であって、前記目標データは、前記水平／上下スペクトル比の絶対値と、前記水平／上下スペクトル比の位相に基づいて算出されるレシーバーファンクションとを含むことが望ましい。これにより、解の一義性を確保でき、逆解析の解の精度を高めることができる。
【００１０】
かかる層構造の推定方法であって、前記逆解析では、前記目標データを用いて個体を評価する遺伝的アルゴリズムが用いられることが望ましい。これにより、少ない計算量で、目標データに比較的近い解を安定的に算出することができる。また、前記逆解析では、遺伝的アルゴリズムによる解の探索が複数回行われ、得られた複数の解に基づいて探索範囲が設定され、設定された探索範囲に基づいて、再度遺伝的アルゴリズムによる解の探索が行われ、前記未知数の値が求められることが望ましい。これにより、より最適な解を安定的に算出することができる。
【００１１】
層構造の解析を行う解析装置であって、
地震観測点における水平成分及び上下成分の観測データを記憶するメモリと、
前記観測データに基づいて水平／上下スペクトル比を算出し、水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、地震観測点における層構造を未知数として前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定するコントローラと、
を有することを特徴とする解析装置。
このような解析装置によれば、精度の高い層構造解析装置を提供できる。
【００１２】
＝＝＝本実施形態の層構造推定方法の概要＝＝＝
本実施形態では、地震観測点における地震観測記録を活用して、地盤の層構造を推定する。地震観測点の地震観測記録は例えばＫ−ＮＥＴ等において公開されており、本実施形態における層構造の推定の際には、このような公開された地震観測記録を利用することができる。
【００１３】
地震観測点は通常地表に設けられているので、地震観測記録に基づいて層構造の推定ができれば、ボーリング等を必要とせずに深部の層構造を探査できることになる。そして、地震観測点が層構造の探査点にもなるので、全国各地の地震観測点の地震観測記録があれば、各地震観測点において層構造を探査することができることになる。本実施形態では、ある地震観測点における層構造の探査について説明を行なう。
【００１４】
図１は、ある地震観測点における観測データである。観測データには、地震観測点における水平成分の時間変化データと、上下成分の時間変化データが含まれている。図中の点線で囲まれた部分は、初期微動部である。この初期微動部は、地震基盤を伝播するＰ波に起因して形成されるものであり、地震基盤を伝播するＳ波の影響は受けていないと考えられている。但し、観測データの水平成分の波形を見て分かる通り、水平成分にも初期微動部が現れている。
【００１５】
図２は、初期微動部における水平成分の生成理由の説明図である。以下、初期微動部において水平成分が生成される理由について説明する。震源からのＰ波が地震基盤中を伝播し、このＰ波が堆積層と地震基盤との境界面に到達すると、堆積層にＰ波が透過するだけでなく、Ｐ波のエネルギーの一部がＳ波に変換されて、堆積層にＳ波が伝播する。つまり、震源からのＳ波が境界面に到達する前に、堆積層においてＳ波が伝播する状態になる。この堆積層のＳ波が、初期微動部の水平成分として現れるのである。（なお、震源からのＰ波のエネルギーの一部は、境界面で反射もしている。）
地震観測点では堆積層を伝播した波を観測しているので、地震観測点で観測された観測データには、地震観測点における堆積層の波動伝播速度（Ｐ波及びＳ波の波動伝播速度）が反映されている。そこで、観測データに基づいて堆積層の波動伝播速度を算出し、地震観測点における堆積層の構造（速度構造）を推定することが考えられる。
【００１６】
但し、観測データの初期微動部には、地震基盤を伝播するＳ波の影響はないものの、地震基盤を伝播したＰ波の影響が当然含まれている。このため、観測データの初期微動部に基づいて層構造を推定するためには、地震基盤を伝播したＰ波の特性を除去する必要がある。しかし、地震基盤を伝播するＰ波を直接的に測定することは困難である。そこで、本実施形態では、観測データから振動源の特性の影響を除去するため、「水平／上下スペクトル振幅比」を用いている。
【００１７】
図３は、水平／上下スペクトル振幅比の特性の説明図である。同図において、Ｓｐ（ω）は、境界面に入射したＰ波のスペクトルである。Ｈ（ω）は、地震観測点における水平成分のスペクトルである。Ｖ（ω）は、地震観測点における上下成分のスペクトルである。Ｓｐ（ω）は未知であるが、水平成分スペクトルＨ（ω）と上下成分スペクトルＶ（ω）は、初期微動部に基づいて算出できることが知られている。
【００１８】
ここで、水平成分スペクトルＨ（ω）や上下成分スペクトルＶ（ω）は、Ｓｐ（ω）を用いて以下のように表現できる。
Ｈ（ω）＝Ｔｈ（ω）×Ｓｐ（ω）
Ｖ（ω）＝Ｔｖ（ω）×Ｓｐ（ω）
なお、上式において、Ｔｈ（ω）は水平成分の増幅率（水平成分の伝達関数）であり、Ｔｖ（ω）は上下成分の増幅率（上下成分の伝達関数）である。そして、これらＴｈ（ω）及びＴｖ（ω）は、堆積層の構造に依存している。
【００１９】
ここで、水平成分スペクトルＨ（ω）と上下成分スペクトルＶ（ω）の比である水平／上下スペクトル比は、以下のようになる。
Ｈ（ω）／Ｖ（ω）＝Ｔｈ（ω）／Ｔｖ（ω）
つまり、水平／上下スペクトル比は、Ｓｐ（ω）の影響を受けずに、堆積層の構造に依存するものとなる。そこで、本実施形態では、堆積層の波動伝播速度等を未知変数として、地震観測点の観測データから算出された水平／上下スペクトル比Ｈ（ω）／Ｖ（ω）の逆解析により、堆積層の構造の推定を行っている。
【００２０】
加えて、本実施形態では、逆解析の精度を高めるため、水平／上下スペクトル比の位相情報も利用している。具体的には、水平／上下スペクトル比の絶対値（振幅比）と位相とに基づいて、レシーバーファンクションＲＦ（後述）を算出でき、このレシーバーファンクションＲＦの逆解析も行っている。
【００２１】
図４は、本実施形態の層構造推定方法のフロー図である。まず、本実施形態では、逆解析の目標データを作成する（Ｓ００１）。なお、この目標データは、水平／上下スペクトル振幅比とレシーバーファンクションＲＦの２つである。次に、逆解析の対象となる未知数の探索範囲（未知数が取り得る値の範囲）を設定する（Ｓ００２）。そして、この探索範囲において未知数の逆解析を行う（Ｓ００３）。本実施形態では、逆解析に遺伝的アルゴリズムが用いられる。以下、本実施形態の層構造推定方法の各工程を詳しく説明する。
【００２２】
＝＝＝目標データの作成（Ｓ００１）＝＝＝
図５は、目標データ作成（Ｓ００１）の処理のフロー図である。解析装置であるコンピュータは、解析プログラムに従って各処理を実行する。
【００２３】
まず、コンピュータは、所定の地震観測所（ここでは成田）の複数の観測データ（観測データ群）を取得する（Ｓ１０１）。例えば、コンピュータは、Ｋ−ＮＥＴ等において公開されている観測データ群をダウンロードし、コンピュータのハードディスクに格納する。この観測データには、前述の水平成分及び上下成分の時間変化データだけでなく、震源の深さや震央距離（震源から地震観測点までの間の水平方向距離）に関するデータも含まれている。
【００２４】
図６は、成田にある地震観測所（地震観測点）の観測データ群の説明図である。図中の縦軸は震源の深さを示し、横軸は震央距離を示している。また、図中には、観測データ群の各観測データの震源の深さと震央距離に応じた座標に、点が記されている。
【００２５】
次に、コンピュータは、各観測データの見かけ入射角を算出する（Ｓ１０２）。図７は、見かけ入射角の説明図である。図に示す通り、見かけ入射角は、地震観測点から見た震源の角度であり、震源−地震観測点間の方向と地震観測点の鉛直方向との角度である。このため、コンピュータは、観測データに含まれる震源の深さのデータと震央距離のデータとに基づいて、見かけ入射角を算出することができる。図６では、各観測データの見かけ入射角に応じて、点の形状を変えて示している。図中に示す通り、同じ地震観測点の観測データであっても、見かけ入射角は様々であることが分かる。
【００２６】
次に、コンピュータは、見かけ入射角と震央距離が類似する複数の観測データを抽出する（Ｓ１０３）。ここでは、図６に示されるグループＡに属する複数の観測データを抽出する。なお、類似する複数の観測データのみを抽出する理由については、後述する。
【００２７】
次に、コンピュータは、抽出された複数の観測データに基づいて、水平／上下スペクトル振幅比の平均値を算出する（Ｓ１０４）。具体的には、まず、コンピュータは、観測データ毎に、水平成分の時間変化データに基づいて水平成分スペクトルＨ（ω）を算出し、上下成分の時間変化データに基づいて上下成分スペクトルＶ（ω）を算出する。そして、コンピュータは、観測データ毎に、この水平成分スペクトルＨ（ω）と上下成分スペクトルＶ（ω）との比を算出することにより、水平／上下スペクトル比を算出する。そして、コンピュータは、観測データ毎に、複素数を用いて示される水平／上下スペクトル比の絶対値を算出し、これを水平／上下スペクトル振幅比として算出する（図８の左上図参照）。そして、コンピュータは、各観測データの水平／上下スペクトル振幅比の平均値を算出する。水平／上下スペクトル振幅比の平均値は、図８の左上図に示すように、振動数に対する振幅比として表される。
【００２８】
次に、コンピュータは、レシーバーファンクションの平均値を算出する（Ｓ１０５）。具体的には、コンピュータは、前述のＳ１０４において算出された水平／上下スペクトル比をフーリエ逆変換することにより、観測データ毎のレシーバーファンクションＲＦ（図８の左下図参照）を算出する。又は、コンピュータは、前述のＳ１０４において算出された水平／上下スペクトル振幅比（図８の左上図参照）とその位相情報（図８の右上図参照）とに基づいて振幅の時間変化を算出し、観測データ毎のレシーバーファンクションＲＦ（図８の左下図参照）を求めることもできる。したがって、レシーバーファンクションは、水平／上下スペクトル比の位相を反映した情報である。そして、コンピュータは、各観測データのレシーバーファンクションの平均値を算出する。レシーバーファンクションの平均値は、図８の左下図に示すように、振幅の時間変化として表される。なお、このレシーバーファンクションは、境界面に単位パルスのＰ波が見かけ入射角で入射したときの地震観測点における応答と考えることができる。
【００２９】
そして、コンピュータは、上記のＳ１０４及びＳ１０５において算出された「水平／上下スペクトル振幅比の平均値」と「レシーバーファンクションの平均値」を、後述する逆解析の際の目標データとして設定する（Ｓ１０６）。
【００３０】
ところで、目標データを作成する際に、図６に示される全ての観測データを用いずに、一部の観測データのみを抽出する理由について説明する。
【００３１】
図８は、Ａグループの観測データに基づく水平／上下スペクトル比（振幅比と位相）及びレシーバーファンクションの説明図である。図９は、Ｂグループの観測データに基づく水平／上下スペクトル比（振幅比と位相）及びレシーバーファンクションの説明図である。Ａグループの観測データは、見かけ入射角が２０〜３０度であり、震央距離が２５〜３０ｋｍの範囲の観測データである。一方、Ｂグループの観測データは、見かけの入射角が４０〜５０度であり、震央距離が４０〜５０ｋｍの範囲の観測データである。図８と図９とを比較して分かる通り、見かけ入射角や震央距離が異なると、水平／上下スペクトル比（振幅比と位相）やレシーバーファンクションの特性が異なっている。
【００３２】
仮に、目標データの作成の際に、一部の観測データのみを抽出するのではなく、全ての観測データを用いるとすると、特性の大きく異なる複数の水平／上下スペクトル振幅比を平均化したり、特性の大きく異なる複数のレシーバーファンクションを平均化したりすることになる。そうすると、算出された目標データに大きな誤差が含まれるおそれがあり、逆解析の目標データには不向きである。
【００３３】
一方、本実施形態のように、見かけ入射角と震央距離が類似する複数の観測データを抽出するようにすれば、この観測データにより求められる水平／上下スペクトル振幅比やレシーバーファンクションは同様な特性を示すので、これらの平均値には大きな誤差が含まれることはなく、逆解析の目標データに適している。また、観測データにノイズが含まれていても、平均化された水平／上下スペクトル振幅比やレシーバーファンクションにはノイズの影響が低減されている。
【００３４】
＝＝＝逆解析の探索範囲を設定（Ｓ００２）＝＝＝
図１０は、逆解析の探索範囲の設定（Ｓ００２）の処理のフロー図である。解析装置であるコンピュータは、解析プログラムに従って各処理を実行する。
【００３５】
まず、コンピュータは、目標データのレシーバーファンクション（Ａグループの各観測データに基づいて算出されたレシーバーファンクションの平均）に基づいて、ＰＳ−Ｐ時間を算出する（Ｓ２０１）。具体的には、コンピュータは、目標データのレシーバーファンクションの最大ピーク（図８の左下図の三角印）の時間を検出し、この時間をＰＳ−Ｐ時間とする。ここでは、ＰＳ−Ｐ時間は１．１８（ｓ）と算出されたものとする。
【００３６】
なお、ＰＳ−Ｐ時間は、Ｐ波と境界面でＰ波から変換されたＳ波との到達時間差を示す。レシーバーファンクションの特性上、レシーバーファンクションの最大値発生時刻がＰＳ−Ｐ時間として現れる。このＰＳ−Ｐ時間は、堆積層の構造に関係する値であり、経験式（例えばＤｂ＝１．４（ｋｍ／ｓ）×［ＰＳ−Ｐ時間］）に基づいて、地震基盤までの概略深さＤｂを推定することができる。
【００３７】
次に、コンピュータは、堆積層の層モデルの設定を行い（Ｓ２０２）、各層の探索範囲の設定を行う（Ｓ２０３）。図１１は、堆積層の層モデル及び探索範囲の説明図である。
【００３８】
通常、地表面から深い層ほどＳ波速度Ｖｓが速くなるので、０．４ｍ／ｓ〜２．２ｍ／ｓまでの範囲を所定の速度Ｖｓで区切って、層モデルを設定している。このため、層Ｎｏ．４〜Ｎｏ．１１の各層のＳ波速度Ｖｓの探索範囲は、重複することなく連続的な範囲になっている。一方、地表面から比較的浅いＶｓ＝０．４ｋｍ／ｓ以下の層では、既存のボーリング資料などを利用できるので、このボーリング資料が利用できるように層Ｎｏ．１〜Ｎｏ．３が設定され、これらの層のＳ波速度Ｖｓの探索範囲がボーリグ資料に基づいて設定される。このため、層Ｎｏ．１〜Ｎｏ．３のＳ波速度Ｖｓの探索範囲は、一部重複している。
【００３９】
各層のＳ波速度Ｖｓの探索範囲の設定後、各層のＰ波速度Ｖｐの探索範囲が設定される。Ｐ波速度ＶｐはＳ波速度Ｖｓに拘束され、経験式（例えばＶｐ＝１．５５Ｖｓ＋０．８１５）から計算される。そこで、各層のＰ波速度Ｖｐの探索範囲は、Ｓ波速度に基づいて経験式から計算されるＰ波速度Ｖｐに対して、経験式の誤差を考慮した付加範囲を追加して、設定される。この付加範囲は、層の深さに応じて異なっても良い。なお、層Ｎｏ．１〜Ｎｏ．３のＰ波速度Ｖｐの探索範囲の設定の際には、ボーリング資料が利用される。
【００４０】
各層のＳ波速度Ｖｓの探索範囲の設定後、各層の層厚の探索範囲が設定される。一番深い層である層Ｎｏ．１１では、層厚の探索範囲は無限大までに設定される。層Ｎｏ．８〜１０の探索範囲は、レシーバーファンクションから推定される概略深さＤｂに基づいて、０．０１〜Ｄｂ／１．５（ｋｍ）として設定される。同様に、層Ｎｏ．７の探索範囲は０．０１〜Ｄｂ／２．０（ｋｍ）、層Ｎｏ．６の探索範囲は０．０１〜Ｄｂ／３．０（ｋｍ）として、概略深さＤｂに基づいて設定される。これに対し、層Ｎｏ．１〜５では、概略深さＤｂに基づかずに、経験的な範囲又はボーリング資料等を参考にして、探索範囲が設定される。
【００４１】
また、各層のＳ波速度Ｖｓの探索範囲の設定後、各層の減衰特性を示すＱ値の探索範囲が設定される。Ｑ値には周波数ｆ（Ｈｚ）に対する依存性があり、以下のようにモデル化される。
１／Ｑｓ＝１／（Ｑ０・ｆ）＋１／Ｑｉ
そこで、本実施形態では、Ｑ０及びＱｉの探索範囲を各層に設定している。
【００４２】
＝＝＝逆解析（Ｓ００３）＝＝＝
＜逆解析の手法について＞
本実施形態では、逆解析を行う際に、探索範囲の中から、目標データに比較的近い状態になるパラメータ（Ｖｓ、Ｖｐ、Ｈ、Ｑ０、Ｑｉ）の値を探索することになる。このようなパラメータの探索方法として、例えば「全探索法」、「ランダムサーチ」、「最急勾配法（山登り法）」等がある。
【００４３】
全探索法は、各パラメータの探索範囲の全ての可能性について網羅的に調べるものであり、その中で最も目標データに近づいたパラメータを解とするものである。しかし、本実施形態のように多数のパラメータが存在する場合、探索時間が長くなり、現実的ではない。
ランダムサーチでは、探索範囲の中からランダムに解を選び出し、これを所定回数繰り返し、その中で最も目標データに近づいたパラメータを解とするものである。但し、ランダムサーチでは、ランダムに解が選び出されるため、解の精度が運次第である。
最急勾配法では、ある状態から目標データに近づく方向へ変化するようにパラメータの値を変化させ、どのパラメータをどの方向へ変化させても目標データから遠ざかるとき、そのパラメータの値を解とするものである。山登りにたとえると、現在の地点から隣り合う地点の中で最も標高の高い地点へ移動を繰り返し、標高の高い地点に移動できなくなった地点を頂上（解）とするものである。但し、最急勾配法では、初期値の設定によって、求められる解が変わってくる。山登りにたとえると、登り始めの地点によって、辿り着く頂上（解）が変わってくる。このため、初期値の設定によっては、低い頂上（目標データに遠い解）にしか辿り着けないことがある（局所安定）。
そこで、本実施形態では、「遺伝的アルゴリズム」と呼ばれる手法を利用して、探索範囲の中から、目標データに比較的近い状態になるパラメータ（Ｖｓ、Ｖｐ、Ｈ、Ｑ０、Ｑｉ）を探索する。
【００４４】
＜遺伝的アルゴリズムについて＞
図１２は、遺伝的アルゴリズムの説明図である。遺伝的アルゴリズムは、環境に適合するように生物が進化するように、目標データに徐々に近づくようにパラメータの値を変化させて、最適解を得ようとするものである。以下に説明するように、遺伝的アルゴリズムでは、確率論的なデータ処理と適応度による選択などの操作で最適解を得ようとしている。
【００４５】
まず、コンピュータは、初期集団を生成する（Ｓ４０１）。初期集団には、多数の個体がある（本実施形態では５００個）。各個体は２進数で表された遺伝子を持っており、この遺伝子は各パラメータの値を示している。本実施形態では、ある個体の遺伝子は、各層の各パラメータ（Ｖｓ、Ｖｐ、Ｈ、Ｑ０、Ｑｉ）の特定の値を示している。初期集団には多様性のある個体を用意するため、各パラメータの探索範囲に一様に分布するような遺伝子になるように多数の個体を用意する。
【００４６】
次に、コンピュータは、各個体に対して評価を行う（Ｓ４０２）。環境に適合する個体ほど高い評価が与えられる。本実施形態では、ある個体について評価を行うとき、まず、コンピュータは、その個体の遺伝子の示すパラメータ（各層のＶｓ、Ｖｐ、Ｈ、Ｑ０、Ｑｉ）の値に基づいて、水平／上下スペクトル振幅比と、レシーバーファンクションを算出する。そして、コンピュータは、算出された水平／上下スペクトル振幅比及びレシーバーファンクションと、目標データの水平／上下スペクトル振幅比及びレシーバーファンクションとを比較し、目標データに近くなる個体ほど高い評価を与える。本実施形態では、遺伝子を基にして算出された水平／上下スペクトル振幅比と目標データの水平／上下スペクトル振幅比との差分を１〜１０Ｈｚの範囲で積分し、また、遺伝子を基にして算出されたレシーバーファンクションと目標データのレシーバーファンクションとの差分を０〜３秒の範囲で積分し、２つの積分値の合計に基づいて評価が行われる。なお、水平／上下スペクトル振幅比に対しては、１〜５Ｈｚの範囲には５．０の重み付けがされ、５〜１０Ｈｚの範囲には１．０の重み付けがされる。また、レシーバーファンクションに対しては、０．５の重み付けがされる。
【００４７】
次に、コンピュータは、各個体に対して選択・交叉・突然変異を行う。選択により、評価の低い個体は淘汰され、評価の高い個体がより多くの子孫を残すようになる。交叉では、複数の個体（親）から遺伝子を受け継ぐ新しい個体（子）を生成する。交叉では、個体の評価に応じて親となる確率を定めても良い。突然変異では、低い確率で遺伝子の一部を変化させる（例えば、２進数の遺伝子の「０」を「１」に、又は「１」を「０」に変更する）。突然変異により、交叉だけでは生成できない子を生成できるので、集団内の個体の多様性を維持することができる。但し、最も評価の高い個体に突然変異を起こすと、その時点での最良解が消滅することになるので、評価の高い個体をそのまま複製したり、突然変異の対象から外したりしても良い。なお、交叉や突然変異により生成された個体であっても、その遺伝子の示すパラメータの値が探索範囲（図１１参照）を外れないようにする必要がある。
【００４８】
コンピュータは、Ｓ４０２〜Ｓ４０４の処理を、設定世代数まで繰り返し行う。本実施形態では、第１５０世代まで繰り返す。そして、集団の中の個体の中で最も評価の高い個体が抽出され、その個体の遺伝子の示す各パラメータの値が解となる。
【００４９】
遺伝的アルゴリズムによれば、少ない計算量で、目標データに比較的近い解を安定的に算出することができる。
【００５０】
＜逆解析について＞
図１３は、逆解析（Ｓ００３）の処理のフロー図である。解析装置であるコンピュータは、解析プログラムに従って各処理を実行する。
【００５１】
まず、コンピュータは、遺伝的アルゴリズムを用いて、第１次逆解析を行う（Ｓ３０１）。なお、コンピュータは、初期集団を変えて、遺伝的アルゴリズムによる解の探索を１０回行う。これにより、コンピュータは、１０個の解を取得することになる。但し、１０個の解の中から最良な解を抽出するわけではない。
【００５２】
図１４は、Ｓ３０１の第１次逆解析によって得られた１０個の解の説明図である。同図は、それぞれの解の示す各層のパラメータ（ここではＶｓ、Ｖｐ、Ｈ）をグラフにしたものである。図から分かる通り、初期集団のランダムな設定や、ランダムな突然変異の影響により、遺伝的アルゴリズムにより探索された解は、それぞれ異なる値になっている。そこで、次に、コンピュータは、第１次逆解析の結果に基づいて、各パラメータの値の平均値及び標準偏差を算出する（Ｓ３０２）。
【００５３】
図１５は、第１次逆解析結果から算出された平均値及び標準偏差の説明図である。
次に、コンピュータは、平均値及び標準偏差に基づいて、各層の探索範囲を設定する（Ｓ３０３）。図１６は、第２次逆解析における各層の探索範囲の説明図である。各パラメータの探索範囲の下限値は、平均値から標準偏差を引いた値である。また、各パラメータの探索範囲の上限値は、平均値から標準偏差を足した値である。なお、本実施形態において２段階で逆解析を行うのは、速度構造の絞込みを目的としている。このため、第２次逆解析のＱ値の探索範囲は、第１次逆解析のＱ値の探索範囲（図１１参照）と同じである。
【００５４】
次に、コンピュータは、遺伝的アルゴリズムを用いて、第２次逆解析を行う（Ｓ３０１）。遺伝的アルゴリズムの手法は、第１次逆解析の場合と同様である。なお、第２次逆解析でも遺伝的アルゴリズムによる解の探索を１０回行う。これにより、コンピュータは、１０個の解を取得することになる。但し、第２次逆解析では、この１０個の解を平均化し、その値を逆解析結果とする。
【００５５】
図１７は、目標データと最適データの比較の説明図である。同図では、本実施形態の逆解析結果による水平／上下スペクトル振幅比とレシーバーファンクションが実線で示されている。また、同図では、観測データに基づく水平／上下スペクトル振幅比とレシーバーファンクションが点線で示されている（つまり、点線は、目標データである）。図に示される通り、逆解析結果は観測結果の特性に近似している。
【００５６】
図１８は、第２次解析の結果の説明図である。同図では、本実施形態の解析結果（実線）とともに、検層データ（点線）も示している。解析結果は、地震基盤までの深さや、堆積層における波動伝播速度（Ｖｓ、Ｖｐ）等を精度良く検出できることを示している。
【００５７】
以上説明した本実施形態によれば、コンピュータは、ハードディスクにダウンロードした水平成分及び上下成分の観測データに基づいて、地震基盤を伝播したＰ波の特性が除去された水平／上下スペクトルを算出し、この水平／上下スペクトルに基づいて逆解析のための目標データが設定される。目標データには水平／上下スペクトル比の絶対値と位相を反映したデータが含まれており、具体的には、水平／上下スペクトル比の絶対値と、水平／上下スペクトル比の絶対値及び位相に基づいて算出されるレシーバーファンクションとが目標データとなる。このように、本実施形態では、位相を反映したデータが目標データに含まれることにより、単に水平／上下スペクトル比の絶対値だけを目標データとして逆解析した場合よりも、解の一義性を確保でき、逆解析結果の精度を高めることができる。そして、逆解析結果の精度を高めることにより、堆積層の層構造（波動伝播速度、層厚、減衰特性等）を精度良く推定できる。
【００５８】
また、本実施形態では、目標データを用いて個体を評価する遺伝的アルゴリズムが逆解析に用いられている。つまり、本実施形態では、目標データである環境に適合する個体が生存するように、遺伝的アルゴリズムによる逆解析が行われる。本実施形態では、目標データに水平／上下スペクトル比の位相の影響が含まれているので、単に水平／上下スペクトル比の絶対値だけを目標データとして遺伝的アルゴリズムによる逆解析を行った場合よりも、逆解析結果の精度を高めることができる。また、本実施形態では、パラメータの数が多いため、他の探索方法と比べて、少ない計算量で、目標データに比較的近い解を安定的に算出することができる。
【００５９】
遺伝的アルゴリズムでは、その性質上、最終的に得られる解が、ランダムに生成される初期集団の影響を受ける。このため、前述の最急勾配法と同様に、得られた解が局所安定の解の可能性がある。そこで、本実施形態では、まず第１次逆解析として遺伝的アルゴリズムを複数実施して複数の解を取得し、この解の平均値と標準偏差に基づいて第２次逆解析用の探索範囲が設定され、この探索範囲に基づいて第２次逆解析として遺伝的アルゴリズムが実施される。これにより、より最適な解を安定的に算出することができる。
【００６０】
なお、特に図示していないが、本実施形態は、ハードディスクを備えたコンピュータにより実現されている。メモリ（記憶装置）として機能するハードディスクには、解析プログラムと、ダウンロードした観測データが記憶される。コンピュータ（詳しくは、コンピュータのＣＰＵ）は、解析プログラムに従って上記の各処理を行い、観測データに基づいて堆積層の層構造を解析する。このようなコンピュータによれば、精度の高い層構造解析装置を提供できる。
【図面の簡単な説明】
【００６１】
【図１】ある地震観測点における観測データである。
【図２】初期微動部における水平成分の生成理由の説明図である。
【図３】水平／上下スペクトル振幅比の特性の説明図である。
【図４】本実施形態の層構造推定方法のフロー図である。
【図５】目標データ作成の処理のフロー図である。
【図６】成田にある地震観測所（地震観測点）の観測データ群の説明図である。
【図７】見かけ入射角の説明図である。
【図８】Ａグループの観測データに基づく水平／上下スペクトル比（振幅比と位相）及びレシーバーファンクションの説明図である。
【図９】Ｂグループの観測データに基づく水平／上下スペクトル比（振幅比と位相）及びレシーバーファンクションの説明図である。
【図１０】探索範囲の設定の処理のフロー図である。
【図１１】堆積層の層モデル及び探索範囲の説明図である。
【図１２】遺伝的アルゴリズムの説明図である。
【図１３】逆解析の処理のフロー図である。
【図１４】第１次逆解析によって得られた１０個の解の説明図である。
【図１５】第１次逆解析結果から算出された平均値及び標準偏差の説明図である。
【図１６】各層の探索範囲の説明図である。
【図１７】目標データと最適データの比較の説明図である。
【図１８】第２次解析の結果の説明図である。
【特許請求の範囲】
【請求項１】
地震観測点における水平成分及び上下成分の観測データに基づいて、水平／上下スペクトル比を算出し、
水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、
地震観測点における層構造を未知数として、前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定する
ことを特徴とする層構造の推定方法。
【請求項２】
請求項１に記載の層構造の推定方法であって、
前記目標データは、前記水平／上下スペクトル比の絶対値と、前記水平／上下スペクトル比の位相に基づいて算出されるレシーバーファンクションとを含む
ことを特徴とする層構造の推定方法。
【請求項３】
請求項１又は２に記載の層構造の推定方法であって、
前記逆解析では、前記目標データを用いて個体を評価する遺伝的アルゴリズムが用いられることを特徴とする層構造の推定方法。
【請求項４】
請求項３に記載の層構造の推定方法であって、
前記逆解析では、
遺伝的アルゴリズムによる解の探索が複数回行われ、
得られた複数の解に基づいて探索範囲が設定され、
設定された探索範囲に基づいて、再度遺伝的アルゴリズムによる解の探索が行われ、前記未知数の値が求められる
ことを特徴とする層構造の推定方法。
【請求項５】
層構造の解析を行う解析装置であって、
地震観測点における水平成分及び上下成分の観測データを記憶するメモリと、
前記観測データに基づいて水平／上下スペクトル比を算出し、水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、地震観測点における層構造を未知数として前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定するコントローラと、
を有することを特徴とする解析装置。
【請求項１】
地震観測点における水平成分及び上下成分の観測データに基づいて、水平／上下スペクトル比を算出し、
水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、
地震観測点における層構造を未知数として、前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定する
ことを特徴とする層構造の推定方法。
【請求項２】
請求項１に記載の層構造の推定方法であって、
前記目標データは、前記水平／上下スペクトル比の絶対値と、前記水平／上下スペクトル比の位相に基づいて算出されるレシーバーファンクションとを含む
ことを特徴とする層構造の推定方法。
【請求項３】
請求項１又は２に記載の層構造の推定方法であって、
前記逆解析では、前記目標データを用いて個体を評価する遺伝的アルゴリズムが用いられることを特徴とする層構造の推定方法。
【請求項４】
請求項３に記載の層構造の推定方法であって、
前記逆解析では、
遺伝的アルゴリズムによる解の探索が複数回行われ、
得られた複数の解に基づいて探索範囲が設定され、
設定された探索範囲に基づいて、再度遺伝的アルゴリズムによる解の探索が行われ、前記未知数の値が求められる
ことを特徴とする層構造の推定方法。
【請求項５】
層構造の解析を行う解析装置であって、
地震観測点における水平成分及び上下成分の観測データを記憶するメモリと、
前記観測データに基づいて水平／上下スペクトル比を算出し、水平／上下スペクトル比の絶対値と位相とに基づく目標データを設定し、地震観測点における層構造を未知数として前記目標データに基づいて逆解析により前記未知数を求めて、地震観測点における層構造を推定するコントローラと、
を有することを特徴とする解析装置。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図１４】
【図１５】
【図１６】
【図１７】
【図１８】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図１４】
【図１５】
【図１６】
【図１７】
【図１８】
【公開番号】特開２００６−３３７１９１（Ｐ２００６−３３７１９１Ａ）
【公開日】平成１８年１２月１４日（２００６．１２．１４）
【国際特許分類】
【出願番号】特願２００５−１６２６６６（Ｐ２００５−１６２６６６）
【出願日】平成１７年６月２日（２００５．６．２）
【出願人】（０００００３６８７）東京電力株式会社 (2,580)
【Ｆターム（参考）】
【公開日】平成１８年１２月１４日（２００６．１２．１４）
【国際特許分類】
【出願日】平成１７年６月２日（２００５．６．２）
【出願人】（０００００３６８７）東京電力株式会社 (2,580)
【Ｆターム（参考）】
[ Back to top ]