JP6453892B2 - 二次尤度汎関数を用いてデータの統計モデリングを行う方法およびシステム - Google Patents

Description

関連出願の相互参照
本出願は、２０１３年９月１８日出願の米国仮特許出願第６１／８７９，６２９号明細書を優先権主張するものである。

本発明は、対数尤度関数を代替すべくデータの予測、特に推定された係数の個数がデータ点の個数と同等かまたはより大きい非漸近的統計的推定問題において統計的推定を行う方法、より具体的には二次尤度汎関数に関する。

検出画像、受信信号、世論調査、検索、または他の任意のデータ収集方法からの測定データはノイズまたは不確実性の影響下にあり得る。当該方法は、データを過度に解釈することなく、なるべく多くの情報をデータから抽出する能力を向上させる必要がある。

統計的推定は、ノイズが多いデータ、すなわちデータの実現値毎に変化する未知のランダムな成分を有し、且つ統計的にしか特徴付けできないデータから情報を抽出する必要がある場合に用いられる。統計的推定の目的は、データ内の信号、すなわちデータの実現値毎に変化しない再現可能な成分をモデリングすることである。信号は、データに対するモデルの適合度により値が決定されるパラメータの関数としてモデリングされる。

ノイズが多いデータをモデリングする処理には広範な用途があり、いくつか例を挙げれば対象追跡、信号処理、画像処理（ＣＴ、ＳＰＥＣＴ、ＰＥＴ、Ｘ線等の医療画像処理を含む）、市場調査、サプライチェーン管理、在庫管理、および金融市場が含まれるがこれらに限定されない。履歴データから抽出された情報は往々にして将来的な動向および／または付随するリスクを予測するために用いられる。基本的なモデルが正しい限り、予測および／またはリスク解析の品質は推定の精度に応じて決定される。

統計的推定が広義の最適化の分野に関係し、これには非統計的な方法が含まれ、ノイズが多いデータが関与せず、最適化は単に最良のパラメータの発見する問題に過ぎないことを指摘せねばならない。

統計的推定の起源は、最小二乗（ガウス、１８０９年）法まで２世紀遡り、後に最大尤度（ＭＬ）法（フィッシャー１９１２、１９２２年）に発展した。統計的に独立したデータの確率分布関数（ｐｄｆ）を与えられたならば、ＭＬ法は、モデルが前提として、当該データの条件付き確率を最大化するかまたは、同等的に、対数尤度関数（ＬＬＦ）を最小化する。

ここに、ｘはｎ次元空間であっていずれの次元も連続的、離散的、更には圏論的であってよく、θはモデルパラメータ、ｘ_ｉは観測値、ｆ（ｘ，θ）はｐｄｆ、δはｎ次元のディラック（１９５８年）デルタ関数である。積分記号は、連続的な次元についての積分、および離散的および圏論的次元にわたる和を表し、但し実用的には積分を和で近似してよいものと理解されたい。

多くの用途において、ｐｄｆの積分的正規化は典型的には単位元に固定されている。この場合、式（１）の右辺第１項は定数であって省略可能であり、次式が得られる。

追加的な項が式（１）に含まれているのは、例えば検出器により観測された事象の発生率または製品の販売率等、観測値の正規化自体が関心対象の場合を考慮に入れるためである。その場合、式（１）はポアソン（１８３７年）分布の非ビン化ＬＬＦである。

ＭＬ法には３個の明確な利点がある。
１．データのビン化を必要とせずに確率分布関数（ｐｄｆ）を直接推定する。
２．データ点の個数がパラメータの個数を大幅に上回る漸近的極限において、ＭＬ法により推定パラメータの分散は、競合統計値の分散以下である。当該漸近的極限において、推定パラメータの分散は、その推定に用いたデータ点の個数に反比例する。所定の精度に対して、ＭＬ法は従って、他の方法よりも小さいサンプルからパラメータを推定することができるため、サンプリング効率がより高い。代替的な推定量の効率は実際、ＭＬ分散と他の推定量の分散との比として定義される。これにより、ＭＬ推定量の効率が定義から単位元に設定される一方、競合者の効率は単位元以下の分数に設定される。
３．推定パラメータの不確実性の共分散行列は、情報行列（フィッシャー、１９２２年）からの漸近的極限において容易に計算され、ＬＬＦの二次偏導関数のヘッセン行列は最小値をとる。

非漸近領域において、パラメータの個数がデータ点の個数と同等か上回る場合、モデルがランダムな統計的ノイズを再現可能な信号として扱うことを回避すべく解を制約する必要がある（例：Ｐｕｅｔｔｅｒ、Ｇｏｓｎｅｌｌ＆Ｙａｈｉｌ、２００５参照）。一般的な方式は、信号を基底関数の一般的な線形「ノンパラメトリック」結合として表してその係数を推定するものである。（基底関数を特徴付ける更なる非線形パラメータがあってもよい）。推定により、有意な係数値を得ると共に、有意でない係数をゼロにするかまたは少なくとも最小化することにより制約することが目標である。このようにして、信号をノイズから分離することが望まれる。

最も信頼性が高いパラメータ化は最も保守的なものであり、最小限の複雑さまたはオッカムの剃刀としても知られる、入力データと整合する最も単純な基本的パラメータを求めるものである。単純性は文脈依存であるが、大多数の連続的なアプリケーションの場合、最も単純な解は最も滑らかなものである。ＰＩＸＯＮ（登録商標）法は、データにより許容される最大且つ空間的に適合可能な平滑化を課すことによりこの解を与える（本明細書に引用するＰｉｎｙａ＆Ｐｕｅｔｔｅｒ、１９９３；Ｐｕｅｔｔｅｒ＆Ｙａｈｉｌ １９９９；Ｐｕｅｔｔｅｒ ｅｔ ａｌ２００５；米国特許第５，９１２，９９３号、第６，３５３，６８８号、第６，４９０，３７４号、第６，８９５，１２５号、第６，９９３，２０４号、第８，０１４，５８０号、第８，０８６，０１１号、第８，０９０，１７９号、第８，１６０，３４０号、第８，３９６，３１３号；米国公開特許第２０１２／０２６３３９３号を参照）。ＡＬＧＥＢＲＯＮ（商標）法は、金融システムにおける予測およびリスク査定（米国特許第７，３２８，１８２号）等、連続的な空間に制約されない離散的な問題に対して設計された同等な技術である。

フィッシャー（１９１２、１９２２年）の先駆的な業績以来、統計的推定の一般的傾向はＭＬおよびそのＬＬＦ推定量を利用するものであった。しかし、ＭＬには多くの重大な短所があって、その有用性が制約される場合がある。
１．ＭＬは漸近的に効率的なだけである。非漸近領域の解に追加的な制約が適用された場合、ＭＬ法の効率はもはや保証されない。
２．パラメータの共分散行列は、非漸近領域の情報行列から推定できない。無制約ＭＬ法では実際、信号として扱われるノイズを増幅することにより顕著なアーチファクトが生じることが多い。解を制約することで上述のようなアーチファクトを減らすことができるが、残余の精度は従って情報行列ではなく制約条件により決定される度合が高い。
３．ＬＬＦは一般に、パラメータθの二次関数ではなく、パラメータの決定に計算の労力を掛けることは、特に大規模な問題の場合に、漸近的サンプリング効率が多少上がるにせよ引き合わない。
４．パラメータに関するＬＬＦの勾配は、ｐｄｆ ｆ（ｘ，θ）に反比例する項を有している。
低ｐｄｆ領域のデータは、恐らく異常値（赤色データ）を含み、従ってパラメータの推定に大きなバイアスおよび／または変動が生じる恐れがある。

上述の制約に鑑みて、ＭＬ法は、非漸近領域のノンパラメトリック推定の他の推定量に比べて特段の利点がない。精度および計算効率の方が、もはや存在しないＭＬの「利点」よりも重要である。

本発明によれば、画像または信号の再構築、データ解析等において、特にパラメータの個数がデータ点の個数と同等か上回る場合に、ノイズまたは他の不確実性を有する入力データをモデリングする方法を提供する。特に、二次尤度汎関数（ＱＬＦ）は、統計的推定問題における対数尤度関数を代替すべく用いられる。（ａ）ＱＬＦがフィッシャー一致を示し、従って統計的推定に適していること、および（ｂ）大多数の関心対象問題について漸近的極限における効率が対数尤度関数よりも１０〜１５％小さいだけであることが示されている。従って、対数尤度関数が、非漸近領域において本質的な利点を有しておらず、当該分野でＱＬＦの方が精度および計算効率並びに低確率事象に対する低感度により理想的な推定量であることが指摘される。以下の詳細な記述により、データ空間における関数の効率的な積分を含む、ＱＬＦを構築および最適化する明示的なプロシージャを提供する。この新技術の潜在的用途の例を示す。

ＬＬＦと同様に、ＱＬＦはデータのビン化を必要とせず、二次形式であることで解の計算および制約の適用が簡素化され、低ｐｄｆ領域においてデータに対する感度が低い。漸近的極限においても、一般に用いられる基底関数に対するサンプリング効率においてＱＬＦはＬＬＦに大きく見劣りしない。線形非パラメータｐｄｆを考慮する場合、ＱＬＦの利点が明らかになる。

本開示は統計的推定に焦点を合わせているが、本明細書に示すいくつかの技術が一般的な最適化問題に適用できることは当業者には自明であろう。

本発明の一態様において、対象物体を記述するデータおよび複数のパラメータを含んでいると共にノイズ部を有する入力信号をソースから受信し、初期パラメータのグループを選択し、二乗可積分な基底関数の組の線形結合を含むノンパラメトリック確率分布関数（ｐｄｆ）を推定し、二次尤度汎関数（ＱＬＦ）を計算し、ＱＬＦに基づいてデータに対する初期パラメータの適合度を評価し、所定の条件が満たされるまで新たなパラメータのグループを選択して当該新たなパラメータのグループの適合度を評価することによりＱＬＦを反復的に最適化することにより、コンピュータプロセッサ内に対象物体のモデルを構築する方法およびシステムを提供する。受容可能な適合が達成されたならば、最適化されたパラメータを用いて構築された対象物体のモデルの出力を表示することができる。

本発明の別の態様において、データ成分およびノイズ成分を有する入力信号から対象物体の再構築された画像を生成する方法は、コンピュータに、複数のパラメータを含む当該入力信号を画像ソースから受信し、複数のパラメータから初期パラメータのグループを選択し、二乗可積分な基底関数の組の線形結合を含むノンパラメトリック確率分布関数（ｐｄｆ）を推定し、
の形式（θはパラメータを表し、ｘは観測値を表し、ｆ（ｘ、θ）はｐｄｆ）で二次尤度汎関数（ＱＬＦ）を計算し、当該ＱＬＦに基づいて当該データに対する当該初期パラメータの適合度を評価し、所定の条件が満たされるまで新たなパラメータのグループを選択して当該新たなパラメータのグループの適合度を評価することによりＱＬＦを反復的に最適化して、当該最適化されたパラメータに基づいて当該対象物体の再構築された画像の表示を含む出力を生成する命令を実行させるステップを含んでいる。いくつかの実施形態において、入力信号は第１の平面画像データおよび第２の平面画像データを含み、ＱＬＦを選択、推定、計算、評価して反復的に最適化する各ステップは第１の平面画像データおよび第２の平面画像データの各々について実行され、出力を生成するステップは対象物体の３次元画像を表示するステップを含んでいる。

