(57)【特許請求の範囲】
【請求項１】
地下領域を表す地球物理学的データ又はデータ属性の１つ又は２つ以上の２Ｄ又は３Ｄ離散化セット（かかる各データセットを「原データボリューム」という）中の地質学的特徴部を識別する方法であって、
（ａ）データ窓の形状及びサイズを選択するステップを有し、
（ｂ）各原データボリュームに関して、前記窓を前記原データ中の複数個のオーバーラップ位置又は非オーバーラップ位置に動かして各データボクセルが少なくとも１つの窓内に含まれるようにし、各窓に関して、成分がその窓内からのボクセル値から成るデータ窓ベクトルＩを形成するステップを有し、
（ｃ）前記データ窓ベクトルを用いて統計的分析を実施し、データ値に関する分布をコンピュータ計算するステップを有し、前記統計的分析は、複数個の原データボリュームの場合に共同して実施され、
（ｄ）前記データ値分布を用いて前記データ中の外れ値又は異常を識別するステップを有し、
（ｅ）前記外れ値又は異常を用いて前記地下領域の地質学的特徴を予測するステップを有する、方法。
【請求項２】
前記データ値分布は、
（ｉ）全てのデータ窓ベクトルの平均行列及び共分散行列をコンピュータ計算する方法、
（ｉｉ）独立成分分析法、
（ｉｉｉ）クラスター化法を用いて前記データをクラスター化する方法、及び
（ｉｖ）別の統計的分析法
から成る一群の統計的分析技術のうちの１つを用いてコンピュータ計算される、請求項１記載の方法。
【請求項３】
前記統計的分析は、（ｉ）を用いて実施され、更に、主成分分析を用いるステップを有する、請求項２記載の方法。
【請求項４】
共分散行列の固有値及び固有ベクトルをコンピュータ計算し、前記固有ベクトルは、対応の原データボリュームの１組の主成分であり、
前記ステップ（ｄ）及び前記ステップ（ｅ）は、原データボリュームを前記固有ベクトルの選択されたサブセット上に投射して投射された部分データボリュームを生成するステップを有し、前記固有ベクトルのサブセットは、これらの対応の固有値に基づいて選択され、前記ステップ（ｄ）及び前記ステップ（ｅ）は、前記投射部分データボリューム中に捕捉されない前記原データボリュームの部分である残差データボリュームと、前記残差データボリューム中の異常な特徴を識別し、前記異常な特徴を用いて前記地下領域の物理的特徴を予測するステップとを更に有する、請求項３記載の方法。
【請求項５】
前記データ窓は、Ｎ次元であり、ここで、Ｎは、１≦Ｎ≦Ｍであるような整数であり、Ｍは、データセットの次元数である、請求項１記載の方法。
【請求項６】
前記選択された窓サイズ及び形状に関する前記平均行列及び共分散行列は、相補窓を用いてコンピュータ計算され、ステップ（ａ）で選択された前記窓中の各場所に相当する相補窓は、前記窓を原データボリューム中を動かしているときにその場所に現われる１組のデータ値を表している、請求項３記載の方法。
【請求項７】
前記選択されたサブセットは、テキスチャ、カオス又は他のデータ若しくは幾何学的属性によって識別されたパターンの内部類似性に基づいて選択される、請求項４記載の方法。
【請求項８】
前記固有ベクトルの前記選択されたサブセットは、最も大きい値から最も小さい値まで順番付けられた固有値を合計して最も大きなＮ個の固有値の和を全ての固有値の和で除算して得られた値がＲのあらかじめ選択された値よりも大きくなるようにし、この場合、０＜Ｒ＜１であり、次に、前記Ｎ個の最も大きな固有値と関連した前記Ｎ個の固有ベクトルを選択することによって求められる、請求項４記載の方法。
【請求項９】
地下領域を表す地球物理学的データ又はデータ属性の２Ｄ又は３Ｄ離散化セット（「原データボリューム」）中の地質学的特徴部を識別する方法であって、
（ａ）データ窓の形状及びサイズを選択するステップを有し、
（ｂ）前記窓を前記原データ中の複数個のオーバーラップ位置又は非オーバーラップ位置に動かして各データボクセルが少なくとも１つの窓内に含まれるようにし、各窓に関して、成分がその窓内からのボクセル値から成るデータ窓ベクトルＩを形成するステップを有し、
（ｃ）全てのデータ窓ベクトルの共分散行列をコンピュータ計算するステップを有し
（ｄ）前記共分散行列の固有ベクトルをコンピュータ計算するステップを有し、
（ｅ）前記原データボリュームを前記固有ベクトルの選択されたサブセット上に投射して投射された部分データボリュームを生成するステップを有し、
（ｆ）前記投射部分データボリューム中の外れ値又は異常を識別し、前記外れ値又は異常を用いて前記地下領域の地質学的特徴を予測するステップを有する、方法。
【請求項１０】
投射された部分データボリュームを生成するための前記固有ベクトルの前記選択されたサブセットは、固有ベクトルをこれらの関連の固有値に基づいてなくすことにより求められる、請求項９記載の方法。
【請求項１１】
前記固有ベクトルの前記選択されたサブセットは、ユーザにより対話形式で選択されるか自動的に識別されたノイズ又は幾何学的特性に基づくかのいずれかである、請求項９記載の方法。
【請求項１２】
前記固有ベクトルの前記選択されたサブセットは、前記原データボリューム中の自明な異常を突き止める基準を案出し、前記基準を用いて１つ又は２つ以上の自明な異常を選択し、関連データ成分（固有ベクトル上への前記原データボリュームの投射）が前記選択された自明な異常に寄与し又は所定量を超えるバックグラウンド信号の原因となる１つ又は２つ以上の固有ベクトルを識別し、次に、残りの固有ベクトルのうちの幾つか又は全てを選択することによって求められ、前記ステップ（ｆ）は、前記固有ベクトルの前記選択されたサブセットを求めるために用いられた前記自明な異常よりも微妙な異常の発見を可能にする、請求項９記載の方法。
【請求項１３】
前記ステップ（ｅ）の次に、前記原データボリュームに代えて前記投射部分データボリュームを用いて前記ステップ（ａ）〜（ｅ）を繰り返し実施するステップと、次に前記ステップ（ｆ）で用いられる更新された投射部分データボリュームを生成するステップとを更に有する、請求項１２記載の方法。
【請求項１４】
地下領域を表す地球物理学的データ又はデータ属性の２Ｄ又は３Ｄ離散化セット（「原データボリューム」）中の地質学的特徴部を識別する方法であって、
（ａ）データ窓の形状及びサイズを選択するステップを有し、
（ｂ）前記窓を前記原データ中の複数個のオーバーラップ位置又は非オーバーラップ位置に動かして各データボクセルが少なくとも１つの窓内に含まれるようにし、各窓に関して、成分がその窓内からのボクセル値から成るデータ窓ベクトルＩを形成するステップを有し、
（ｃ）全てのデータ窓ベクトルの共分散行列をコンピュータ計算するステップを有し、  （ｄ）前記共分散行列の固有値及び固有ベクトルをコンピュータ計算するステップを有し、
（ｅ）ボクセルの異常の程度をコンピュータ計算する方法を選択し、この方法を用いて所定の閾値よりも異常であるものとコンピュータ計算されたボクセルから成る部分データボリュームを定めるステップを有し、
（ｆ）前記投射部分データボリューム中の１つ又は２つ以上の外れ値又は異常を識別し、前記外れ値又は異常を用いて前記地下領域の地質学的特徴を予測するステップを有する、方法。
【請求項１５】
添字ｘ，ｙ，ｚ，ｉ，ｊ，ｋにより示されたボクセルの異常の程度Ｒ′は、次式、即ち、
0005530452.tif 000013
（上式において、Ｉi,j,kは、ボクセルｉ，ｊ，ｋを含む前記ステップ（ｂ）からのデータ窓ベクトルの成分であり）、
0005530452.tif 000014
（前記離散化原データボリュームは、Ｎx×Ｎy×Ｎz個のボクセルから成る）からコンピュータ計算され、前記選択された窓の形状及びサイズは、ｎx×ｎy×ｎz個のボクセルであり、Ｎ＝（Ｎx−ｎx）×（Ｎy−ｎy）−（Ｎz−ｎz）である、請求項１４記載の方法。
【請求項１６】
前記異常の程度を求めるのに、前記データボリュームを前記固有ベクトルの選択されたサブセット上に投射して投射部分データボリュームを生成し、前記固有ベクトルのサブセットは、これらの対応の固有値に基づいて選択され、前記投射データボリュームに捕捉されない前記原データボリューム部分である残差データボリュームを求め、前記残差データボリュームは、前記ステップ（ｆ）において前記地下領域の物理的特徴を予測するために用いられる前記部分データボリュームである、請求項１４記載の方法。
【請求項１７】
前記異常の程度は、前記原データボリュームを前記固有ベクトルの選択されたサブセット上に投射して前記ステップ（ｆ）で用いられる前記部分データボリュームを生成することによって求められる、請求項１４記載の方法。
【請求項１８】
前記地下領域の前記予測地質学的特徴を用いてその石油が潜在的に存在しているか又は存在していないかを推定するステップを更に有する、請求項１記載の方法。
【請求項１９】
前記データ中の外れ値又は異常の識別は、（ｉ）前記データ値分布中の各データ窓の発生確率（又は同等の基準）をコンピュータ計算するステップと、（ｉｉ）低い確率データ領域を考えられる外れ値又は異常として識別するステップとを含む、請求項１記載の方法。
【請求項２０】
地下領域から炭化水素を産出する方法であって、
（ａ）前記地下領域の地球物理学的探査の結果を得るステップと、
（ｂ）少なくとも１つには、請求項１記載の方法を用いて識別された前記領域の物理的特徴に基づいて前記地下領域の潜在的な石油の存在の予測を得るステップと、
（ｃ）前記潜在的な石油の存在の確実な予測に応答して、坑井を前記地下領域中に掘削し、炭化水素を産出するステップとを有する、方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、一般に主に地球物理学的探査の分野に関し、特に地球物理学的データを処理する方法に関する。具体的には、本発明は、実世界の地質学的特徴（潜在的な炭化水素集積を含む）を表す、１つ又は２つ以上の地質学的又は地球物理学的データセット、例えば地震探査データセット中の領域を、事前トレーニングデータを使用することなく、強調表示する方法であって、所望の物理的特徴が未処理データ内に微妙な形でしか現れておらず、目立った異常によって見えにくくなっている場合に用いられる方法である。
【０００２】
〔関連出願の説明〕
本願は、２００８年１１月１４日に出願された米国特許仮出願第６１／１１４，８０６号及び２００９年７月３１日に出願された米国特許仮出願第６１／２３０，４７８号の権益主張出願である。
【背景技術】
【０００３】
地震探査データセットには、しばしば、複雑なパターンが含まれ、かかるパターンは、微妙であり、しかも多数の地震探査又は属性／派生ボリューム及び多数の空間的尺度で示される。ここ数十年間、地質学者及び地球物理学者は、炭化水素の存在を示す多くの重要パターンを抽出するために、様々な手法を開発してきた。しかしながら、これらの方法の殆どでは、１データボリューム又は多くて２ボリューム内での既定の特徴による既知のパターンか大まかに定義されたパターンかのいずれかの検索が必要である。これらの「テンプレート利用」又は「モデル利用」アプローチでは、かかる仕様に適合しない微妙な又は予想外の異常は、しばしば見落とされる。これらのアプローチについては、同じ技術的課題を取り扱うという点を除き、本発明とほとんど共通性がないため、ここではこれ以上検討しない。
【発明の概要】
【発明が解決しようとする課題】
【０００４】
これらの既知の方法の殆どは、解釈者が、事前に指定された特徴により、１データボリューム又は多くて２ボリューム内で既知のパターン又は大まかに定義されたパターンを検索することを伴う。これらの「テンプレート利用」又は「モデル利用」アプローチでは、かかる仕様に適合しない微妙な又は予想外の異常は、しばしば見落とされる。従って、地震探査データの１つ又は２つ以上のボリューム内の異常な領域を、それが何であるか又はどこにあるかに関する演繹的知識なしに、多数の空間的尺度にわたって自動的に強調表示することのできる統計的分析方法を開発することが望ましい。本発明は、この要望を満たす。
【課題を解決するための手段】
【０００５】
一実施形態では、本発明の要旨は、地下領域を表す地球物理学的データ又はデータ属性の１つ又は２つ以上の２Ｄ又は３Ｄ離散化セット（かかる各データセットを「原データボリューム」という）中の地質学的特徴部を識別する方法であって、（ａ）データ窓の形状及びサイズを選択するステップを有し、（ｂ）各原データボリュームに関して、窓を原データ中の複数個のオーバーラップ位置又は非オーバーラップ位置に動かして各データボクセルが少なくとも１つの窓内に含まれるようにし、各窓に関して、成分がその窓内からのボクセル値から成るデータ窓ベクトルＩを形成するステップを有し、（ｃ）データ窓ベクトルを用いて統計的分析を実施し、データボリュームに関する分布をコンピュータ計算するステップを有し、統計的分析は、複数個の原データボリュームの場合に共同して実施され、（ｄ）データ値分布を用いてデータ中の外れ値又は異常を識別するステップを有し、（ｅ）外れ値又は異常を用いて地下領域の地質学的特徴を予測するステップを有することを特徴とする方法にある。
【０００６】
次に、本発明の方法を用いて識別される地質学的特徴を用いると、炭化水素蓄積の存在を予測することができる。
【０００７】
本発明の内容及びその利点は、以下の詳細な説明及び添付の図面を参照すると良好に理解されよう。
【図面の簡単な説明】
【０００８】
【図１Ａ】本発明の方法の試験例の用途として、合成地震探査データの３Ｄボリュームからのイメージ（２Ｄタイムスライス）を示す図である。
【図１Ｂ】本発明の方法の試験例の用途として、情報の９０％を占める最初の１６個の主成分により定められる本発明により生成された原イメージの残差を示す図である。
【図１Ｃ】本発明の方法の試験例の用途として、上位（すなわち最初）の１６個の主成分を３０×３０窓の形で示す図である。
【図２】残差分析を用いた本発明の方法の一実施形態における基本的ステップの略図である。
【図３】単一の窓サイズを用いて本発明の窓ＰＣＡ実施形態を多くのデータボリュームに適用する際の基本ステップを示すフローチャートである。
【図４Ａ】窓中の互いに異なる画素に関するデータボリューム（大きな長方形）及びそのデータのサンプル（小さな長方形）の２Ｄスライスの略図であり、画素（１，１）に関するデータサンプルを示す図である。
【図４Ｂ】窓中の互いに異なる画素に関するデータボリューム（大きな長方形）及びそのデータのサンプル（小さな長方形）の２Ｄスライスの略図であり、ｉ番目の画素に関するデータサンプルを示す図である。
【図５Ａ】共分散行列の効率的なコンピュータ計算を行なうために図４Ａ及び図４Ｂの２Ｄデータセットに関してサンプルに含まれないデータの細分結果を示す図である。
【図５Ｂ】共分散行列の効率的なコンピュータ計算を行なうために図４Ａ及び図４Ｂの２Ｄデータセットに関してサンプルに含まれないデータの細分結果を示す図である。である。
【発明を実施するための形態】
【０００９】
図１Ａ〜図１Ｃ及び図２は、カラーディスプレイの黒色及び白色再現結果を示す図である。
【００１０】
本発明を例示の関連実施形態と関連して説明する。以下の説明が本発明の特定の実施形態又は特定の使用に特有である程度まで、これは、例示であるに過ぎず、本発明の範囲を制限するものと解されてはならない。これとは対照的に、本発明は、特許請求の範囲に記載された本発明の範囲に含まれる全ての変形例、改造例及び均等例を含むものである。
【００１１】
本発明は、事前トレーニングデータを使用することなく、マルチボリュームの地震探査又はその他の地球物理学的データ（例えば電磁データ）内の異常なパターンを多数の空間的尺度にわたって検出する方法である。本発明の方法は、窓型統計分析に基づくものであり、本発明のある実施形態では、次のような基本ステップを含む。
１．ユーザが指定したサイズ及び形状の窓の範囲内のデータについて統計的分布を抽出するステップ。主成分分析（Principal Component Analysis：ＰＣＡ）、独立成分分析（Independent Component Analysis：ＩＣＡ）、クラスター化分析（Clustering Analysis）のような標準的な統計手法を用いてもよい。
２．（ａ）抽出された分布での各データ窓の発生確率（又は同等の測定基準）をコンピュータ計算し、（ｂ）低確率データ領域を予想される異常として識別することにより、データ内の異常な領域を抽出するステップ。
【００１２】
本発明の特に簡便な実施形態は、窓主成分分析（「Window Principal Component Analysis ：ＷＰＣＡ」）、残差分析（Residual Analysis ）及びクラスター化分析の組み合わせ（下記に詳細に説明する）を伴うものである。しかしながら、当該技術分野の当業者は、その他の統計的分析手法をどのように使用し又は好適に適合させれば同じ目標を達成できるか、その方法を容易に正しく理解するはずである。
【００１３】
主成分分析（「ＰＣＡ」）の有用な一般化は、独立成分分析（「ＩＣＡ」）として知られ、データが標準的な多次元ガウス分布と大きく異なる場合に好ましい方法である。この場合、本発明の方法は、それに応じて窓ＩＣＡ（「ＷＩＣＡ」）を用いるために一般化され、次いで残差分析が一般化される（「外れ値検出（Outlier Detection ）」と呼ばれる）。一実施形態では、本発明は、動作窓に対してＰＣＡを用い、次いで主成分（「ＰＣ」）から内積及びデータ残差をコンピュータ計算する。これは、地震探査用途だけでなく、より広範な多次元データ処理分野にわたっても有利に適用できると考えられる。これは、画像、音声及び信号処理の分野を含む。
【００１４】
主成分分析（“ＰＣＡ”）は、最初にペアソン（Pearson）（「オン・ラインズ・アンド・プランズ・オブ・クローゼスト・フィット・トュー・システムズ・オブ・ポインツ・イン・スペース（On Lines and Planes of Closest Fit to Systems of Points in Space）」，フィロソフィカル・マガジン（Philos. Magazine），第２巻，１９０１年，ｐ．５５９‐５７２）により提案され、更にホテリング（Hotelling）（「アナリシズ・オブ・ア・コンプレックス・オブ・スタティスティカル・ヴァリアブルズ・イントュー・プリンシパル・コンポーネンツ（Analysis of a Complex of Statistical Variable Into Principal Components）」，ジャーナル・オブ・エデュケーション・フィコロジー（Journal of Education Psychology），第２４巻，１９３３年，ｐ．４１７‐４４１）により発展したデータ分析のための周知な古典的な技術である。地震探査データへの主成分分析の最初の公知の適用であると考えられているものは、カリ・カヒューネン（Kari Karhunen）及びマイケル・ローブ（Michel Loeve）（（ワタナベ（Watanabe），「カーヒューネン‐ローブ・エクスパンション・アンド・ファクター・アナリシス（Watanabe "Karhunenh‐Loeve Expansion and Factor Analysis"）」，トランザクションズ・オブ・ザ・フォース・プラーク・コンフェレンス（Transactions of the Fourth Prague Conference），ジェイ・コゼンスニク（J. Kozensnik）編集，プラーク（Prague），チェコスロヴァキア・アカデミー・オブ・サイエンス（Czechoslovakia Academy of Science），１９６７年）に因んでカーヒューネン‐ローブ変換の形態で得られたものである。この方法は、ＰＣＡを用いて１組の地震探査トレース中の情報内容を説明しており、入力データセットの形態は、可変サイズの多次元窓ではなく、地震探査トレース全体である。ワタナベ（Watanabe）の主要な用途は、地震探査トレース全体を分解し、最初の７つの主成分トレースを用いて殆どのコヒーレントエネルギーを再構成し、それにより非地質学的ノイズを除去することであった。
【００１５】
ＰＣＡは、大抵の場合、多くの測定特性を統計的に独立した１組の属性にセットに減少させる地震探査分析において一般的に用いられている（これについては、例えば、フォーニア（Fournier）及びディレイン（Derain），「ア・スタティスティカル・メソドロジー・フォー・デライヴィング・リザーバ・プロパティーズ・フロム・セイスミック・データ（A Statistical Methodology for Deriving Reservoir Properties from Seismic Data）」，ゲオフィジクス（Geophysics），第６０巻，１９９５年，ｐ．１４３７‐１４５０及びハーゲン（Hagen），「ジ・アプリケーション・オブ・プリンシパル・コンポーネンツ・アナリシス・トゥー・セイスミック・データ・セッツ（The Application of Principal Components Analysis to Seismic Data Sets）」，ジオエクスプロレーション（Geoexploration），第２０巻，１９８２年，ｐ．９３‐１１１を参照されたい）。
【００１６】
地震解釈プロセスは、原（オリジナル）データからの多くの派生生成物を生じさせる場合が多い。これら属性は、様々な程度に相関するので、ＰＣＡは、多量の情報を保持しながら多くの特性を減少させる巧妙な手法であった。
【００１７】
今日まで、地質学的又は地球物理学的データにおけるスコーピング及び踏査に関して関心のある地質学的特徴を見出すのに専用に用いられる非移動式窓利用型統計的外れ値検出技術が存在すると考えられている。しかしながら、かかる技術は、専用信号処理又は貯留層特徴付け用途のための地震探査データの特定のサブセット又はドメインに利用されている。ケイ（Key）及びスミソン（Smithson），（「ニュー・アプローチ・トゥー・セイスミック・リフレクション・イヴェント・ディテクション・アンド・ヴェロシティ・デタミネーション（New Approach to Seismic Reflection Event Detection and Velocity Determination）」，ジオフィジクス（Geophysics），第５５巻，１９９０年，ｐ．１０５７‐１０６９）は、ＰＣＡをあらかじめプレスタック地震探査データ中の２Ｄ移動窓に適用し、信号コヒーレンシー可干渉性の尺度として合成固有値の比を取っている。プレスタック地震探査データ中の特徴を検出するのに主成分それ自体は用いられていない。シーベル（Sheevel）及びペイラジアン（Payrazyan），（「プリンシパル・コンポーネント・アナリシス・アプライド・トゥー・スリーディー・セイスミック・データ・フォー・リザボイア・プロパティ・エスティメーション（Principal Component Analysis Applied to 3D Seismic Data for Reservoir Property Estimation）」，ソサイエティ・オブ・ペトロレイウム・エンジニアズ・アニュアル・コンファレンス・アンド・エキジビション（Society of Petroleum Engineers Annual Conference and Exhibition），１９９９年）は、小さな１Ｄ移動垂直窓を用いてトレースに基づく主成分を利用した主成分をコンピュータ計算し、坑井較正から遠ざかる貯留層特性を予測する分類アルゴリズム中に大抵の地質学的特徴を見るＰＣを入力する。この場合も又、この１Ｄ単一データセットアプローチは、そのデータ中の異常又は外れ値を自動的には識別するようにはなっていない。チョー（Cho）及びスペンサー（Spencer）（「エスティメーション・オブ・ポララリゼーション・アンド・スロウネス・イン・ミックスド・ウェーブフィールズ（Estimation of Polarization and Slowness in Mixed Wavefields）」ジオフィジクス（Geophysics），第５７巻，１９９２年，ｐ．８０５‐８１４）及びリッチウォルスキー等（Richwalski et al.）（「プラクティカル・アスペクツ・オブ・ウェーブフィールド・セパレーション・オブ・ツー・コンポーネント・サーフェイス・セイスミック・データ・ベースド・オン・ポララリゼーション・アンド・スロウネス・エスティメーツ（Practical Aspects of Wavefield Separation of Two-Component Surface Seismic Data Based on Polarization and Slowness Estimates）」，ジオフィジカル・プロスペクティング（Geophysical Prospecting），第４８巻，２０００年，ｐ．６９７‐７２２）は、２Ｄ窓ＰＣＡを周波数ドメイン内で用いて既定数のＰ波及びＳ波の伝搬をモデル化している。
【００１８】
ウー等（Wu et al.）（「エスタブリッシング・スパティアル・パターン・コリレーションズ・ビトゥイーン・ウォータ・サチュレーション・タイム‐ラプス・アンド・セイスミック・アンプリチュード・タイム‐ラプス（Establishing Spatial Pattern Correlations Between Water Saturation Time-Lapse and seismic Amplitude Time-Lapse）」，ペトロレウム・ソサイアエティーズ・シックスス・アニュアル・カナディアン・インターナショナル・ペトロレウム・カンファレンス（Petroleum Society's 6th Annual Canadian International Petroleum Conference），（フィフティーシックスス・アニュアル・テクニカル・ミーティング（56th Annual Technical Meeting）），２００５年）の目的は、単一地震探査又は時間経過地震探査ボリュームを貯留層モデル中のフローシミュレーションデータに最適に相関させて空間パターンの実際の飽和時間経過値を推定することにある。このアプローチは、原データボリュームに関してではなく、ＰＣＡ分析からの第１の主固有ベクトル上へのこれらデータの投射に関して２点間比較を実施することにある。かくして、これら目的は、地震探査データ中の異常パターンの識別ではなく、既知のモデルに対する地震探査データの相関である。
【００１９】
ビショップ（Bishop）に付与された米国特許第５，８４８，３７９号明細書（発明の名称：（Method for Characterizing Subsurface Petrophysical Properties Using Linear Shape Attributes 」，１９９８年）は、本発明によって取り組まれる技術的課題であるスコーピング及び踏査方式で関心のある地質学的特徴を識別するのではなく、地下岩石特性を予測してファシイ又はテキスチャ分析のための地震探査データを分類する方法を開示している。ビショップは、地震探査トレースをあらかじめ指定された時間又は深さの間隔内でリニア・シェイプス（Linear Shapes）と呼ばれている直交波形ベースの一次結合に分解するためにＰＣＡを用いて統計的分析を実施する。線形属性（Linear Shape Attribute：ＬＳＡ）は、特定のトレース形状を再構成するために用いられる重み（又は固有値）のサブセットとして定義される。また、ビショップは、窓をオーアーラップさせ、多くのデータボリュームを同時に分析し又は統計的分布を用いて異常のあるデータ領域を検出することを開示しているわけではない。
【００２０】
地質学的データ及び地球物理学的データを統計的に分析する他のアプローチは、例えば人工神経回路網、遺伝的アルゴリズム及びマルチポイント統計学的手法のような方法を用いるが、異常のあるパターンの自動的検出を目的としているわけではない。加うるに、これら方法は、一般に、成功の度合いが制限されている。というのは、これらの内部作業は、不明瞭であり、しかも、これら方法は、多量のトレーニングデータに非常に依存しているからである。
【００２１】
上述のように、ＰＣＡ及びＩＣＡは、高次元（即ち、多変量又は多属性）信号を分離し統計的に無相関の（すなわち独立した）成分にするために一般的に用いられる方法である。本発明の窓ＰＣＡ及び窓ＩＣＡでは、原データの各点をその近傍の点集合（即ち、窓）として表すことにより原データから導き出されたデータセットに対して成分分析が適用される。図３のフローチャートを参照してこの概念を説明するために、単一の３次元データボリュームに対する固定窓サイズを用いたＷＰＣＡの実施を、以下に概要説明する。同じ手順又はそのＩＣＡ相当物を、２Ｄデータに対して、又は多数の２Ｄもしくは３Ｄデータボリュームに対して同時に適用できるはずである（図３のステップ３１参照）。Ｎx×Ｎy×Ｎzのサイズの３Ｄ地震探査ボリュームについて考察する。
【００２２】
（ステップ３２）窓の形状（例えば楕円体又は直方体）及びサイズ（例えば、半径ｒ、ｎx×ｎy×ｎz）を選択する。
【００２３】
３Ｄ地震探査ボリュームＩi,j,k内の各ボクセルは、ｎx×ｎy×ｎzの次元ベクトルベクトルＩi,j,kとして表され、かかるベクトルは、各ボクセルの窓近傍内のボクセル値を包含する。
【００２４】
（ステップ３３）全てのｎ次元ベクトル
0005530452.tif 000002
の平均及び共分散行列を次のようにコンピュータ計算する。
0005530452.tif 000003
【００２５】
相関行列を
0005530452.tif 000004
としてコンピュータ計算する。上式において、ｔ及びｋは、ベクトルＩの２つの指数であり、従って、３次元における２つの異なる空間座標セットを表す。
【００２６】
（ステップ３４）固有値（主値）｛λ1＞λ2＞…λn｝及び固有ベクトル（主成分）
0005530452.tif 000005
の｛ν1，ν2…，νn｝
をコンピュータ計算する。変形例として、共分散行列の固有値をコンピュータ計算してもよく、かかる数値は、尺度係数分だけ相関行列の固有値と異なる。これらの固有ベクトルは、ｎx×ｎy×ｎzのサイズとなり、ベクトル型から再び窓型に戻すと、データ内の様々な（独立した）空間パターンを、共通性が最大のものから最小のものまで順番に表す。対応する固有値は、どれほど多くの原データを各固有値が占めているか（すなわち分散量）を表す。
【００２７】
下記の地震探査又は属性データの部分ボリュームの１つ又は２つ以上を生成する。その後、原データボリュームからは明らかでなかった異常があるかどうかについてかかる部分ボリュームを検査する。
（ａ）（ステップ３５）投射：各主成分又は主成分群（例えば、クラスター化分析から選択されるもの）を用いて再現できる原データの一部。これは、各主成分又は主成分群について平均中心を求めて標準化した地震探査ボリュームの内積を取ることにより達成される。従って、ベクトルＢに対するベクトルＡの投射は、ｐｒｏｊ（Ａ）＝（Ａ・Ｂ）Ｂ／｜Ｂ｜2を意味し、Ｂ方向のベクトルである。
（ｂ）（ステップ３６）残差：最初のｋ−１（すなわち共通性が最大の）主成分により捕捉されていない原ボリューム内の残差信号。本発明の好ましい実施形態では、これは、平均中心を求めて標準化した地震探査ボリュームを、｛ｖk，ｖk+1，・・・，ｖn｝により張られた部分空間に投射することによって達成され、その結果、次の通りである。
0005530452.tif 000006
上式において、Ｒは、０〜１の範囲内でユーザが定義した閾値である。変形例として、あるいは、予測値をボトムアップ方式で加えることもできるはずであるが、これは、殆どの場合、コンピュータ計算をより厄介にするはずである。
（ｃ）外れ値：項目（ｂ）の残差分析は、本発明のある実施形態において各ボクセルの「異常の程度」を割り出すやり方である。各ボクセルの「異常の程度」をコンピュータ計算する代替的なやり方では、（ａ）及び（ｂ）の属性データボリュームは必要とされず、（上記の残差Ｒと関連するものの同じではないため）Ｒ′で示され、次式により求められる。
0005530452.tif 000007
【００２８】
この異常の程度の測定を用いて、部分データボリュームが展開される。また、この測定は、最初の少数の固有値により張られた空間内にある「外れ値」も拾い上げるが、場合によっては、上記の２つのステップに比べて、コンピュータ計算がより集約的になり得る。しかしながら、この場合、上記のステップ３４は、これを省くか、又は相関行列のコレスキー分解と簡単に置き換えることができ、より迅速にＲ′を求められるようになることに注目されたい。
【００２９】
上記の基本的アプローチには、異なるデータ標準化方式を用いる変形が存在する。この方法は、任意の数の地震探査ボリュームに拡大適用することができる。ユーザが試みることのできる調整可能なパラメータは、（１）窓形状、（２）窓サイズ、及び（３）残差予測値の閾値Ｒである。
【００３０】
地震探査データの２次元スライスに対して３×３ＷＰＣＡを適用した結果が図１Ａ〜Ｃに示されている。図１Ａは、合成３Ｄ地震探査データボリュームからのイメージ（２Ｄタイムスライス）を示す。実地では、この表示は通常、色で行われ、地震反射波の振幅を表示する（例えば、青＝プラス、赤＝マイナス）。図１Ｂは、最初の１６個の主成分が情報の９０％を占めた後の原データイメージの残差を示す。留数は、異常パターンでは高数値で、この場合には断層である。図１Ｂのカラー版では、少量の残差を青で示してもよく、図１Ｂの残差表示に今や明確に見て取れる異常な断層系を暖色で強調表示してもよい。図１Ｃには、上位（すなわち最初）の１６個の主成分１４が、３０×３０窓の形で示されている。断層が下２列内の数個の主成分に捕捉されていることが見て取れる。
【００３１】
２次元合成地震探査断面に対して９×９ＷＰＣＡを適用した結果が図２の概略的フローチャートに示されている。合成３Ｄ地震探査データボリュームからの２Ｄ断面が符号２１のところに表示されている。通常は、色を用いて地震反射波の振幅を表す。８ｍ秒と小さく、微妙すぎて目では検出できない背斜が、背景の水平反射率内に埋め込まれている。入力イメージの最初の４つの主成分（固有ベクトル）を、２２に表示する。表示２３は、情報の９９％を占める最初の４つの主成分に対する原イメージの投射又は予測値を示す。表示２４は、原イメージから投射イメージを取り去った後の残差を示す。埋め込まれていた微妙な特徴が、今や、約４４０ｍ秒の深さ（往復走時）のトレース番号（横方向位置を１方向に測定して）３０〜５０の間に現れている。カラー表示では、「暑い」色を用いて、埋め込まれた微妙な特徴の位置を明らかにしてもよい。
【００３２】
図３のフローチャートは、本発明の方法の実施形態の概要を説明したものであり、単一の窓サイズを用いて多数のデータボリュームにＷＰＣＡが適用されている。
【００３３】
標準的パターンの構築における一般化及び効率
【００３４】
以下の項では、本発明の窓主成分分析の改良について説明する。かかる改良により、コンピュータ計算の減少を通じて、適用をより簡便にできるとともに、主成分又は独立成分の解釈並びに選択的な保持又は除去を通じて、結果の使用を改善できる。
【００３５】
コンピュータ計算効率：上記の共分散行列をコンピュータ計算する単純な方法は、大きいデータセットについてはメモリ要求事項とプロセッサ要求事項との両面で、コンピュータ計算が厄介である。従って、本願では、次のような事実を活用する代替方法を開示する。すなわち、ＰＣＡの個々のベクトルは、データを横切って移動する窓である。例えば、｛Ｉ1，Ｉ2，…，ＩN｝の数値を有する１Ｄデータセットについて考えてみる。Ｋ＜Ｎのサイズの窓の共分散行列を求めるためには、入力値の平均及び二次モーメントを、次のようにコンピュータ計算することができる。
0005530452.tif 000008
【００３６】
この方法が、データのサブベクトル（より高次元の部分行列）の平均及び内積を取ることだけを伴うものであり、よって、原データから導き出される数々の小サイズの窓の保存及び操作を行わずに済むことに注目されたい。従って、このコンピュータ計算方法の修正により、オブジェクト指向のソフトウェアは、効率的なアレイインデキシング（array indexing）（マットラボ（Matlab）及び「サムド−エリア・テーブルズ・フォア・テキスチャ・マッピング（Summed-Area Tables for Texture Mapping）」コンピュータ・グラフィックス１８（Computer Graphics 18），２０７（１９８４年）においてクロウ（Crow）が説明しているデータ構造であるエリア総和テーブルの使用）により、最小限の記憶努力及びコンピュータ計算努力で共分散行列をコンピュータ計算できるようになる。
【００３７】
あるいは、共分散行列のコンピュータ計算を、徐々に小さくなる領域に対する一連の相互相関演算として表すことにより、コンピュータ計算効率を確保してもよい。このアプローチを説明するために、図４Ａ〜Ｂに記載の２次元データセットをｎ＝ｎx×ｎyのサイズ、また、二次元窓をｍ＝ｍx×ｍyのサイズと考えてみる。その時、相関行列Ｗ（ｔ，ｋ）は、最初に各標本データの平均をコンピュータ計算し次に内積行列をコンピュータ計算しその後かかる行列を標準化して平均を差し引くことにより、得られる。
【００３８】
まず、いずれも１／（ＤＳ１内のピクセル数）に等しい入力値から成る標本データ（例えばＤＳ１）のサイズのカーネルで、データボリュームを畳み込むことにより、平均をコンピュータ計算することができる。この演算の結果は大きい行列になるが、平均は、当該出力の左上の角に配置されたサイズｍの窓に位置付けられた数値である。一般に、この種の演算は、ｃｏｒｒＷ（カーネル，データ）で示され、その結果が、上記のように読み出されたサイズｍの窓である。高速フーリエ変換（ＦＦＴ）を用いる演算の実行は、ｎ×ｌｏｇ（ｎ）に比例する時間を要し、標本化窓のサイズとは無関係である。ｍがｌｏｇ（ｎ）より十分大きい場合、このＦＦＴアプローチは、陽関数表示アプローチより迅速である。
【００３９】
次に、データセットの副標本に対して一連のｃｏｒｒＷ演算を行うことにより、内積行列Ｕ（ｔ，ｋ）をコンピュータ計算する。この行列の行ｉは、Ｕ（ｉ，：）で示され、Ｕ（ｉ，：）＝ｃｏｒｒＷ（ＤＳｉ，データ）としてコンピュータ計算することができる。よって、このやり方による行列のポピュレートは、ｍ×ｎｌｏｇ（ｎ）に比例して又はそれ以上に時間を要する。ただし、データの様々な部分領域に対して数回のｃｏｒｒＷ演算を行うことによりＵ（ｔ，ｋ）をコンピュータ計算するほうが有利である。具体的には、次式のように書き直せる。
0005530452.tif 000009
上式において、ｃｏｒｒＷ（データ，ＤＮＳｉ）は、ＤＮＳｉとＤＮＳｉの近傍のデータとの相互相関であって、ＤＮＳｉの位置がｍx又はｍyの範囲内のものを示す。演算ｃｏｒｒＷ（データ，データ）を、全ての行について１回だけ行う必要があり、その後、ｃｏｒｒＷ（データ，ＤＮＳｉ）をｍ回、コンピュータ計算する必要がある。利点は、次のような事実からもたらされる。すなわち、ＤＮＳｉは通常、データセットのサイズよりはるかに小さく、そのため、ｃｏｒｒＷ（データ，ＤＮＳｉ）は、ｃｏｒｒＷ（データ，データ）に比べて、はるかに小さい入力に関する相互相関である。同様に、ｃｏｒｒＷ（データ，ＤＮＳｉ）のコンピュータ計算は、いっそう小さい部分領域に対する数回のｃｏｒｒＷ演算に分けることができる。
【００４０】
ＤＮＳｉの大部分は、異なる標本についても同じであり、１回につき標本化窓の１方向に沿った部分が異なるだけである。例えば、図５Ａ〜Ｂの説明図について考えてみる。図５ＡにＡ，Ｂ及びＣで示した領域は、１つにまとめると、ピクセル１により標本化されていないデータボリュームの全エリアを形成する。これは、さらに細分してコンピュータ計算数を減らすことのできるエリアである。Ａ及びＣにより張られた「垂直方向」エリアについて考え、図５Ｂに記載の異なる標本化領域ＤＳｉと比較してみる。類似した垂直方向エリアが、数個の小領域の和集合、Ｃ１＋Ｃ２＋Ｃ３＋Ｃ４＋Ａ１＋Ａ２により張られている。（図５Ａの領域Ｂに相当する分割は、図５Ｂでは和集合Ｂ１＋Ｂ２である。）一般に、かかるあり得るエリアは、ｍxだけ異なる形で、各々がＤＳｉの一意的な横方向位置に対応している。換言すると、Ａ＋Ｃ内に包含されるデータは、異なる多くの標本データＤＳｉについて同じであり、従って、ｍx回操作するだけでよく、このエリアに関してｍy回のコンピュータ計算の節約となる。従って、ｃｏｒｒＷ（データ，ＤＮＳｉ）のコンピュータ計算は、このやり方で最適化して次式に従ってコンピュータ計算してもよく、
0005530452.tif 000010
上式において、文字で示した領域は、当該文字及び数字で表示した全ての領域の和集合を意味する。例えば、式中のＣは、図５Ａでは領域Ｃを、図５ＢではＣ１＋Ｃ２＋Ｃ３＋Ｃ４を指し、そのため、Ａ＋Ｃは、図５ＢではＡ１＋Ａ２＋Ｃ１＋Ｃ２＋Ｃ３＋Ｃ４を表す。ｃｏｒｒＷ（データ，Ａ＋Ｃ）のコンピュータ計算は、Ｕ（ｔ，ｋ）のｍy行について１回だけ、また、ｃｏｒｒＷ（データ，Ｂ＋Ｃ）についても同様に行えば済むため、各行についてコンピュータ計算する必要のある唯一の部分は、ｃｏｒｒＷ（データ，Ｃ）である。効率確保は、Ｃで示した領域が通常、その他の領域に比べて著しく小さいという事実からもたらされる。このやり方を進めると、アルゴリズムは、３Ｄデータセット及び窓まで（そして、実際には何次元までであろうと）拡大適用される。
【００４１】
最後に、行列Ｕを適切に標準化し平均を差し引くことにより、相互相関行列Ｗ（ｔ，ｋ）を得る。
0005530452.tif 000011
上式において、ｎＤＳは、各標本データ内の要素の数である
【００４２】
マスクの使用：非常に大きいデータセットについては、上記のコンピュータ計算効率でさえ、利用可能なコンピュータ計算資源により時宜を得たやり方で結果を得るためには十分でない可能性がある。このような場合には、事前定義マスクに対して（ａ）固有ベクトルによる内積コンピュータ計算、又は（ｂ）主成分コンピュータ計算を適用することができる。マスクは、コンピュータ計算が行われる対象のデータの空間サブセットである。マスクは、（ａ）ユーザが対話形式で生成してもよく、又は（ｂ）派生属性を用いて自動的に生成してもよい。（ｂ）の実施例は、勾配推定アルゴリズムを用いての、局所勾配の高いデータ領域の事前選択になるはずである。内積コンピュータ計算は、主成分コンピュータ計算より厄介であり、これにより、必要に応じて一方又は両方のコンピュータ計算への、マスクの適用が動機付けられる。
【００４３】
カノニカルパターンの用途
【００４４】
さらに、コンピュータ計算された主／独立成分は、テキスチャ、カオス又は他の特性により測定された類似のパターンを表示するグループにクラスター化されるのが良い。残差ボリュームと一緒に、個々の又はグループをなす主成分上への原地震探査データの投射又は予測値により、異常のあるパターンが強調表示された多くの地震探査ボリュームが導き出される。本発明のこれら実施形態について以下に詳細に説明する。
【００４５】
多数の窓／空間スケール：更に、多数の入れ子型窓サイズのための共分散行列を一度にコンピュータ計算する安易な仕方と比較して、かかる共分散行列を階層的順序でコンピュータ計算する技術的努力を合理化することが可能である。この場合も又、２つの窓サイズＫ1＜Ｋ2の一次元例を考察する。Ｋ2についての平均及び第２のモーメントを先ず最初に上述した方法の使用によってコンピュータ計算し、その後、Ｋ1について同じ量を次のようにコンピュータ計算するのが良い。
0005530452.tif 000012
【００４６】
上述の公式により、増分的労力により小さな窓に関する量をコンピュータ計算することができるということに注目されたい。この方法を高次元の窓の入れ子型シリーズに拡張することが容易である。
【００４７】
主成分及び投射の利用：本発明により生成される主成分及び投射を利用し、組み合わせ、そして視覚化することができる多くの方法が考えられる。好ましい一具体化例では、上述した残差を用いて異常の識別を行なう。同様に有効なアプローチは、ＰＣの選択されたサブセット上への原データの選択された投射を実施することである。サブセットは、（ａ）ユーザにより対話的にか（ｂ）ＰＣに対するコンピュータ計算基準を自動的に用いるかのいずれかによって選択されるのが良い。（ｂ）の一例は、自動幾何学的アルゴリズムを用いて「チャネル」又は管状構造物に似た特徴を有するＰＣの選択である。別の例は、ノイズ検出アルゴリズム又は分散基準を用いて「ノイズのある」ＰＣを排除する投射を作ることによって入力データ中のノイズを減少させることである。当業者であれば、この説明から他の例を認識されよう。
【００４８】
種々の窓サイズでの投射結果を視覚化する別の有用な仕方は、（ａ）ＰＣ投射のユーザにより又は自動的に選択される組み合わせ、（ｂ）種々の残差閾値のところでの残差又は（ｃ）ノイズ成分の視覚化を含む。もう１つの有用な変形例は、「分類ボリューム」の視覚化を含み、これは、各データ場所をこの場所での最も高い値を持つＰＣ投射がどれであるかを一義的に定める色で色分けすることを含む。
【００４９】
反復ＷＰＣＡ：図３に概要説明された仕事の流れにより作られる残差ボリュームがより異常のあるパターンを含む領域中の大きな値を示すということが判明した。その結果として、入力データ中の微妙なパターンは、残差ボリューム中の自明な異常によって隠される場合が多い。ＷＰＣＡの感度を極めて微妙なパターンに増大させるため、２つの択一的な反復アプローチを用いるのが良い。
【００５０】
反復固有ベクトル除去：この第１の択一的な手順は、以下のステップを含むのが良い。  １．図３のフローチャートの最初の４つのステップを実施する（固有ベクトル及び固有値生成により）。
２．投射が多量のバックグラウンド信号及び大抵の自明な異常を再構成する固有ベクトルを識別する。
３．データを先のステップにおいて識別されなかった固有ベクトルのサブセット上にのみ投射する（バックグラウンド信号及び大抵の自明な異常の信号は、この投射像中で減衰されるはずである）。
４．ＷＰＣＡを先のステップで生成した投射像に対して実施する。
５．必要に応じてステップ１〜３を繰り返す。
【００５１】
反復マスキング又はデータ除去：この第２の択一的手順は、次のステップを含むのが良い。
１．図３の最初の４つのステップを実施する（固有ベクトル及び固有値生成により）。  ２．種々の残差ボリュームを検査することにより、入力データ中の大抵の自明な異常に相当する領域を識別する。
３．ＷＰＣＡをデータに対して実施し、識別された領域を
ａ．ＷＰＣＡによる分析に先立ってこれら領域中の全ての属性値をゼロに設定し、又は
ｂ．ＷＰＣＡへの入力としてこれら領域を含まないようにする。
４．ＷＰＣＡを新たなデータセットに対して実施する。
５．必要に応じてステップ１〜３を繰り返す。
【００５２】
ＷＰＣＡ分類：主成分は、投影の強度に基づいて像を分類するために用いられるのが良い。かかる分類は、特に原データが多数のボリュームから成る場合、簡便な視覚化により、選択された主成分内に示された特定のパターンを有する領域を識別するのを助ける。この変形例は、以下のステップを有するのが良い。
１．図３の最初の４つのステップを実施する（固有ベクトル及び固有値生成により）。  ２．データ中の各点にその点周りに窓内の大抵の信号を再構成する固有ベクトルに対応した番号を割り当てる。これは、各点が１（即ち、最初の固有ベクトル）とＮ＝ｎx×ｎy×ｎz（即ち、最後の固有ベクトル）との間の数を含む分類されたボリュームを構成する。
３．次に、各値（又は値の群）を１−Ｎから独特の色又は透明度（これらの組み合わせ）を割り当てることにより分類結果を視覚化する。この手順は、Ｎ次元画像のパターン利用型分類の一形態である。連続スペクトル残差又は投射又は予測値ではなく、依然として投射画像中の信号の大きさに基づいてカテゴリを出力することにより、この手順は、微妙な特徴に対する感度の低下という影響を受けにくくなる。
【００５３】
かくして、本発明の方法は、大きな高次元データセット、例えば地震探査データから特徴を抽出するのに有利である。例えば、ＰＣＡを地震探査データに適用する大抵の文献記載方法は、これらがデータ窓に関する固有モード分解を実施するという点においてのみ本発明の方法に似ている。一例は、上述したウー等（Wu et al.）の方法である。これらのアプローチは、幾つかの基本的な仕方において本発明とは異なっている。先ず最初に、これらは、小さな１Ｄの垂直に動く窓だけをＰＣＡに対する入力として地震探査データに適用する。３Ｄ移動窓は、フローシミュレーションデータにのみ用いられる。第２に、第１のＰＣは、時間経過地震探査データとフローシミュレーションデータの両方を再構成するために用いられるに過ぎない。他の投影又は数学的組み合わせ、例えば残差ボリュームの再構成を実施する。最後に、多数の地震探査ボリュームを同時に検査し、地震探査データに固有のパターンを抽出させる（即ち、既に存在している地質学的モデルには制約されない）試みはなされない。
【００５４】
上述の開示内容は、本発明を例示する目的で本発明の特定の実施形態に関する。しかしながら、本明細書において説明した実施形態の多くの改造例及び変形例が可能であることは当業者には明らかであろう。かかる改造例及び変形例は全て、特許請求の範囲に記載された本発明の範囲に含まれるものである。
【図１Ａ】
0005530452.tif 000015
【図１Ｂ】
0005530452.tif 000016
【図１Ｃ】
0005530452.tif 000017
【図２】
0005530452.tif 000018
【図３】
0005530452.tif 000019
【図４Ａ】
0005530452.tif 000020
【図４Ｂ】
0005530452.tif 000021
【図５Ａ】
0005530452.tif 000022
【図５Ｂ】
0005530452.tif 000023