本発明の更に別の態様において 入力信号に含まれる対象物体を記述するデータをモデリングするシステムは、コンピュータ可読媒体と、当該コンピュータ可読媒体に接続されたパラメータ最適化プロセッサと、当該パラメータ最適化プロセッサに接続されていて前記パラメータ最適化プロセッサとの間で再構築されたモデルの電子表現を受送信すべく適合された通信インターフェースとを含み、当該コンピュータ可読媒体が、前記パラメータ最適化プロセッサにより実行されたならば前記パラメータ最適化プロセッサに以下の動作、すなわち被写体データを収集すべく構成されたソースから入力信号を受信し、当該被写体データに対応する初期パラメータのグループを生成し、二乗可積分な基底関数の組の線形結合を含むノンパラメトリック確率分布関数を推定し、
の形式（θはパラメータを表し、ｘは観測値を表し、ｆ（ｘ、θ）はｐｄｆ）で二次尤度汎関数（ＱＬＦ）を計算し、当該ＱＬＦに基づいて当該データに対する当該初期パラメータの適合度を評価し、所定の条件が満たされるまで新たなパラメータのグループを選択して当該新たなパラメータのグループの適合度を評価することによりＱＬＦを反復的に最適化して、最適化されたパラメータを用いて構築された当該対象物体のモデルを含む出力を生成するステップを含む動作を実行させるソフトウェア命令を保存している。いくつかの実施形態において、当該データは重みｗ_ｉを含み、ＱＬＦは、
の形式を有している。
本システムは、データおよび基底関数を用いてソース項を計算するソフトウェア命令を更に含んでいてよい。ＱＬＦは、基底関数を用いるグラム行列を計算して、グラム行列、パラメータおよびソース項を組み合せてＱＬＦを求めることにより得られる。いくつかの実施形態において、入力信号は画像データであり、出力は、グラフィカルユーザーインターフェースに表示された対象物体の２次元、３次元または４次元表現を含んでいる。画像データは、Ｘ線、ＣＴ、放射断層撮影、ＳＰＥＣＴおよびＰＥＴからなるグループから選択されてよく、対象物体は患者の身体部位である。画像データは、出力が３次元表現を含むように少なくとも２平面で取得されてよい。いくつかの実施形態において、画像データは、出力が４次元表現を含むように少なくとも２平面で取得され、更に時間を含んでいる。

本発明の更に別の態様において、統計的推定を用いて入力データからモデルを生成する方法を提供し、改良点として、モデルを生成するためのパラメータを最適化すべく対数尤度関数（ＬＬＦ）を二次尤度汎関数（ＱＬＦ）で代替することが含まれる。

図１は、本発明を実施可能な汎用計算環境のブロック図である。 図２は、従来の反復的パラメータ推定処理の模式的フロー図である。 図３は、本発明による反復的パラメータ推定処理の模式的フロー図である。 図４ａは、黒い点で印を付けた７個の部位を有する円形領域の通常のボロノイ分割および小円で印を付けたボロノイセルの重心と、７個の部位を有する同じ領域のＣＶＴとを示す。図４ｂは、黒い点で印を付けた７個の部位を有する円形領域の通常のボロノイ分割および小円で印を付けたボロノイセルの重心と、７個の部位を有する同じ領域のＣＶＴとを示す。 図５は、本発明の少なくとも１個の実施形態による例示的な統計的推定システムのブロック図である。 図６は、強い周期的成分、若干の白色ノイズ、および定常的な背景からなるシミュレーションされた信号を示す。 図７は、単一のピークおよび若干の白色ノイズを示す図６の信号のパワースペクトル（ゼロ周波数は省略）を示す。 図８は、図６のソースからの「遠距離で検出された」２００個の光子のヒストグラムである。 図９は、図８の２００個の光子のヒストグラムの周波数（ゼロ周波数は省略）を有し、周期性を示さないパワースペクトルのプロット図である。 図１０は、正しい周波数で明らかな周期性を示す、検出された２００個の光子のＱＬＦパラメータ化により決定された周波数（ゼロ周波数は省略）を有するパワースペクトルのプロット図である。 図１１は、実時間定位生検の実行を支援すべくＸ線データをモデリングする本発明のＱＬＦ法の例示的アプリケーションを示すフロー図である。 図１２は、血管造影プロシージャの実行中に実時間撮像用のＸ線データをモデリングする本発明のＱＬＦ法の例示的アプリケーションを示すフロー図である。 図１３は、関心対象領域の３次元または４次元画像を生成すべくＣＴまたは放射断層撮影データをモデリングする本発明のＱＬＦ法の例示的アプリケーションを示すフロー図である。

本発明の態様について述べる前に、本発明を実装できる適当なコンピュータシステム環境１００（図１）を簡単に説明することは有用であろう。コンピュータシステム環境１００は、適当な計算環境の一例に過ぎず、本発明の利用範囲または機能に関して何らの限定を示唆するものではない。また、計算環境１００が、例示的な動作環境１００に示す構成要素の任意の１個または組み合わせに関して何らかの依存性または必要性を有しているものと解釈すべきでない。

本発明は、多数の他の汎用または専用コンピュータシステム環境または構成で動作する。本発明での利用に適した公知のコンピュータシステム、環境、および／または構成の例として、パーソナルコンピュータ、サーバコンピュータ、携帯またはラップトップ機器、マルチプロセッサシステム、マイクロプロセッサ利用システム、セットトップボックス、プログラマブル家電、ネットワークＰＣ、ミニコンピュータ、メインフレームコンピュータ、電話システム、上述のシステムまたは装置の任意のものを含む分散計算環境等が含まれるがこれらに限定されない。

本発明は、コンピュータにより実行されるプログラムモジュール等、コンピュータ実行可能な命令の一般的な文脈で記述できる。一般に、プログラムモジュールには、特定のタスクを実行する、または特定の抽象的なデータ型を実装するルーティン、プログラム、オブジェクト、コンポーネント、データ構造等が含まれる。当業者は、本明細書の記述および／または図面を以下に述べる任意の形式の計算機可読媒体に保存できるコンピュータ実行可能な命令として実装することができる。

本発明はまた、通信ネットワークを介して接続されたリモート処理装置によりタスクが実行される分散型計算環境で実施されてもよい。分散計算環境において、プログラムモジュールは、メモリストレージ装置を含むローカルおよびリモートコンピュータ記憶媒体の両方に搭載することができる。

図１を参照するに、本発明を実装する例示的なシステムは、コンピュータ１１０の形式で汎用計算装置を含んでいる。コンピュータ１１０の構成要素は、処理部１２０、システムメモリ１３０、およびシステムメモリを含む各種のシステム要素を処理部１２０に接続するシステムバス１２１を含んでいてよいがこれらに限定されない。システムバス１２１は、メモリバスまたはメモリコントローラ、周辺バス、および各種バスアーキテクチャの任意のものを用いるローカルバスを含む数種類のバス構造の任意のものであってよい。非限定的な例として、このようなアーキテクチャは、業界標準アーキテクチャ（ＩＳＡ）バス、マイクロチャネルアーキテクチャ（ＭＣＡ）バス、拡張ＩＳＡ（ＥＩＳＡ）バス、ビデオ電子標準協会（ＶＥＳＡ）ローカルバス、およびメザニンバスとしても知られる周辺装置相互接続（ＰＣＩ）バスを含んでいる。

コンピュータ１１０は典型的に、各種の計算機可読媒体を含んでいる。計算機可読媒体は、コンピュータ１１０からアクセス可能な任意の入手可能な媒体であってよく、揮発性および不揮発性媒体、着脱可能および着脱不可能媒体のいずれも含まれる。非限定的な例として、計算機可読媒体は、コンピュータ記憶媒体および通信媒体を含んでいてよい。コンピュータ記憶媒体は、計算機可読命令、データ構造、プログラムモジュールその他のデータ等の情報を保存する任意の方法または技術で実装された揮発性および不揮発性の着脱可能および着脱不可能媒体を含んでいる。コンピュータ記憶媒体は、ＲＡＭ、ＲＯＭ、ＥＥＰＲＯＭ、フラッシュメモリその他メモリ技術、ＣＤ−ＲＯＭ、デジタル多用途ディスク（ＤＶＤ）その他の光ディスクストレージ、磁気カセット、磁気テープ、磁気ディスクストレージその他の磁気ストレージ装置、または所望の情報を保存すべく利用可能であってコンピュータ１１０からアクセス可能な他の任意の媒体を含んでいるがこれらに限定されない。通信媒体は典型的に、計算機可読命令、データ構造、プログラムモジュールその他のデータを、搬送波その他の搬送機構等の変調されたデータ信号で表現し、且つ任意の情報配信媒体を含んでいる。「変調されたデータ信号」という用語は、信号内の情報を符号化すべく自身の１個以上の特徴が設定または変更された信号を意味する。非限定的な例として、通信媒体は、有線ネットワークまたは直接的有線接続等の有線媒体、および音響、ＲＦ、赤外線その他の無線媒体等の無線媒体を含んでいる。上述の任意の要素の組み合わせもまた、計算機可読媒体の範囲内に含めるべきである。

システムメモリ１３０は、読出し専用メモリ（ＲＯＭ）１３１およびランダムアクセスメモリ（ＲＡＭ）１３２等の揮発性および／または不揮発性メモリの形式でコンピュータ記憶媒体を含んでいる。起動中等にコンピュータ１１０内の要素間での情報転送を支援する基本ルーティンを含む基本入出力システム１３３（ＢＩＯＳ）は典型的にはＲＯＭ１３１に保存されている。ＲＡＭ１３２は典型的に、処理部１２０から直ちにアクセス可能な、および／または現在処理部１２０上で動作中であるデータおよび／またはプログラムモジュールを含んでいる。非限定的な例として、図１に、オペレーティングシステム１３４、アプリケーションプログラム１３５、他のプログラムモジュール１３６、およびプログラムデータ１３７を示す。

コンピュータ１１０はまた、他の着脱可能／着脱不可能な揮発性／不揮発性コンピュータ記憶媒体（「非一時的コンピュータ可読媒体」）を含んでいてよい。一例に過ぎないが、図１に、着脱不可能な不揮発性磁気媒体との間で読み出しまたは書き込みを行うハードディスク装置１４１、着脱可能な不揮発性磁気ディスク１５２との間で読み出しまたは書き込みを行う磁気ディスクドライブ１５１、およびＣＤＲＯＭその他の光媒体等の着脱可能な不揮発性光ディスク１５６との間で読み出しまたは書き込みを行う光ディスクドライブ１５５を示している。例示的な動作環境で利用可能な他の着脱可能／着脱不可能な揮発性／不揮発性コンピュータ記憶媒体は、磁気テープカセット、フラッシュメモリカード、デジタル多用途ディスク、デジタルビデオ・テープ、ソリッドステートＲＡＭ、ソリッドステートＲＯＭ等を含むがこれらに限定されない。ハードディスク装置１４１は典型的にインターフェース１４０等の着脱不可能なメモリインターフェースを介してシステムバス１２１に接続され、磁気ディスクドライブ１５１および光ディスクドライブ１５５は典型的にインターフェース１５０等の着脱可能なメモリインターフェースを介してシステムバス１２１に接続されている。

上述の、および図１に示すドライブおよび付随するコンピュータ記憶媒体は、計算機可読命令、データ構造、プログラムモジュールおよびコンピュータ１１０用の他のデータを保存する。例えば図１において、ハードディスク装置１４１を、オペレーティングシステム１４４、アプリケーションプログラム１４５、他のプログラムモジュール１４６、およびプログラムデータ１４７を保存しているものとして示す。これらの構成要素は、オペレーティングシステム１３４、アプリケーションプログラム１３５、他のプログラムモジュール１３６、およびプログラムデータ１３７と同一であっても、または異なっていてもよいことに注意されたい。オペレーティングシステム１４４、アプリケーションプログラム１４５、他のプログラムモジュール１４６、およびプログラムデータ１４７は、少なくとも異なる複製であることを示すために本明細書では異なる番号を付与されている。

ユーザーは、キーボード１６２等の入力装置、マイクロホン１６３（電話を介して与えられた入力を表す場合もある）、およびマウス、トラックボールまたはタッチパッド等のポインティング装置１６１を介してコンピュータ１１０に命令および情報を入力することができる。他の入力装置（図示せず）としてジョイスティック、ゲームパッド、衛星パラボラアンテナ、スキャナ等を含んでいてよい。これらおよび他の入力装置は多くの場合、システムバスに接続されたユーザー入力インターフェース１６０を介して処理部１２０に接続されているが、パラレルポート、ゲームポートまたは汎用シリアルバス（ＵＳＢ）等の他のインターフェースおよびバス構造により接続されていてもよい。モニタ１９１その他の種類の表示装置もまた、ビデオインターフェース１９０等のインターフェースを介してシステムバス１２１に接続されている。モニタに加え、コンピュータはまた、出力周辺装置インターフェース１９５を介して接続可能なスピーカー１９７およびプリンタ１９６等の他の周辺出力装置を含んでいてよい。

コンピュータ１１０は、リモートコンピュータ１８０等の１個以上のリモートコンピュータへの論理接続を用いてネットワーク化された環境で動作可能である。リモートコンピュータ１８０は、パーソナルコンピュータ、携帯機器、サーバ、ルータ、ネットワークＰＣ、ピア装置その他の一般的なネットワークノードであってよく、典型的にはコンピュータ１１０との関連で上に述べた要素の多くまたは全てを含んでいる。図１に示す論理接続は、ローカルエリアネットワーク（ＬＡＮ）１７１および広域ネットワーク（ＷＡＮ）１７３を含んでいるが、他のネットワークをも含んでいてよい。そのようなネットワーク環境は、オフィス、企業内コンピュータネットワーク、イントラネットおよびインターネットで普及している。

ＬＡＮネットワーク環境で用いる場合、コンピュータ１１０はネットワークインターフェースまたはアダプタ１７０を介してＬＡＮ１７１に接続される。ＷＡＮネットワーク環境で用いる場合、コンピュータ１１０は典型的に、インターネット等のＷＡＮ１７３上で通信を確立するモデム１７２その他の手段を含んでいる。モデム１７２は、内蔵または外付けのいずれにせよ、ユーザー入力インターフェース１６０その他の適当な機構を介してシステムバス１２１に接続されていてよい。ネットワーク化された環境において、コンピュータ１１０またはその一部との関連で示すプログラムモジュールは、リモートメモリストレージ装置に保存されていてよい。非限定的な例として、図１に、リモートアプリケーションプログラム１８５をリモートコンピュータ１８０に常駐するものとして示している。図に示すネットワーク接続が例示的であって、コンピュータ間で通信リンクを確立する他の手段が用いてもよいことが理解されよう。

図５は、少なくとも１個の実施形態による、入力画像、信号その他のデータの正確なモデルを構築すべく利用可能な例示的な統計的推定システム５００の簡略ブロック図である。システム５００は、プロセッサ５０２、オペレーティングシステム５０４、メモリ５０６、および入出力インターフェース５０８を含んでいる。メモリ５０６はＱＬＦエンハンサ５１０、およびパラメータデータ５１２およびＱＬＦエンハンサ５１０を用いるパラメータの最適化に用いるデータを保存する区画を含んでいてよい。

動作時においてプロセッサ５０２は、メモリ５０６に保存されたパラメータ推定部５１０を実行することができる。ＱＬＦエンハンサ５１０は、プロセッサにより実行されたならば本開示に従いパラメータを推定する動作をプロセッサに実行させるソフトウェア命令を含んでいてよい。ＱＬＦエンハンサ５１０は、オペレーティングシステム５０４と協働的に動作可能であって、本明細書に記述するようにパラメータデータ５１２を利用して、入出力インターフェース５０８で表示すべく得られたモデルの出力、すなわち再構築されたデータを生成することができる。再構築されたデータの例として画像、プロット、グラフ等が含まれるがこれらに限定されない。具体的な説明例を以下の第１２節に示す。

以下の詳細な記述は、二次尤度汎関数が統計的に整合性を有しているため統計的推定に統計的推定に適していることを最初に示すために編集されている。対数尤度関数が非漸近領域で本質的に利点を有しておらず、二次尤度汎関数が精度および計算効率の面で理想的な推定ツールとなることを示す。データ空間における関数の効率的な積分を含む、二次尤度汎関数を構築および最適化する明示的なプロシージャを提供する。

以下の記述では、読者の便宜を図って表１に掲載する多くの略語を用いる。

本明細書で用いる「信号」とは、コンピュータプロセッサにアップロードされて解釈されデジタルフォーマットで構成されたデータを意味する。信号の例として、デジタル画像データ、通信伝送、工業測定、市場調査データ、サプライチェーン管理データ、在庫管理、および金融市場等が含まれるがこれらに限定されず、いずれも信号内にノイズ成分が含まれており、ノイズの除去はデータの解析結果を最適化するために望ましい。

本明細書で用いる「ノイズ」とは、測定／受信設備の限界、対象物体と信号検出器／受信器との干渉（例：大気、電磁気、機械的、その他）、および非強調信号に混入または導入されて、入力に含まれるデータを正確に反映する所望の品質の出力表現を生成する能力を損なう恐れがある他のあらゆる不確実性の結果、信号に含まれる望ましくない成分を意味する。

本明細書で用いる「対象」または代替的に「対象物体」とは、データにより記述される項目を指す。データが画像データである場合、対象は、関心対象の特徴、例えば身体の内臓、地質学的構造、天体等、検出された信号を用いて情報が得られる任意の構造である。データが市場調査データである場合、対象は、当該データに関する、または当該データが収集された母集団であってよい。

本発明の技術の詳細な記述を明快にすべく章に分割する。以下に各章の目次を示す。
１．二次尤度関数（ＱＬＦ）
２．フィッシャー一致の実証
３．漸近的効率
４．線形ノンパラメトリック確率密度関数
５．重み付き二次尤度汎関数
６．不均一性補正
７．バックグラウンド
８．ＱＬＦの最適化
９．前処理部
１０．非線形ｐｄｆ
１１．数値積分
１２．例
ｉ．疎な信号検知
ｉｉ．デジタルマンモグラフィ
ｉｉｉ．血管造影
ｉｖ．計算断層撮影
ｖ．放射断層撮影

１．二次尤度汎関数（ＱＬＦ）
ＬＬＦには漸近的効率の利点があるが、一般に最適化が困難である。対照的に、二次尤度汎関数（ＱＬＦ）は、常に二次であるがＬＬＦではない。ＱＬＦの利点は、漸近的効率が若干低下する代償として計算が便利な点である。ＬＬＦの場合と同様に、ＱＬＦもデータのビン化を必要とせず、次式で与えられる。
ここに、θはデータのモデル構築に用いるパラメータを表し、ｘは、信号を表す。

ノイズがガウスｐｄｆに従うデータのＱＬＦとＬＬＦとの差異を強調することが重要である。ＱＬＦはｐｄｆにおいて二次形式のであるのに対し、ＬＬＦは対数形式である。また、ガウスｐｄｆの対数の平均値が二次形式となり得る。（これらは、ベル状のガウスｐｄｆの標準偏差、すなわち幅の二次形式ではない）。関心対象パラメータがガウス平均である場合、ガウスＬＬＦはこれらのパラメータの二次関数であり、ＭＬ法は最小二乗法に帰着する。対照的に、ガウスｐｄｆのＱＬＦは平均の指数関数である。

一方、ｐｄｆがパラメータの一次形式である、すなわち振幅が求めるパラメータである（３章および４章を参照）基底関数の線形結合である場合、ｐｄｆが二次形式であるＱＬＦはまた、パラメータの二次形式でもある。一方、線形ｐｄｆのＬＬＦは振幅の対数関数である。

従って、線形ｐｄｆのガウスＬＬＦおよびＱＬＦは共に各々のパラメータの二次関数に帰着し、類似の技術を適用してパラメータを求めることができる。しかし、基本的なノイズパターンは全く異なっていて互いに排他的である。ガウスｐｄｆを有する問題のパラメータ推定は、ＬＬＦを用いて求めることができ、且つそのように求めるべきであるのに対し、線形ｐｄｆはＱＬＦを用いる方が扱いやすい。

ガウスＬＬＦと線形ＱＬＦの差異を別の仕方で述べると、ＬＬＦは何らかの物理量を計量する際に用いられるのに対し、ＱＬＦは事象の発生率を測定するために用いられる。無論、事象をビン単位でカウントし、次いでポアソンＬＬＦ推定量、すなわち式（１）を用いてもよい。実際、カウント値が高ければポアソンｐｄｆはガウス分布を近似し、従ってカウントを測定したい物理量と見なすことができるが、このアプローチは低いカウント値において顕著に問題となる。このときＱＬＦが非ビン化データを用いる代替的な統計的推定法を提供する。

２．フィッシャー一致
ＬＬＦおよびＱＬＦは共にｐｄｆの同一の完全無制約な推定を与える。式（１）および（４）をｆに関して変化させることで次式が各々得られる。

確率変数δｆは任意であるため、角括弧内の２項は等しくゼロでなければならず、両式は同一の解を与える。

完全無制約ｐｄｆは、観測毎に変化する全てのランダムな統計的変動を含むため実用面で殆ど関心対象外である。従って次に区分的定数ｐｄｆを考察する。同様の導出により、両方の推定量を最小化することで再び同一の推定ｐｄｆが得られ、ビン内の周波数密度は定数ｐｄｆの領域に対応することを示す。
ここに、Ｎ_ＫはビンＫ内の観測値の個数（カウント値）、Ｖ_Ｋはその量である。

上述の考察はまた、ＬＬＦおよびＱＬＦがフィッシャー（１９２２）一致、すなわち、サンプルが無限に増加するに従い正しいｐｄｆを与えることを示す。これは、両方の推定量が正当であることを意味するため、重要な結論である。両者の間の選択は、各々の利点と短所を加味するユーザーに委ねられる。

３．漸近的効率
ＭＬの大きな利点はその漸近的効率にある。本章は、ＱＬＦの効率がどれくらいＭＬを下回るのかを評価する。

本発明は、二乗可積分な基底関数の組の線形結合として表現可能なノンパラメトリックｐｄｆに焦点を合わせる。

ここに、Ｂ_α（ｘ）は基底関数であり、係数θ^＊ _αはｐｄｆを特徴付ける未知パラメータθ_αの複素共役であり、上式の第２項は簡潔なベクトル行列表記を用いる。

ｐｄｆは、実数且つ非負でなければならないが、ｐｄｆが実数且つ非負のままであるように制約されている限り、式（９）は複素係数を有する複素基底関数、例えばフーリエ変換を許す。実基底関数の場合、パラメータは実数である。

最も簡単な場合は、単一の実非負基底関数であり、その唯一の未知数が非負パラメータである。
ｆ（ｘ）＝θＢ（ｘ）、Ｂ（ｘ）≧０（単一実基底関数）（１１）

１個の実基底関数しか存在しないため複素数およびベクトル表記は省略しており、以下で具体的には基底関数の正規化は単位元に設定される。

ＬＬＦ推定は、
であって、Ｎは観測されたカウント値である一方、ＱＬＦ推定は次式で与えられる。

最初に、ＬＬＦおよびＱＬＦ推定は共に偏りが無い点に注意されたい。所定個数の観測値について基底関数Ｂ（ｘ）にわたり式（１３）および（１４）の平均を求めて、式（１２）の正規化を用いれば、いずれの場合も次式が得られる。

観測値の個数が、平均λを有するランダムなポアソン変数であれば期待値に偏りが無い。

次に分散を見るに、ＬＬＦ推定の分散は単にポアソン分散である。

ＱＬＦの場合、式（１４）の和の項を二乗し、Ｂ（ｘ）にわたり平均を求め、次いでＮにわたり平均を求める必要がある。結果的に得られる分散は次式で与えられる。
ここで式（１８）が式（１２）の単位積分だけでなく任意の正規化について正しくなるように一般的な正規化を再導入する。ＱＬＦ推定量の効率は従って次式で与えられる。

また式（１９）から、効率が並進およびスケーリングに対して不変である点にも注意されたい。

ドット積不等式（コーシー、１８２１年；シュヴァルツ、１８８８年）を適用すれば、式（１９）から、ＱＬＦ効率が必然的に単位元を上界とし、当該上界に達するのはＢ（ｘ）が自身の非ゼロ領域にわたり定数である場合に限られる。また、Ｂが三乗可積分でなければ効率がゼロまで減少し続けることが分かる。興味深い問題は、実際に関心対象の基底関数に対して効率が単位元からどこまで低下するかである。従って次に、一般的に用いられる基底関数を考える。（簡潔を旨として議論を一次元に限定するが、より高い次元に計算を拡張することは容易である）。

べき法則基底関数に対して次式が成り立つ。

式（１９）からＱＬＦ効率が得られる。
ｅ＝（３α＋１）（α+１）／（２α＋１）^２（べき法則ＱＬＦ効率）（２１）

予想されるように、効率は、定数基底関数α＝０に対して最大値に達し、α→−１／３に近づくにつれてＢ（ｘ）^３にわたる積分が発散する場合ゼロに低下する。より現実的な線形すなわちα＝１の場合、効率はｅ＝８／９≒０．８９である。任意の三角基底関数についても２個の線形部分に分割にすることで分かるように同じ効率が得られる。次式のように非ゼロの切片を有する一般的な線形基底関数の場合、
効率は、上述の限界の間、０．８９≒８／９＜ｅ＜１に存在することが示される。

関心対象の他の場合では次式のガウス基底関数
Ｂ（ｘ）∝ｅｘｐ（−ｘ^２）（ガウス基底関数）、（２３）
に対して
および、
次式の逆放物線
に対してｅ＝１４／１５≒０．９３が含まれる。

最後に、２個の確率ｐおよびｑ＝１−ｐを有する圏論的二項の場合を考える。式（１９）は次式となり、
ｅ＝（ｐ^２＋ｑ^２）^２／（ｐ^３＋ｑ^３）（二項ＱＬＦ効率）、（２５）
ｐ＝１、ｑ＝０、ｐ＝０、ｑ＝１、およびｐ＝０．５、ｑ＝０．５に対してｅ＝１であり、その中間で最小値ｅ＝８／９≒０．８９に低下する。

従って、実用的な基底関数に対して、ＱＬＦの漸近的効率の損失は１０〜１５％の範囲にあり、これは以下に詳述するようにＱＬＦの利点を前提として許容可能であると結論付けることができる。しかし、ＬＬＦおよびＱＬＦの効率はいずれも非漸近領域では知られていないことに注意されたい。

４．線形ノンパラメトリックｐｄｆ
ＱＬＦの利点は、複数の二乗可積分基底関数（式）（９）の組み合わせである線形ノンパラメトリックｐｄｆを考慮することにより直ちに分かる。
θ^＊ _αに関する各推定量の式（１）および（４）の勾配をゼロに設定することで
が得られる。

３章で議論した単一の基底関数の場合を除いてＬＬＦ勾配の式はパラメータに関して極めて非線形式に見えるが、ＱＬＦ勾配は線形である。

ＱＬＦおよびその勾配を簡潔なベクトル行列表記で表すことで次式が得られる。

ここに、Ｑはグラム行列
であり、ｂはソース項
である。

後述する５章のＱＬＦの最適化は、グラム行列Ｑが常に正半定値であり、基底関数が線形独立ならば正定値であるという事実により極めて容易になる。また、Ｑの要素が有限であるために必然的に基底関数が二乗可積分でなければならない点にも注意されたい。

式（２８）におけるＱＬＦの第２の利点は、Ｑが使用する基底関数系、すなわち式（３０）だけに依存するのに対し、全てのデータの依存性がソース項ｂ、すなわち式（３１）に現れる点である。

５．重み付き二次尤度汎関数
上述のＱＬＦの変型例は、データ点に重み付けするものであり、各事象ｉは位置ｘ_ｉおよび既知の重みｗ_ｉにより特徴付けられる。例えば、観測される事象は既知だが異なるエネルギーを有する光子の検出であってよく、関心対象の分布関数は堆積エネルギーの分布関数であって光子の個数の分布関数ではない。別の例において、事象は異なる価格での販売であってよく、所望の分布関数はドル建ての販売高の分布関数であって販売事象の個数の分布関数ではない。このような場合、ＱＬＦは次式のように変更できる。

式（３２）の完全無制約ｐｄｆおよび区分的定数ｐｄｆの解は、式（７）および（８）を反映し、
推定されたｐｄｆの整合性を同様に示す。

３章の基底関数の線形結合の解析は上と同様に続く。グラム行列Ｑは不変であって、ソース項は次式のようになる。

特に、Ｑはここでも、使用する基底関数系だけに依存するのに対し、重み付けを含む全てのデータの依存性はソース項ｂに現れる。

以下の議論において、ｂが式（３５）により与えられるものと仮定し、重み付けの無い特別なケースは次式の通りである。
ｗ_ｉ≡１（重み付けの無い特別なケース）（３６）

６．不均一性の補正
恐らくはピクセル毎に急速に変化する、または「無効ピクセル」を含んでいるために、検出器全体で検出確率が均一でない場合、軽微な混乱が生じる。

この問題に対処する標準的な方法は、強い一定のソースについてカウント数Ｕ（ｘ）を測定するものである。（慣習的に、Ｕ（ｘ）は通常、単位元を平均または中央値として正規化される。）本発明の状況に当該較正を含める最も簡単な方法は、基底関数を次式で与えられる新たな基底関数で代替するものである。

解析はこれまで通り続く。

７．背景
「背景」はｘに対して緩慢に変化するｐｄｆの部分であるのに対し、「信号」はより急速に変化し得る。背景は往々にして異なる処理に起因するため、推定されたモデルにおける信号と背景の差異を維持する意味がある。

いくつかのアプリケーションにおいて、背景は未知のパラメータ、例えば１個の定数およびｘに対して緩慢に変化する項を有する数個の基底関数により特徴付けることができる。背景基底関数は次いで推定処理に含められ、最終的には単に背景としてラベル付けされて個別に報告される。

他のアプリケーションでは、背景は既知または他のデータから推定可能であって、推定パラメータの一部ではない個別の既知項Ａ（ｘ）としてｐｄｆに追加すべきである。信号および背景用の組み合わせｐｄｆは従って次式で与えられる。

信号が基底関数すなわち式（９）の線形結合として記述される場合、ｐｄｆは次式となる。

定数項を除けば、推定量は
但し
となる。

固定された外部背景を導入した結果、単に解空間のソース項ｂから背景項ａを減算するだけで済む。但し、ａは全てのｘにわたる積分すなわち式（４１）で与えられるのに対し、ｂは観測値ｘ_ｉについての和すなわち式（３１）である点に注意されたい。

８．最適化
式（４）の一般的なＱＬＦの最適化は、ｐｄｆがどのようにパラメータ化されるかに依存する。主な関心対象である線形ノンパラメトリックｐｄｆの場合、ＱＬＦは凸二次形式すなわち（式）（２８）であって、Ｑおよびｂは式（３０）および（３１）により与えられる。解きたい勾配方程式である式（２９）は従ってθの一次方程式である。しかし、殆ど全ての場合において、図２の基本フロー図２００に示すように、一次方程式は反復的に解く必要がある。本フロー図は当分野で公知の従来方式の反復的処理に従うため、「従来技術」とラベル付けされている。

フロー図２００に簡単に示すように、従来方式反復的処理は、モデリングしたい初期パラメータ２０２およびデータ２０３を提供する入力２０１から始まる。初期パラメータ２０１を用いてデータモデル２０７を構築する。ステップ２０４でデータ２０３に対するデータモデル２０７の適合度をパラメータ推定部により計算し、次いでステップ２０５で適合度の良否を評価する。適合度が所定の範囲または閾値に収まっている場合、推定パラメータは表示（ステップ２０８）および／またはメモリに保存（ステップ２０９）されてよい。指定された適合度条件に達していない場合、所望の条件が満たされるまでステップ２０６でパラメータを更新し、更新されたパラメータに基づいて変更されたデータモデルを構築して、適合度を評価することにより本処理を反復する。

図２に示す基本タイプの反復的スキームは通常、パラメータの推定には効果的であるが、同時に以下を含む重大な問題点があり得る。
１．大規模な問題には数百万個のデータ点および数百万個のパラメータが含まれる場合がある。従ってＢ、Ｂ^†およびＱを数兆個の要素を有する完全な行列としてメモリに保存するのは現実的でない。疎な行列表現は大抵可能であるが、多くのシステムにおいて疎な行列でも大き過ぎる。そのような場合、Ｂ、Ｂ^†およびＱを演算子として扱う必要があり、各行列乗算が膨大なパラメータ化を用いて新たに計算される。
２．必然的に、Ｂ、Ｂ^†およびＱ、または任意の関連行列の膨大な再配置を必要とするいかなる方法も用いることができない。特に、逆グラム行列Ｑ^−１またはＱのコレスキー分解の計算は論外である（例：Ｐｒｅｓｓ他、２００７年）。
３．多くの問題においてＱは扱い難い、すなわち、最大と最小の固有値の比の絶対値が極めて大きい。これは、最適化においてソース項ｂの統計的ノイズが、制御されなければ大幅に強調され、解における顕著なアーチファクトが生じることを意味する。最適化方法は従って、良好な制約機構を有していなければならない。

上述の考察は、反復的な共役勾配（ＣＧ）法の利用を示唆する（Ｆｏｘ、Ｈｕｓｋｅｙ＆Ｗｉｌｋｉｎｓｏｎ、１９４８年；Ｈｅｓｔｅｎｅｓ＆Ｓｔｉｅｆｅｌ、１９５２年）。ＣＧには、行列を明示的に計算される必要が一切なく、Ｂ、Ｂ^†およびＱを演算子として適用できる利点がある。ＱＬＦ勾配すなわち式（２９）のＱθ項は２段階のステップ、すなわちパラメータ空間からデータ空間への前方射影、およびデータ空間からパラメータ空間への後方射影Ｂに分解される。

推定パラメータの個数がデータ点の個数と同等か上回る非漸近的問題では、ノイズが信号として認識されて解にアーチファクトが生じるのを回避すべく追加的な制約が必要である。制約方法（例：Ｐｕｅｔｔｅｒ他による考察、２００５年）を、他の統計推定量に適用する場合と同様にＱＬＦ推定量に適用できるが、追加的な統計推定量が必要な場合がある。

大多数のＰＩＸＯＮ（登録商標）法の適用では、例えば、パラメータ推定部とは異なる適合度（ＧＯＦ）統計量を用いて平滑化の程度を判定する。例外として、パラメータ推定部の勾配にＰＩＸＯＮ（登録商標）平滑化を適用するが追加的なＧＯＦ統計量を必要としない新たな技術（米国特許第８，６７５，９８３号）がある。

ＣＧ法もまた、凸集合への射影法（Ｂｉｅｍｏｎｄ、ＬａｇｅｎｄｉｊｋおよびＭｅｒｓｅｒｅａｕ、１９９０年）により制約を課すものである。例えば、パラメータは多くの場合負でなくなければならない。
θ≧０（非負条件）（４３）

上の条件は、各反復の後で単に負の成分をゼロに切り捨て、あたかも切り捨てが生じなかったように継続することで課せられる（Ｐｕｅｔｔｅｒ他、２００５年）。

９．前処理部
ＣＧ法は、前処理部を追加して、次式によりＱＬＦの勾配を代替することにより大幅に加速できる（例：ＧｏｌｕｂおよびＶａｎ Ｌｏａｎ、１９９６年：Ｓａａｄ、２００３年）。
Ｇ＝Ｐ^−１∇Ｌ＝Ｐ^−１（Ｑθ−ｂ）（前処理部を有する勾配）、（４４）
従って、一次方程式
Ｐθ＝ｃ（前処理方程式）（４５）
を容易に解くことができる。結果的に得られる方法は、前処理済共役勾配（ＰＣＧ）法として知られる。

知的な問題依存型前処理部は、ＰＣＧ法の収束速度を大幅に向上させることができる。フェイルセーフな一般的選択は、ＰをＱすなわち次式で要素が与えられる対角行列の対角部に設定する単純な再スケーリングである。

スケーリング前処理部が機能するのは半定値行列Ｑの対角支配に起因する。

前処理により収束が加速する程度は、Ｑの非対角項が式（４７）内の極限からの距離およびＱの条件数に依存する。

しかし、式（４６）は、各々の基底関数について個別に評価する必要があり、演算子としてＢ^＊ _α次いでＢ_αの適用を含むため、大規模な問題で実行的に計算するのは困難であろう。代替策（ＡｈｎおよびＦｅｓｓｌｅｒ、２００３年）は、式（４６）を次式で代替するものである。

計算上の利点は、式（４８）における和を単位パラメータの配列への単一の前方射影として計算ができる点である。

１０．非線形の場合
式（９）における解の線形ノンパラメトリック表現の簡便性は強調して余りある。にもかかわらず、基底関数を特徴付けるために追加的な非線形パラメータを導入することが有用な場合が多い。一例として幅が調節可能な基底関数が挙げられる。この場合非線形パラメータに関してＱＬＦの偏導関数を追加して、非線形ＣＧ最適化（例：Ｐｒｅｓｓ他、２００７年）を求めることが必要である。

１１．数値積分
上の１０章までにおいてグラム行列がシステムだけに依存し、且つ全てのデータ依存性がソース項ｂに現れると述べた。この差異の実用性に現れるのは、Ｑがデータ空間にわたる空間積分すなわち式（３０）であるのに対し、ｂはデータ点にわたる和すなわち式（３１）であることである。

固定されたメッシュにわたる空間積分は、空間解像度が各次元において相応の数のメッシュ点を必要とする（「次元の呪い」）ため、問題の次元が増えるにつれて計算量が膨大になる。この難点に対処する２通りの方法がある。一つの方法は、データ空間内で適応型の階層的直交格子を形成するものである。代替的な方法は、「自然な」無構造メッシュはデータ点自体の位置により与えられるものと認識することである。データ位置のボロノイ（１９０８年）分割によりデータ空間を分割することができ、その内部で（通常は有限な）データボリュームＶがセルに分割され、空間内の各点が自身に最も近いデータ位置周辺のセルに割り当てられる。
Ｖ_ｉ＝｛ｘ∈Ｖ｜｜ｘ−ｘ_ｉ｜＜｜ｘ−ｘ_ｋ｜、∀ｋ≠ｉ｝、ｉ＝１、．．．、Ｎ（ボロノイセル）（５０）

ボロノイ分割に関する議論についてはＰｒｅｓｓ他（２００７年）を参照されたい。ボロノイセルのボリュームＶへの制約は、無境界ボロノイセルを切り捨てることにより実現される。

データ空間が分割されたならば、データ点の位置でｐｄｆを計算するだけでよい。関数Ｆ（ｘ）にわたる空間積分は、各データ点における被積分関数の値とそれらのボロノイセルのボリュームとの積にわたる和で近似される。

ボロノイセルは確かに、データ点の位置の統計的変動に起因して変動するボリュームを有するが、積分された関数は多数のデータ点にわたり滑らかでなければならないため、ボロノイボリュームの変動により生じた余分の分散は、計算された積分に重大な影響を及ぼさない筈である。

別の可能な方法は、分割のノードが自身の重心でもある（Ｄｕ、ＦａｂｅｒおよびＧｕｎｚｂｕｒｇｅｒ、 １９９９年）重心ボロノイ分割（ＣＶＴ）を構築するものである。ｐｄｆは次いでＣＶＴのノードで計算されて、より滑らかなセルボリュームにわたり積分される。ソース項ｂに対して必要なデータ点の位置におけるｐｄｆ値は、内挿により計算することができる。ＣＶＴの決定には多少の計算労力が必要であるが、これを行うために並列化された確率的方法（Ｊｕ、ＤｕおよびＧｕｎｚｂｕｒｇｅｒ、２００２年）およびＣＶＴエネルギー関数のＣ^２平滑度（Ｌｉｕ他、２００９年）を利用するより新たなニュートン型最適化がある。

本明細書に示す全てのアプリケーションの中心に、反復的パラメータ推定におけるＱＬＦの計算があり、その基本ステップを図２に示す。４〜１１章で述べたように、当該計算の本質を図３のフロー図に示す。データと重み３０１、および基底関数３０２を用いてソース項３０４およびグラム行列３０５を計算し、次いでこれらをパラメータ３０３と組み合わせてＱＬＦ３０７を与える。グラム行列３０５は通常、特に大規模な問題においては明示的に計算されないが、基底関数３０２およびそれらの共役はステップ３０６で演算子として適用される。追加的なオプションは、１１章に詳述するようにデータのボロノイ分割３０９に基づいて、データ空間にわたる積分をボロノイ積分３０８で近似するものである。このオプションのステップを、当該フロー図の破線接続で示す。４〜１１章に詳述するいくつかの他のオプション的特徴は、明快さのため省略しているが、当業者には明らかな筈である。

図４ａ、４ｂは通常のボロノイ分割をＣＶＴと比較しており、図４ａは黒点で印を付けた７個の部位を有する円形領域の通常のボロノイ分割、および小さい（開いた）円で印を付けたボロノイセルの重心を示し、図４ｂは７個の部位を有する同じ領域のＣＶＴを示す。

１２．事例
以下に、統計的推定に本発明のアルゴリズムを適用した例を示す。

事例ｉ：疎な信号検出
強い信号であっても、過度に疎にサンプリングされたならば検出が困難な場合がある。例えば、軍事信号は、少数の光子しか検出されていなくても遠距離から観測することができる。図６に、強い周期的信号に加えて非ゼロの定数背景を上回る若干の白色ノイズを発するソースのシミュレーションを示す。図７に、異なる周波数（ゼロ周波数は省略）の相対的な強度を示す図７の信号のパワースペクトル（Ｐｒｅｓｓ他、２００７年）を示す。

検出器が「ソースから非常に遠い」ため、比較的少数の光子、例えば２００個の光子が検出される。図８の光子の到着時間のヒストグラムは若干の高周波挙動を示すが、図９のヒストグラムのパワースペクトルに周期性は検出されない。一方、ＱＬＦ推定は図１０に示すように、正しい周波数でパワースペクトルに周期的信号を明瞭に検出する。

事例ｉｉ：デジタルマンモグラフィ
デジタルマンモグラフィは基本的に従来のマンモグラフィと同じＸ線マンモグラフィシステムを用いるが、当該システムはフィルムカセットの代わりにデジタルレセプタおよびコンピュータを備えている。標準的マンモグラフィにおいて、画像はＸ線カセットを用いてフィルムに記録され、当該フィルムは放射線科医により光ボックスを用いて確認されてから施設のアーカイブのジャケットに保存される。デジタルマンモグラフィにおいて、胸部画像は特別な電子Ｘ線検出器を用いて撮像され、当該検出器により当該画像をデジタル画像に変換してコンピュータモニタで確認してデジタルアーカイブシステムに保存する。デジタルマンモグラフィにおいて、画像の拡大、向き、明度、およびコントラストを変更して放射線科医が特定の領域をより明確に視認できるよう支援することができる。

デジタル全視野マンモグラフィは乳癌のスクリーニングおよび診断を可能にするものであり、いずれは従来のマンモグラフィを代替するであろう。デジタルスポット視野マンモグラフィは、より高速且つより正確な定位生検を可能にするものであり、異なる二方向から胸部を撮像するため、より正確に生検の３次元誘導を行う。その結果、検査時間が短縮されると共に、患者が静止姿勢をとり続けなければならない時間が大幅に短縮されるため患者の快適性および利便性が大幅に向上する。デジタルスポット視野マンモグラフィにおいて画像はデジタル的に撮像されて直ちにシステムモニタに表示されるのに対し、従来の定位生検ではマンモグラムフィルムを露光、現像、次いで視認する必要があるため、胸部生検が完了までの時間が極めて長くなる。

デジタルマンモグラフィは、標準的フィルムマンモグラフィに比べて多くの利点をもたらす潜在的可能性がある。すなわち（ａ）密な胸部組織と疎な胸部組織とのコントラストが向上し、（ｂ）画像取得が高速化（１分未満）され、（ｃ）試験時間が短縮（フィルムベースのマンモグラフィの約半分）され、（ｄ）画像の保存が容易になり、（ｅ）医師が画像を操作することで乳癌の診断がより正確になり、（ｆ）マンモグラムを反復する必要なしにフィルムの露出不足／過剰を修正可能であり、（ｇ）他の医師による遠隔診断用に画像のデジタル送信が可能になる。

ＱＬＦは、ＰＩＸＯＮ（登録商標）法等のノイズ低減技術と協働して、収集したデジタルデータをよりうまく利用してマンモグラフィ画像の鮮明さを向上させることでこれらの目標の実現を支援することができる。これにより医師は、より少ないＸ線カウント数で同じ画質を得ることが可能になり、より少ない線量、より短い走査時間、またはこれら二つの組み合わせが可能になる。

具体的には、ＱＬＦは勾配が線形であるため、勾配が非線形式であるポアソンＬＬＦよりも極めて高速なノイズ低減処理を行うことができる。これは特に、医師が生検針を誘導できるように画像を実時間で見る必要がある定位生検法で重要である。図１１に、ステップ１１０２でＱＬＦを用いて単一平面マンモグラムを強調し（ステップ１１０４）、その結果をステップ１１０６で医師に提示する例示的なフロー図を示す。第２の平面マンモグラムの誘導（ステップ１１１２）を用いて任意選択の定位生検を実行し、ＱＬＦ（ステップ１１１４）を用いて当該マンモグラムを強調して、生検の実時間誘導を行うべくステップ１１１６で医師に提示する。重心ボロノイ分割は、特に低カウント領域で更なるノイズ低減をもたらす。

また、より良質な画像を用いて、潜在的に異常な領域へ医師を自動的に案内する機械学習ツールを開発することができる。このようなツールは、異常の領域を見逃がす（偽陰性）ことがないが、一旦証明されたならば、試験結果の読み込み効率を大幅に向上させ、読み込みの完了までに要する時間を短縮して、診断に対する医師の信頼度が増すことを保証すべく徹底的に試験されなければならない。

ＱＬＦは、より良質な画像の生成だけでなく、基底関数の振幅の形式で画像の信頼性が高い数値表現を与えることによっても上記処理を支援することができる。経験的に、スペクトル表現として知られるこれらの振幅は多くの場合、元のピクセル画像強度よりも有用な機械学習ツール用の入力データ（より優れた機械学習言語）であることが分かっている。

事例ｉｉｉ：血管造影法
血管造影法は、特に動脈、静脈および心室に注目して、身体の血管および器官の内部を視覚化すべく用いる医療画像技術である。これは従来、放射線不透過性造影剤を血管内に注入してＸ線検出器を用いて撮像することにより行われる。血管造影法の種類に応じて、心臓の左側および動脈系を観察する場合は大抵大腿動脈を、または心臓の右側および静脈系を観察する場合は頸静脈または大腿動脈を通して血管内に進入する。ガイドワイヤおよびカテーテルシステムを用いて、ある種の造影剤（Ｘ線を吸収することで現れる）を血液に添加してＸ線画像で視認可能にする。今日用いられている血管造影システムは一般に、定位生検と同様に、より良質な３次元位置調整を可能にすべく二つの異なる視点からの画像を与える。

撮像するＸ線画像は、静止画像または動画像のいずれでもよい。心臓以外の全ての構造について、画像は通常、デジタル減算血管造影法（ＤＳＡ）と呼ばれる技術を用いて撮像される。この場合画像は通常、毎秒２〜３フレームで撮像されるため、介在する放射線科医が１個または複数の血管を通る血液流を評価することができる。ＤＳＡは、造影剤で満たされた血管だけを観察できるように骨および他の器官を減ずる。減算技術は患者が静止したままでいる必要があり、拍動および動いている心臓に用いることができないため、心臓の画像はＤＳＡを使用せずに毎秒１５〜３０フレームで撮像される。これら両方の技術により、介在する放射線科医または心臓外科医が血流を阻害して痛みの原因になり得る血管内部の塞がりまたは狭まり（狭窄）を観察することができる。

従来の血管造影Ｘ線検出器は、蛍光透視法、すなわち蛍光媒体を用いてＸ線光子を吸収し、それらのエネルギーを画像増倍管により光学光子として再放射して電荷結合素子（ＣＣＤ）ビデオカメラにより検出する処理を利用することにより、画像を記録してモニタで再生することができる。当該装置は、各Ｘ線光子から大量且つ可変的な数の光学光子が生じるため、厳密には光子カウント装置ではない。また、膨大な放射線量が生じるため、患者に影響を及ぼすだけでなく、臨床医が診療を実施できる時間も制限される。（防護服を着用しているが、完全に保護される訳でないため、週次および通算診療時間を制限する必要がある）。この放射線量を減らすべく、光子カウント装置として直接Ｘ線を検出する固体Ｘ線検出器を構築する努力を現在継続中である。導入予定の光子カウント装置においてＱＬＦは特に、画質を向上させて放射線量を減らす有効なツールであろう。

再び、ＱＬＦはノイズ低減技術と協働して、将来的な光子カウント装置において、画質を向上させるとともに放射線量を減らす特に有効なツールとなろう。その処理は、一般に二つの撮像平面を有し、医師が血管を通して導線を誘導するのと同時に処理を実時間で実行する必要がある点を除いてマンモグラフィで概説したものと同様である。図１２を参照するに、ステップ１２０２Ａ、１２０２ＢにおいてＸ線データを２平面で収集する。各平面は、ステップ１２０４Ａ、１２０４ＢにおいてＱＬＦ強調を用いて別々に処理される。処理された画像はステップ１２０６Ａ、１２０６Ｂにおいて医師に提示され、血管造影法処置を実行している医師を正確に誘導する。臨床医は患者の近くにあって、防護服を着用しているにもかかわらず相当量の放射線に晒されるため、放射線量を低減する潜在能力は本出願において特に重要である。放射線量が大幅に低減することで、臨床医がより長い時間患者に接触できるようになるため、当該医療を実施する方法が変わる。

事例ｉｖ：コンピュータ断層撮影
コンピュータ断層撮影（ＣＴ）は、患者または試験対象に対して侵襲的外科治療を行うことなく当該患者または試験対象の内部構造を調べられることができる医療診断および測定方法並びに試験技術を提供する。この場合、調査対象の物体の多くの射影像が様々な角度から記録され、そこから当該物体の３次元画像を計算することができる。

観測された射影像（データ）を３次元画像に変換することにより断層撮影画像が生成される。例えば、Ｘ線ＣＴ撮影において、Ｘ線ビームが物体に向けられ、当該物体内で変動する構造によりビームは様々な程度に減衰する。物体の反対側で、減衰したビームが検出器により測定される。このような射影像が、物体周辺において様々な角度で生成される。次いでこれらの射影像を用いて、断層撮影画像再構築として知られる処理により３次元画像を計算する。

ＣＴ画像再構築は当初、および今日でも依然として多く、フィルタリングされた後方射影として知られる簡単な変換により実行されている。しかし、射影測定値はノイズが多く、且つ相対ノイズレベルは減衰量に依存する。骨および特に金属等の高密度材料を通る射影は、肉体、水、または他のより低密度の材料を通る射影よりも信号対ノイズ比が低い。検出される光子の個数の大幅且つ空間的な変動に対処することで、画質を向上させるべく統計画像再構築技術の開発に至った。

適切な統計モデリングによりノイズ画像を低減することができるため、患者へのＸ線放射線量を減少させることが可能になる。問題は、（恐らくは非線形な）物理モデルおよび統計的ノイズモデルに従い測定に最も良く適合する画像を見出すことである。再構築は反復的であって、多くのパスを必要とする。ＱＬＦは統計画像再構築、特に個々の光子を検出する現在開発中の装置を用いるもので特に重要な役割を果たし得る。最終的な形式として、データは検出された光子のリストに編成される。各光子は、自身のエネルギーと、検出された時間および位置と共に個々に記録される。これは個々の事象を扱うべく設計されているため、ＱＬＦが理想的に適している検出の微細な量子限界である。

図１３に、例えば心臓、脳、または内臓あるいは身体構造の非侵襲的体積または機能評価を支援すべく３次元または４次元モデルの生成に用いるコンピュータ断層撮影または放射断層撮影データを強調するステップを示すフロー図を示す。処理フローは一般に、図２に示すものに続く。データ（デジタル化されたＣＴまたはＥＴ画像データ１３０１）および初期画像にシステムおよび構築パラメータを加えた値１３０２を用いて初期ｐｄｆ１３０７を計算する。ｐｄｆ１３０７は、光子リストモードデータ１３０３（各光子のエネルギー、時間および位置）に作用してＱＬＦ１３０４を計算する。ＱＬＦ推定値のデータ１３０１、１３０３との適合度をステップ１３０５で評価して適合度の良さを求める。適合度が所定の範囲または閾値に収まる場合、推定パラメータを用いて、表示（ステップ１３０８）および／またはメモリに保存（ステップ１３０９）可能な多次元画像を生成する。指定された適合条件に達しない場合、所望の条件が満たされるまで、ステップ１３０６で画像パラメータを更新し、更新されたパラメータに基づいて、変更されたデータモデルのｐｄｆをステップ１３０７で再計算し、適合度を評価（１３０５）することにより処理を反復して、結果を表示（１３０８）および／またはメモリに保存（１３０９）することができる。

時間情報が利用可能であることで、第４次元として時間を含む４次元ＣＴ撮影を可能にする潜在的能力を有している。これによりＣＴは、画像が時間と共に変化する動的なプロセス、例えば拍動および動いている心臓を測定することができる。制約因子は、各々の短い時間間隔での光子の個数であり、ＱＬＦはリストモードデータの非ビン化再構築によりこれを支援することができる。

事例ｖ：放射断層撮影
放射断層撮影は、ＣＴで用いるものと同様の断層撮影技術を用いて身体の機能工程の３次元画像を生成する核医療画像技術である。違いは、ガンマ線または陽電子を放射する放射性同位元素（放射性核種と呼ばれる）が患者の血流に注入される点である。ガンマ線を放射する放射性核種は単一の光子を放射し、その撮像方法は単一光子放射コンピュータ断層撮影（ＳＰＥＣＴ、または時折ＳＰＥＴ）として知られる。対照的に、放射される陽電子は、体内の電子により消滅して反対方向に移動する２個の光子を形成し、これらは同時に検出される。当該撮像方法は、陽電子放射断層撮影（ＰＥＴ）として知られる。

殆どの時間、放射特性だけが関心対象であるマーカー放射性同位元素が特定のリガンドに取り付けられて放射性リガンドを生成し、これはある種の組織への化学結合特性により関心対象である。この合体により、リガンドと放射性同位元素の組み合わせ（放射性医薬品）を搬送して体内の関心対象部位と結合させることができ、次いでこれにより（同位元素の直接または間接的ガンマ線放射により）リガンド濃度を撮像することができる。

核スキャンは、ＣＴまたは磁気共鳴撮像（ＭＲＩ）走査と合わせて行われるケースが増えており、その組み合わせにより解剖学的および代謝的情報の両方（すなわち、何の構造か、および何を生化学的に行っているか）を与える。核撮像は解剖学的撮像と組み合せて最も有用であるため、最新の核スキャナを一体化された最高級の多検出器列ＣＴスキャナ、またはより最近ではＭＲＩと合わせて利用できる。２種類の走査の間に患者が姿勢を変えることなく２回の走査が同一セッション内で連続的に、または同時に実行できるため、２組の画像がより正確に登録され、核画像上の異常領域を、より完全にＣＴまたはＭＲＩ画像上の解剖学的知見と相関させることができる。これは、動いている器官または脳の外側でより一般的である、より高い解剖学的変動を有する構造の詳細な画像を表示するのに極めて有用である。より最近では、核撮像よりも高い解像度を有する解剖学的モダリティを用いて核画像再構築自体を誘導している。（従来は、ＣＴ画像の役割は、上述の解剖学的相関に加え、核減衰マップの提供に限定されていた）。

コンピュータ断層撮影（ＣＴ）の再構築に極めて類似した技術を３次元画像の生成に用いることは一般的だが、核撮像において収集されたデータセットの光子はＣＴよりもはるかに少ないため、再構築技術はより困難である。次いでＱＬＦ処理を用いて画質を向上させる、および／または患者に浴びせられる放射線量を減らすことができる。放射断層撮影装置の物理的構成はＣＴとは極めて異なるが、これは前方および後方射影演算子（８章参照）だけに影響を及ぼし、実際にＣＴおよび放射断層撮影と極めて良く類似している基本的再構築方法論には影響を及ぼさない。図１３のフロー図は従って、システム構成およびパラメータは異なるものの、放射断層撮影にも同様に適用可能である。

上述の詳細な説明および例に基づいて、データの統計モデリングを行う本発明のＱＬＦ法の他の応用例は、当業者には容易に明らかになろう。従って、本発明の範囲は、本明細書に記述する具体例だけに限定されるものではない。

以下の参考文献のリストは、統計画像再構築の背景および全般的な最先端技術の理解に役立つ。一覧に示す参考文献は、当業者が本発明を完全に理解すると共に実施するために必要な程度に、当分野の技術レベルを示すと共に、本開示の一部として本明細書に引用している。

参考文献
Ｒｅｆｅｒｅｎｃｅｓ
１． Ａｈｎ Ｓ ＆ Ｆｅｓｓｌｅｒ ＪＡ ２００３，“Ｇｌｏｂａｌｌｙ ｃｏｎｖｅｒｇｅｎｔ ｉｍａｇｅ ｒｅｃｏｎｓｔｒｕｃｔｉｏｎ ｆｏｒ ｅｍｉｓｓｉｏｎ ｔｏｍｏｇｒａｐｈｙ ｕｓｉｎｇ ｒｅｌａｘｅｄ ｏｒｄｅｒｅｄ ｓｕｂｓｅｔｓ ａｌｇｏｒｉｔｈｍｓ”，ＩＥＥＥ Ｔｒａｎｓ．Ｍｅｄ．Ｉｍａｇ．，２２，６１３．
２． Ｂｉｅｍｏｎｄ Ｊ，Ｌａｇｅｎｄｉｊｋ ＲＬ ＆ Ｍｅｒｓｅｒｅａｕ ＲＭ １９９０，“Ｉｔｅｒａｔｉｖｅ Ｍｅｔｈｏｄｓ ｆｏｒ Ｉｍａｇｅ Ｄｅｂｌｕｒｒｉｎｇ”，Ｐｒｏｃ．ＩＥＥＥ，７８，８５６．
３． Ｃａｕｃｈｙ Ａ １８２１，Ｏｅｕｖｒｅｓ ２，ＩＩＩ，３７３．
４． Ｄｉｒａｃ Ｐ １９５８，Ｐｒｉｎｃｉｐｌｅｓ ｏｆ Ｑｕａｎｔｕｍ Ｍｅｃｈａｎｉｃｓ，４ｔｈ ｅｄ．，（Ｏｘｆｏｒｄ： Ｃｌａｒｅｎｄｏｎ Ｐｒｅｓｓ），ｐ．５８．
５． Ｄｕ Ｑ，Ｆａｂｅｒ Ｖ ＆ Ｇｕｎｚｂｕｒｇｅｒ Ｍ １９９９，“Ｃｅｎｔｒｏｉｄａｌ Ｖｏｒｏｎｏｉ Ｔｅｓｓｅｌｌａｔｉｏｎｓ：Ａｐｐｌｉｃａｔｉｏｎｓ ａｎｄ Ａｌｇｏｒｉｔｈｍｓ”，ＳＩＡＭ Ｒｅｖ．，４１，６３７．
６． Ｆｉｓｈｅｒ ＲＡ １９１２，“Ｏｎ ａｎ Ａｂｓｏｌｕｔｅ Ｃｒｉｔｅｒｉｏｎ ｆｏｒ Ｆｉｔｔｉｎｇ Ｆｒｅｑｕｅｎｃｙ Ｃｕｒｖｅｓ”，Ｍｅｓｓｅｎｇｅｒ Ｍａｔｈ．，４１，１５５．
７． Ｆｉｓｈｅｒ ＲＡ １９２２，“Ｏｎ ｔｈｅ Ｍａｔｈｅｍａｔｉｃａｌ Ｆｏｕｎｄａｔｉｏｎｓ ｏｆ Ｔｈｅｏｒｅｔｉｃａｌ Ｓｔａｔｉｓｔｉｃｓ”，Ｐｈｉｌ． Ｔｒａｎｓ． Ｒｏｙ． Ｓｏｃ． Ａ，２２２，３０９．
８． Ｆｏｘ Ｌ，Ｈｕｓｋｅｙ ＨＤ ＆ Ｗｉｌｋｉｎｓｏｎ ＪＨ １９４８，“Ｎｏｔｅｓ ｏｎ ｔｈｅ Ｓｏｌｕｔｉｏｎ ｏｆ Ａｌｇｅｂｒａｉｃ Ｌｉｎｅａｒ Ｓｉｍｕｌｔａｎｅｏｕｓ Ｅｑｕａｔｉｏｎｓ”，Ｑｕａｒｔ．Ｊ．Ｍｅｃｈ．Ａｐｐｌ．Ｍａｔｈ．，１，１４９．
９． Ｇａｕｓｓ ＣＦ １８０９，Ｔｈｅｏｒｉａ Ｍｏｔｕｓ Ｃｏｒｐｏｒｕｍ Ｃｏｅｌｅｓｔｉｕｍ ｉｎ ｓｅｃｔｉｏｎｉｂｕｓ ｃｏｎｉｃｉｓ ｓｏｌｅｍ ａｍｂｉｅｎｔｉｕｍ，Ｅｎｇｌｉｓｈ ｔｒａｎｓｌａｔｉｏｎ，１９６３，Ｔｈｅｏｒｙ ｏｆ Ｍｏｔｉｏｎ ｏｆ ｔｈｅ Ｈｅａｖｅｎｌｙ Ｂｏｄｉｅｓ（Ｎｅｗ Ｙｏｒｋ：Ｄｏｖｅｒ）．
１０． Ｇｏｌｕｂ ＧＨ ＆ Ｖａｎ Ｌｏａｎ ＣＦ １９９６，Ｍａｔｒｉｘ Ｃｏｍｐｕｔａｔｉｏｎｓ，３ｒｄ ｅｄ．（Ｂａｌｔｉｍｏｒｅ：Ｊｏｈｎｓ Ｈｏｐｋｉｎｓ Ｕｎｉｖｅｒｓｉｔｙ Ｐｒｅｓｓ）．
１１． Ｈｅｓｔｅｎｅｓ ＭＲ ＆ Ｓｔｉｅｆｅｌ Ｅ １９５２，“Ｍｅｔｈｏｄｓ ｏｆ Ｃｏｎｊｕｇａｔｅ Ｇｒａｄｉｅｎｔｓ ｆｏｒ Ｓｏｌｖｉｎｇ Ｌｉｎｅａｒ Ｓｙｓｔｅｍｓ”，Ｊ．Ｒｅｓ．Ｎａｔ．Ｂｕｒ．Ｓｔａｎｄａｒｄｓ，４９，４０９．
１２． Ｊｕ Ｌ，Ｄｕ Ｑ ＆ Ｇｕｎｚｂｕｒｇｅｒ Ｍ ２００２，“Ｐｒｏｂａｂｉｌｉｓｔｉｃ Ｍｅｔｈｏｄｓ ｆｏｒ Ｃｅｎｔｒｏｉｄａｌ Ｖｏｒｏｎｏｉ Ｔｅｓｓｅｌｌａｔｉｏｎｓ ａｎｄ ｔｈｅｉｒ Ｐａｒａｌｌｅｌ Ｉｍｐｌｅｍｅｎｔａｔｉｏｎｓ”，Ｐａｒａｌｌｅｌ Ｃｏｍｐｕｔｉｎｇ，２８，１４７７．
１３． Ｌｉｕ Ｙ，Ｗａｎｇ Ｗ，Ｌｅｖｙ Ｂ，Ｓｕｎ Ｆ，Ｙａｎ Ｄ−Ｍ，Ｌｕ Ｌ ＆ Ｙａｎｇ Ｃ ２００９，“Ｏｎ Ｃｅｎｔｒｏｉｄａｌ Ｖｏｒｏｎｏｉ Ｔｅｓｓｅｌｌａｔｉｏｎ−Ｅｎｅｒｇｙ Ｓｍｏｏｔｈｎｅｓｓ ａｎｄ Ｆａｓｔ Ｃｏｍｐｕｔａｔｉｏｎ”，ＡＣＭ Ｔｒａｎｓ．Ｇｒａｐｈｉｃｓ ２８，Ａｒｔｉｃｌｅ １０１．
１４．Ｐｉｎｙａ ＲＫ ＆ Ｐｕｅｔｔｅｒ ＲＣ １９９３，“Ｂａｙｅｓｉａｎ Ｉｍａｇｅ Ｒｅｃｏｎｓｔｒｕｃｔｉｏｎ：ｔｈｅ Ｐｉｘｏｎ ａｎｄ Ｏｐｔｉｍａｌ Ｉｍａｇｅ Ｍｏｄｅｌｉｎｇ”Ｐｕｂｌ．Ａｓｔｒｏｎ．Ｓｏｃ．Ｐａｃｉｆｉｃ １０５，６３０．
１５． Ｐｏｉｓｓｏｎ ＳＤ １８３７，Ｐｒｏｂａｂｉｌｉｔｅ ｄｅｓ Ｊｕｇｅｍｅｎｔｓ ｅｎ Ｍａｔｉｅｒｅ Ｃｒｉｍｉｎｅｌｌｅ ｅｔ ｅｎ Ｍａｔｉｅｒｅ Ｃｉｖｉｌｅ，Ｐｒｅｃｅｄｅｅｓ ｄｅｓ Ｒｅｇｌｅｓ Ｇｅｎｅｒａｌｅｓ ｄｕ Ｃａｌｃｕｌ ｄｅｓ Ｐｒｏｂａｂｉｌｉｔｅｓ （Ｐａｒｉｓ： Ｂａｃｈｅｌｉｅｒ），ｐ．２０６．
１６． Ｐｒｅｓｓ ＷＨ，Ｔｅｕｋｏｌｓｋｙ ＳＡ，Ｖｅｔｔｅｒｌｉｎｇ ＷＴ ＆ Ｆｌａｎｎｅｒｙ ＢＰ ２００７，Ｎｕｍｅｒｉｃａｌ Ｒｅｃｉｐｅｓ，３ｒｄ ｅｄ．，（Ｃａｍｂｒｉｄｇｅ：Ｃａｍｂｒｉｄｇｅ Ｕｎｉｖｅｒｓｉｔｙ Ｐｒｅｓｓ）．
１７． Ｐｕｅｔｔｅｒ ＲＣ，Ｇｏｓｎｅｌｌ ＴＲ ＆ Ｙａｈｉｌ Ａ ２００５，“Ｄｉｇｉｔａｌ Ｉｍａｇｅ Ｒｅｃｏｎｓｔｒｕｃｔｉｏｎ：Ｄｅｂｌｕｒｒｉｎｇ ａｎｄ Ｄｅｎｏｉｓｉｎｇ”，Ａｎｎｕ．Ｒｅｖ．Ａｓｔｒｏｎ．Ａｓｔｒｏｐｈｙｓ．，４３，１３９．
１８． Ｐｕｅｔｔｅｒ ＲＣ ＆ Ｙａｈｉｌ Ａ １９９９，“Ｔｈｅ Ｐｉｘｏｎ Ｍｅｔｈｏｄ ｏｆ Ｉｍａｇｅ Ｒｅｃｏｎｓｔｒｕｃｔｉｏｎ”，ｉｎ Ａｓｔｒｏｎｏｍｉｃａｌ Ｄａｔａ Ａｎａｌｙｓｉｓ Ｓｏｆｔｗａｒｅ ａｎｄ Ｓｙｓｔｅｍｓ ＶＩＩＩ，Ｍｅｈｒｉｎｇｅｒ ＤＭ，Ｐｌａｎｔｅ ＲＬ ＆ Ｒｏｂｅｒｔｓ ＤＡ，ｅｄｓ．，ＡＳＰ Ｃｏｎｆｅｒｅｎｃｅ Ｓｅｒｉｅｓ，１７２，３０７．
１９． Ｓａａｄ Ｙ ２００３，Ｉｔｅｒａｔｉｖｅ Ｍｅｔｈｏｄｓ ｆｏｒ Ｓｐａｒｓｅ Ｌｉｎｅａｒ Ｓｙｓｔｅｍｓ，２ｎｄ ｅｄ．（Ｐｈｉｌａｄｅｌｐｈｉａ： ＳＩＡＭ）．
２０． Ｓｃｈｗａｒｚ ＨＡ １８８８，“Ｕｅｂｅｒ ｅｉｎ Ｆｌａｅｃｈｅｎ ｋｌｅｉｎｓｔｅｎ Ｆｌａｅｃｈｅｎｉｎｈａｌｔｓ ｂｅｔｒｅｆｆｅｎｄｅｓ Ｐｒｏｂｌｅｍ ｄｅｒ Ｖａｒｉａｔｉｏｎｓｒｅｃｈｎｕｎｇ”，Ａｃｔａ Ｓｏｃｉｅｔａｔｉｓ Ｓｃｉｅｎｔｉａｒｕｍ Ｆｅｎｎｉｃａｅ ４５，３１８．
２１． Ｖｏｒｏｎｏｉ Ｇ １９０８，“Ｎｏｕｖｅｌｌｅｓ Ａｐｐｌｉｃａｔｉｏｎｓ ｄｅｓ Ｐａｒａｍｅｔｒｅｓ Ｃｏｎｔｉｎｕｓ ａ ｌａ Ｔｈｅｏｒｉｅ ｄｅｓ Ｆｏｒｍｅｓ Ｑｕａｄｒａｔｉｑｕｅｓ”，Ｊｏｕｒｎａｌ ｆｕｅｒ ｄｉｅ Ｒｅｉｎｅ ｕｎｄ Ａｎｇｅｗａｎｄｔｅ Ｍａｔｈｅｍａｔｉｋ １３３，９７．

Claims (21)

1. 対象物体のモデルを構築するコンピュータ実装方法であって、コンピュータに、
   対象物体を記述するデータおよび複数のパラメータを含んでいると共にノイズ部を有する入力信号をソースから受信し、
   前記複数のパラメータから初期パラメータのグループを選択し、
   二乗可積分な基底関数の組の線形結合を含むノンパラメトリック確率分布関数（ｐｄｆ）を推定し、
   の形式（θはパラメータを表し、ｘは観測値を表し、ｆ（ｘ、θ）はｐｄｆ）で二次尤度汎関数（ＱＬＦ）を計算し、
   前記ＱＬＦに基づいて前記データに対する前記初期パラメータの適合度を評価し、
   所定の条件が満たされるまで新たなパラメータのグループを選択して前記新たなパラメータのグループの適合度を評価することにより前記ＱＬＦを反復的に最適化して、
   最適化されたパラメータを用いて構築された前記対象物体のモデルを含む出力を生成する命令を実行させることを特徴とする方法。
2. 請求項１に記載のコンピュータ実装方法において、前記データが重みｗ_ｉを含み、前記ＱＬＦが、
   の形式を有していることを特徴とする方法。
3. 請求項１に記載のコンピュータ実装方法において、前記データおよび基底関数を用いてソース項を計算するステップを更に含んでいることを特徴とする方法。
4. 請求項３に記載のコンピュータ実装方法において、前記ＱＬＦが、
   前記基底関数を用いてグラム行列を計算し、
   前記ＱＬＦを生成するため、前記グラム行列、前記パラメータ、および前記ソース項を組み合せること、により得られることを特徴とする方法。
5. 請求項１に記載のコンピュータ実装方法において、前記入力信号が画像データであり、前記出力が、グラフィカルユーザーインターフェースで表示される前記対象物体の２次元、３次元、または４次元表現を含んでいることを特徴とする方法。
6. 請求項５に記載のコンピュータ実装方法において、前記画像データが、Ｘ線、ＣＴ、放射断層撮影、ＳＰＥＣＴ、およびＰＥＴからなるグループから選択され、前記対象物体が患者の身体部位であることを特徴とする方法。
7. 請求項６に記載のコンピュータ実装方法において、前記画像データが少なくとも２個の平面で取得され、前記出力が３次元表現を含んでいることを特徴とする方法。
8. 請求項６に記載のコンピュータ実装方法において、前記画像データが少なくとも２個の平面で取得されて更に時間を含み、前記出力が４次元表現を含んでいることを特徴とする方法。
9. 入力信号に含まれる対象物体を記述するデータをモデリングするシステムであって、
   コンピュータ可読媒体と、
   前記コンピュータ可読媒体に接続されたパラメータ最適化プロセッサと、
   前記パラメータ最適化プロセッサに接続されていて前記パラメータ最適化プロセッサとの間で再構築されたモデルの電子表現を受送信すべく適合された通信インターフェースとを含み、前記コンピュータ可読媒体に、前記パラメータ最適化プロセッサにより実行されたならば前記パラメータ最適化プロセッサに以下の動作、すなわち、
   被写体データを収集すべく構成されたソースから入力信号を受信し、
   前記被写体データに対応する初期パラメータのグループを生成し、
   二乗可積分な基底関数の組の線形結合を含むノンパラメトリック確率分布関数を推定し、
   の形式（θはパラメータを表し、ｘは観測値を表し、ｆ（ｘ、θ）はｐｄｆ）で二次尤度汎関数（ＱＬＦ）を計算し、
   前記ＱＬＦに基づいて前記データに対する前記初期パラメータの適合度を評価し、
   所定の条件が満たされるまで新たなパラメータのグループを選択して前記新たなパラメータのグループの適合度を評価することにより前記ＱＬＦを反復的に最適化して、
   最適化されたパラメータを用いて構築された前記対象物体のモデルを含む出力を生成するステップを含む動作を実行させるソフトウェア命令を保存していることを特徴とするシステム。
10. 請求項９に記載のシステムにおいて、前記データが重みｗ_ｉを含み、前記ＱＬＦが、
    の形式を有していることを特徴とするシステム。
11. 請求項９に記載のシステムにおいて、前記データおよび基底関数を用いてソース項を計算するステップを更に含んでいることを特徴とするシステム。
12. 請求項１１に記載のシステムにおいて、前記ＱＬＦが、
    前記基底関数を用いてグラム行列を計算し、
    前記ＱＬＦを生成するため、前記グラム行列、前記パラメータ、および前記ソース項を組み合せること、により得られることを特徴とするシステム。
13. 請求項９に記載のシステムにおいて、前記入力信号が画像データであり、前記出力が、グラフィカルユーザーインターフェースで表示される前記対象物体の２次元、３次元、または４次元表現を含んでいることを特徴とするシステム。
14. 請求項１３に記載のシステムにおいて、前記画像データが、Ｘ線、ＣＴ、放射断層撮影、ＳＰＥＣＴ、ＰＥＴからなるグループから選択され、前記対象物体が患者の身体部位であることを特徴とするシステム。
15. 請求項１４に記載のシステムにおいて、前記画像データが少なくとも２個の平面で取得され、前記出力が３次元表現を含んでいることを特徴とするシステム。
16. 請求項１３に記載のシステムにおいて、前記画像データが少なくとも２個の平面で取得されて更に時間を含み、前記出力が４次元表現を含んでいることを特徴とする方法。
17. データ成分およびノイズ成分を有する入力信号から対象物体の再構築された画像を生成する方法であって、コンピュータに、
    複数のパラメータを含む前記入力信号を画像ソースから受信し、
    前記複数のパラメータから初期パラメータのグループを選択し、
    二乗可積分な基底関数の組の線形結合を含むノンパラメトリック確率分布関数（ｐｄｆ）を推定し、
    の形式（θはパラメータを表し、ｘは観測値を表し、ｆ（ｘ、θ）はｐｄｆ）で二次尤度汎関数（ＱＬＦ）を計算し、
    前記ＱＬＦに基づいて前記データ成分に対する当該初期パラメータの適合度を評価し、
    所定の条件が満たされるまで新たなパラメータのグループを選択して前記新たなパラメータのグループの適合度を評価することにより前記ＱＬＦを反復的に最適化して、
    前記最適化されたパラメータに基づいて前記対象物体の再構築された画像の表示を含む出力を生成する命令を実行させることを特徴とする方法。
18. 請求項１７に記載の方法において、前記入力信号が第１の平面画像データおよび第２の平面画像データを含み、前記ＱＬＦを選択、推定、計算、評価して反復的に最適化するステップが、前記第１の平面画像データおよび第２の平面画像データの各々について実行され、出力を生成するステップが前記対象物体の３次元画像を表示するステップを含んでいることを特徴とする方法。
19. 請求項１８に記載の方法において、前記データ成分が重みｗ_ｉを含み、前記ＱＬＦが
    の形式を有していることを特徴とする方法。
20. 請求項１７に記載の方法において、前記データ成分および基底関数を用いてソース項を計算するステップを更に含んでいることを特徴とする方法。
21. 請求項２０に記載の方法において、前記ＱＬＦが、
    前記基底関数を用いてグラム行列を計算し、
    前記ＱＬＦを生成するため、前記グラム行列、前記パラメータ、および前記ソース項を組み合せること、により得られることを特徴とする方法。