(57)【要約】テストサンプルの特徴を予測することのできるシステムの構成要素の部分集合を識別するための方法および装置が記載されている。この方法は、構成要素と構成要素加重との線形結合を発生するステップを含んでおり、この線形結合において、各訓練サンプルが１つの既知の特徴を有している複数の訓練サンプルから発生されたデータから各構成要素に対する値が決定される。ある特徴の確率分布に対するモデルが規定され、このモデルは線形結合を条件とし、また、このモデルは２クラス反応に対する二項分布を、線形結合とこの反応の期待値とをリンクするプロビット関数と組合せたものではない。事前（ｐｒｉｏｒ）分布はゼロに近い高い確率密度を有するハイパープライアを含む線形結合の構成要素加重に対して構成され、また、事前分布とモデルとを組合せることにより、事後分布が生成される。この事後分布を最大にする構成要素加重を有する構成要素の部分集合が識別される。
2005524124.tif 000002
【特許請求の範囲】
【請求項１】
テストサンプルの特徴を予測することのできるシステムの構成要素の部分集合を識別する方法において、
（ａ）構成要素と構成要素加重との線形結合を発生し、この線形結合において、各訓練サンプルが１つの既知の特徴を有している複数の訓練サンプルから発生されたデータから各構成要素に対する値が導かれ、
（ｂ）ある特徴の確率分布に対するモデルを規定し、このモデルは線形結合を条件付けし、また、このモデルは２クラス反応に対する二項分布を、線形結合とこの反応の期待値とをリンクするプロビット関数と組合せたものではなく、
（ｃ）ゼロに近い高い確率密度を有するハイパープライア（ｈｙｐｅｒｐｒｉｏｒ）を含む線形結合の構成要素加重に対して事前（ｐｒｉｏｒ）分布を構成し、
（ｄ）この事前分布とモデルとを組合せて事後（ｐｏｓｔｅｒｉｏｒ）分布を生成し、
（ｅ）この事後分布を最大にする構成要素加重を有する構成要素の部分集合を識別するステップを含んでいる方法。
【請求項２】
モデルは、多項または二項ロジスティック回帰、一般化された線形モデル、コックスの比例ハザードモデル、およびパラメトリック生存モデルを含むグループから選択されたモデルに基づいた尤度関数である請求項１記載の方法。
【請求項３】
モデルは、多項または二項ロジスティック回帰に基づいた尤度関数である項１または２記載の方法。
【請求項４】
ロジスティック回帰は、多項または二項分布を有する特徴をモデル化する請求項２または３記載の方法。
【請求項５】
構成要素の部分集合はロジスティック回帰を規定することによってサンプルを複数の予め規定されたグループの１つに分類することが可能であり、このロジスティック回帰はサンプルを複数のサンプルグループにグループ化することを含んでおり、各サンプルグループは共通のグループ識別子を有している請求項１乃至４のいずれか１項記載の方法。
【請求項６】
ロジスティック回帰は、
【数１】
2005524124.tif 000252
の形態のものであり、ここで、
ｘiTβgは、構成要素加重βgを有する訓練サンプルｉからの入力データから発生された線形結合であり、
ｘiTは、Ｘのｉ番目の行に対する構成要素であり、βgはサンプルクラスｇに対する構成要素加重のセットであり、
訓練サンプルｉがクラスｇのメンバーである場合、ｅig＝１であり、そうでない場合にはｅig＝０であり、
Ｘは、ｐ個の構成要素を含むｎ個の訓練サンプルからのデータである請求項１乃至５のいずれか１項記載の方法。
【請求項７】
構成要素の部分集合は、ロジスティック回帰を規定することによってサンプルをクラスに分類することができ、このクラスは複数の予め規定された順序付けられたクラスの１つであり、前記ロジスティック回帰は、各グループ識別子が順序付けられたクラスのメンバーに対応している一連のグループ識別子を規定し、サンプルを順序付けられたクラスの１つにグループ化することを含んでいる請求項１または２記載の方法。
【請求項８】
ロジスティック回帰は、
【数２】
2005524124.tif 000253
の形態のものであり、ここで、
【数３】
2005524124.tif 000254
は、訓練サンプルｉがｋ以下である識別子を有するクラスに属している確率であり（ここで、順序付けられたクラスの合計はＧである）、
【数４】
2005524124.tif 000255
は、構成要素加重：
【数５】
2005524124.tif 000256
を有する訓練サンプルｉからの入力データから発生された線形結合であり、
Ｘは、ｐ個の構成要素を含むｎ個の訓練サンプルからのデータであり、
ｘiTは、Ｘのｉ番目の行に対する構成要素であり、
ｒijは、
【数６】
2005524124.tif 000257
として定義されており、ここで、
観察ｉがクラスｊ中のものである場合はｃij＝１であり、そうでない場合はｃij＝０である請求項７記載の方法。
【請求項９】
モデルは、一般化された線形モデルに基づいた尤度関数である請求項１または２記載の方法。
【請求項１０】
一般化された線形モデルは、分布の正規指数族として分布された特徴をモデル化する請求項９記載の方法。
【請求項１１】
分布の正規指数族は、正規分布、ガウス分布、ポアソン分布、指数分布、ガンマ分布、ｘ2分布、および逆ガンマ分布からなるグループから選択される請求項１０記載の方法。
【請求項１２】
構成要素の部分集合は、一般化された線形モデルを規定することによってサンプルの予め規定された特性を予測することができ、この一般化された線形モデルは予測される特性をモデル化することを含んでいる請求項１または２記載の方法。
【請求項１３】
一般化された線形モデルは、
【数７】
2005524124.tif 000258
の形態のものであり、ここで、
ｙ＝（ｙ1，・・・，ｙn）Tであり、ｙiはｉ番目のサンプルに関して測定された特性であり、
ａi（φ）＝φ／ｗiであり、ｗiは既知の加重の固定されたセットであり、φは単一の
スケールパラメータであり、
関数ｂ（．）およびｃ（．）はNelderおよびWedderburn氏等による文献（1972）において定義されたものであり、
Ｅ｛ｙi｝＝ｂ’（θi）
Ｖａｒ｛ｙ｝＝ｂ’’（θi）ａi（φ）＝τi2ａi（φ）であり、
ここで、各観察はコバリエイト（ｃｏｖａｒｉａｔｅ）ｘiと線形予測子ηi＝ｘiTβとのセットを有している請求項９または１０記載の方法。
【請求項１４】
モデルは、コックスの比例ハザードモデル、パラメトリック生存モデルおよび加速された生存時間モデルからなるグループから選択されたモデルに基づいた尤度関数である請求項１または２記載の方法。
【請求項１５】
構成要素の部分集合は、コックスの比例ハザードモデル、パラメトリック生存モデルまたは加速された生存時間モデルに基づいて尤度を規定することにより、あるサンプルに対するある事象までの時間を予測することができ、この尤度は、複数のサンプルに対して、そのサンプルが得られた時間からその事象の時間までに経過した時間を測定することを含んでいる請求項１記載の方法。
【請求項１６】
コックスの比例ハザードモデルは、
【数８】
2005524124.tif 000259
の形態のものであり、ここで、
Ｘは、ｐ個の構成要素を含むｎ個の訓練サンプルからのデータであり、
Ｚは、Ｘの行の再配列であるマトリックスであり、ここで、Ｚの行の順序付けは生存時間の順序付けにより誘発された順序付けに対応し、
ｄは、生存時間を順序付けるために必要とされる同じ置換により検閲削除インデックスを順序付けた結果であり、
ＺjはマトリックスＺのｊ番目の行であり、ｄjはｄのｊ番目のエレメントであり、
【数９】
2005524124.tif 000260
Ｒj＝｛ｉ：ｉ＝ｊ，ｊ＋１，・・・，Ｎ｝＝ｊ番目の順序付けられた事象時間ｔ(j)におけるリスクセットである請求項１４記載の方法。
【請求項１７】
パラメトリックハザードモデルは、
【数１０】
2005524124.tif 000261
の形態のものであり、ここで、
【数１１】
2005524124.tif 000262
ｉ番目のサンプルが検閲削除されていない場合はｃi＝１、ｉ番目のサンプルが検閲削除されていない場合はｃi＝０であり、
関数λ（．）およびΛ（．）は、AitkinおよびClayton氏等による文献（1980）において定義されたものであり、
ｘjはＸのｉ番目の行であり、Ｘはｐ個の構成要素を含むｎ個の訓練サンプルからのデータである請求項１４記載の方法。
【請求項１８】
事前分布は、
【数１２】
2005524124.tif 000263
の形態のものであり、
ここで、ｐ（β|ｖ2）はＮ（０，ｄｉａｇ｛ｖ2｝）であり、
ｖはハイパーパラメータであり、
ｐ（ｖ2）はハイパープライア分布である請求項１乃至１７のいずれか１項記載の方法。
【請求項１９】
ハイパープライアは、
【数１３】
2005524124.tif 000264
の形態のジェフェリーのプライア（ｐｒｉｏｒ）である請求項１乃至１８のいずれか１項記載の方法。
【請求項２０】
事後分布は、
【数１４】
2005524124.tif 000265
の形態のものであり、
ここで、
【数１５】
2005524124.tif 000266
は尤度関数である請求項１乃至１９のいずれか１項記載の方法。
【請求項２１】
事後分布は、反復的な処理手続きを使用して最大化される請求項１乃至２０のいずれか１項記載の方法。
【請求項２２】
反復的な処理手続きは、ＥＭアルゴリズムである請求項２１記載の方法。
【請求項２３】
システムは、生物学的システムである請求項１乃至２２のいずれか１項記載の方法。
【請求項２４】
生物学的システムは、バイオテクノロジーアレイである請求項２３記載の方法。
【請求項２５】
バイオテクノロジーアレイは、ＤＮＡアレイ、蛋白質アレイ、抗体アレイ、ＲＮＡアレイ、炭水化物アレイ、化学物質アレイ、脂質アレイからなるグループから選択される請求項２４記載の方法。
【請求項２６】
各グループがテスト治療に対する反応によって規定されている複数の予め規定されたグループの１つに部分集合を分類することのできる被験者の構成要素の部分集合を識別する方法において、
（ｄ）複数の被験者をテスト治療に露出し、そのテスト治療に対する反応に基づいて被験者を反応グループにグループ化し、
（ｅ）被験者の構成要素を測定し、
（ｆ）請求項１乃至２５のいずれか１項記載の方法を使用して被験者を反応グループに分類することのできる構成要素の部分集合を識別するステップを含んでいる方法。
【請求項２７】
構成要素は、遺伝子、小さいヌクレオチド多型（ＳＮＰｓ）、蛋白質、抗体、炭水化物、脂質からなるグループから選択される請求項２６記載の方法。
【請求項２８】
システムから発生されたデータからのテストサンプルの特徴を、そのシステムからの複数のサンプルから予測することのできるシステムの構成要素の部分集合を識別する装置において、
（ａ）各訓練サンプルが１つの既知の特徴を有している複数の訓練サンプルから発生されたデータから各構成要素に対する値が導かれる、構成要素と構成要素加重との線形結合を発生する手段と、
（ｂ）線形結合を条件付けし、また、２クラス反応に対する二項分布を、線形結合とこの反応の期待値とをリンクするプロビット関数と組合せたものではない、ある特徴の確率分布に対するモデルを規定する手段と、
（ｃ）ゼロに近い高い確率密度を有するハイパープライアを含む線形結合の構成要素加重に対して事前分布を構成する手段と、
（ｄ）この事前分布とモデルとを組合せて事後分布を生成する手段と、
（ｅ）この事後分布を最大にする構成要素加重を有する構成要素の部分集合を識別する手段とを具備している装置。
【請求項２９】
計算装置にロードされたときにこの計算装置を制御して、請求項１乃至２７のいずれか１項記載の方法を実行するように構成されているコンピュータプログラム。
【請求項３０】
請求項１乃至２７のいずれか１項記載の方法により実施される請求項２９記載のコンピュータプログラム。
【請求項３１】
請求項２９または３０記載のコンピュータプログラムを提供するコンピュータ読出し可能な媒体。
【請求項３２】
システムからのサンプルをテストしてそのサンプルの特徴を識別する方法において、
特徴の診断に役立つ構成要素の部分集合をテストするステップを含んでおり、
構成要素の部分集合は請求項１乃至２７のいずれか１項記載の方法によって決定されている方法。
【請求項３３】
サンプルの特徴を決定するためにシステムからのサンプルをテストする装置において、
請求項１乃至２７のいずれか１項記載の方法にしたがって識別された構成要素をテストする手段を備えている装置。
【請求項３４】
計算装置上で実行されたときに、システムからのテストサンプルの特徴を予測することのできるそのシステムからの構成要素を識別する方法でこの計算装置を制御するように構成されており、構成要素と構成要素加重との線形結合は、各訓練サンプルが１つの既知の特徴を有する複数の訓練サンプルから発生されたデータから発生され、ゼロに近い高い確率密度を有するハイパープライアを含む構成要素加重に対する事前分布と、線形結合に条件付けするモデルとを組合せることによって事後分布が発生され、それによって事後分布を最大にする構成要素加重を推定し、前記モデルは２つのクラス反応に対する二項分布を、線形結合とこの反応の期待値とをリンクするプロビット関数と組合せたものではないコンピュータプログラム。
【請求項３５】
生物学的システムからのテストサンプルの特徴を予測することのできる生物学的システムの構成要素の部分集合を識別する方法において、
（ａ）構成要素と構成要素加重との線形結合を発生し、この線形結合において、各訓練サンプルが１つの既知の特徴を有している複数の訓練サンプルから発生されたデータから各構成要素に対する値が決定され、
（ｂ）ある特徴の確率分布に対するモデルを規定し、このモデルは線形結合に条件付けし、
（ｃ）ゼロに近い高い確率密度を有するハイパープライアを含む線形結合の構成要素加重に対して事前分布を構成し、
（ｄ）この事前分布とモデルとを組合せて事後分布を生成し、
（ｅ）この事後分布を最大にする構成要素加重を有する構成要素の部分集合を識別するステップを含んでいる方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、システムのサンプルから発生されたデータから、そのシステム内のサンプルの特徴を予測することのできるシステムの構成要素を識別するための方法および装置に関し、とくに、本発明は、それに限定されるものではないが、生物学的方法によって発生されたデータから、その生物学的システムからのサンプルと関連付けられた重要な特徴を予測することのできる生物学的システムの構成要素を識別するための方法および装置に関する。
【背景技術】
【０００２】
種々の重要な特徴に分類されることのできる多くの“システム”が存在する。この“システム”という用語には、本質的に、データが提供されることのできる全てのタイプのシステムが含まれ、それには化学的システム、金融システム（たとえば、個人、グループまたは法人向けクレジットシステム、ローン履歴）、地質学的システムおよびその他多数が含まれる。システムのサンプルから特定の特徴を識別するために［たとえば、金融システムの解析により、金融システム内に存在するグループを識別する（信用が“申し分のない”者および信用貸しは危険である者という非常に簡単な関係で）ことを補助するために］システムから発生されたデータ（たとえば、統計的データ）を使用できることが望ましい。大量の統計的データが存在している場合、そのシステムからのサンプルの特定の特徴を予測するそのデータから構成要素を識別することは困難な作業である。これは一般に、処理すべき大量のデータの大部分が、そのデータが得られた特定のサンプルの重要な特徴を全く、あるいはほとんど示さない可能性があるためである。さらに、テストサンプルデータが訓練サンプルデータに関して高度の可変性を有しているときには、その訓練サンプルデータを使用して識別された構成要素はそのテストサンプルデータに関する特徴の識別には役に立たないことが多い。これは、たとえば、データが多数の異なったソースから得られる状況においてよくあるケースであり、それは、個々のソースのそれぞれからデータが収集される条件を制御することが不可能なことが多いためである。
【０００３】
これらの問題がとくに関係するシステムのタイプの一例は生物学的システムであり、以下の説明はとくに生物学的システムに該当するものである。しかしながら、本発明は生物学的システムによる使用に限定されず、任意のシステムに対して一般的に適用される。
【０００４】
バイオテクノロジーにおける近年の進歩の結果、システムの大規模スクリーニングおよびサンプルの解析に対する生物学的方法が発達してきた。このような方法には、たとえば、ＤＮＡ、ＲＮＡまたは抗体マイクロアレイ解析、プロテオミクス解析、プロテオミクス電気泳動ゲル解析および高度なスループットのスクリーニング技術が含まれる。これらのタイプの方法の結果、テストされたサンプルのそれぞれに対して３０，０００以上の構成要素を有することのできるデータがしばしば発生される。
【０００５】
生物学的システムからのサンプル中の重要な特徴を識別できることが明らかに重要である。たとえば、“病気にかかっている”および“病気にかかっていない”のようなグループを分類すること等が挙げられる。これらの生物学的方法の多くは、生物学的システム中のサンプルの特徴を予測する診断ツール（たとえば、組織または体液をスクリーニングすることにより病気を識別する）として（、あるいは、たとえば、薬剤化合物の効能を決定するツールとして）有用である。
【発明の開示】
【発明が解決しようとする課題】
【０００６】
今日、このような適用におけるバイオテクノロジーアレイのような生物学的方法の使用は、これらのタイプの方法から発生された大量のデータと、および意味のある結果を得るためにデータをスクリーニングする効率的な方法の欠如のために制限されている。その結果、従来技術の方法を使用した生物学的データの解析は、データ中の情報を十分に利用することができないか、あるいは時間を要し、誤った陽性および陰性の結果を出し易く、また、データから意味のある結果が得られなければならない場合には大量のコンピュータメモリを必要とする。これは、迅速で正確なスクリーニングが要求される規模スクリーニングシナリオにおいて問題である。
【０００７】
したがって、とくに生物学的データを解析する改善された方法が必要とされており、さらに一般的には、システムからのサンプルに対する重要な特徴を予測するために任意のシステムからのデータを解析する改善された方法が必要とされている。
【課題を解決するための手段】
【０００８】
第１の特徴において、本発明は、テストサンプルの特徴を予測することのできるシステムの構成要素の部分集合を識別する方法を提供し、この方法は、
（ａ）構成要素と構成要素加重との線形結合を発生し、この線形結合において、各訓練サンプルが１つの既知の特徴を有している複数の訓練サンプルから発生されたデータから各構成要素に対する値が決定され、
（ｂ）ある特徴の確率分布に対するモデルを規定し、このモデルは線形結合を条件付けし、また、このモデルは２クラス反応に対する二項分布を、線形結合とこの反応の期待値とをリンクするプロビット関数と組合せたものではなく、
（ｃ）ゼロに近い高い確率密度を有するハイパープライア（ｈｙｐｅｒｐｒｉｏｒ）を含む線形結合の構成要素加重に対して事前（ｐｒｉｏｒ）分布を構成し、
（ｄ）この事前分布とモデルとを組合せて事後（ｐｏｓｔｅｒｉｏｒ）分布を生成し、
（ｅ）この事後分布を最大にする構成要素加重を有する構成要素の部分集合を識別するステップを含んでいる。
【０００９】
この方法は、訓練サンプルに対する特徴を予測することのできる構成要素の部分集合を識別するために既知の特徴を有する訓練サンプルを使用する。続いて、たとえば、臨床テストのようなテストに対して構成要素の部分集合情報を使用して、組織サンプルが悪性であるかまたは良性であるか、あるいは腫瘍の重量はどれくらいか等の特徴を予測することができ、あるいは特定の症状を有する患者の推定生存時間を提供することができる。ここで使用されているように、“特徴”という用語は、あるサンプルと関連付けられた任意の応答または識別可能な特性または特質のことである。たとえば、ある特徴は、特定のサンプルに対するある事象までの特定の時間であってもよいし、あるいはあるサンプルの寸法または量、もしくはあるサンプルが分類されることのできるクラスまたはグループであってもよい。
【００１０】
本発明の方法は、ベイズの統計的方法を使用して構成要素の加重を推定する。大量の構成要素がシステムから発生された場合（通常、本発明の方法が有効であるケースである）、本発明は、大部分の構成要素は、ある特徴を予測するための構成要素の部分集合の一部を形成している構成要素ではなさそうであるというアプリオリ仮定を行うことが好ましい。したがって、大部分の構成要素の加重はゼロであるようだという仮定がなされる。この仮定を念頭に置いて加重の事後（ｐｏｓｔｅｒｉｏｒ）確率が最大化されるように構成要素加重を設定したあるモデルが構成される。予め定められたしきい値より下の加重を有する構成要素（アプリオリ仮定によるとそれらの大部分となるであろう）はその適用を免除される。このプロセスは、残った診断構成要素が識別されるまで繰返される。この方法は、主に、大部分の構成要素を迅速に除去するアプリオリ仮定のために迅速である。
【００１１】
あるシステムのほとんどの特徴は典型的に確率分布を示し、ある特徴の確率分布は、訓練サンプルから発生されたデータに基づいた統計モデルを使用してモデル化されることができる。本発明の方法は、１つのまたは一連の重要な特徴に対する確率分布をモデル化する統計モデルを使用する。したがって、特定の確率分布を有する１つの重要な特徴に対して、その分布をモデル化した１つの適切なモデルが規定される。この方法は、線形結合を条件とする任意のモデルを使用することができ、また、訓練サンプルから得られたデータに基づいて確率分布を提供する尤度関数の形態の数学的方程式であることが好ましい。この尤度関数は、ある確率分布を記述するために上述されたモデルに基づいていることが好ましい。１実施形態において、モデルは、多項または二項ロジスティック回帰、一般化された線形モデル、コックスの比例ハザードモデル、加速された故障（ｆａｉｌｕｒｅ）モデル、パラメトリック生存モデル、ｘ2分布モデル、または指数分布モデルからなるグループから選択されたモデルに基づいた尤度関数である。
【００１２】
１実施形態において、尤度関数は多項または二項ロジスティック回帰に基づいている。多項または二項ロジスティック回帰は、多項または二項分布を有する特徴をモデル化することが好ましい。二項分布は、オン／オフ状態のような２つの可能性のあるクラスまたはグループを有する統計的分布である。このようなグループの例には、死亡／生存、改善された／改善されていない、衰弱した／衰弱していない、が含まれる。多項分布は、複数のクラスまたはグループが複数のサンプルのそれぞれに対して可能である、換言すると、あるサンプルが複数のクラスまたはグループの１つに分類されることのできる二項分布を一般化したものである。したがって、多項または二項ロジスティック回帰に基づいて尤度関数を規定することにより、サンプルを複数のグループまたはクラスの１つに分類することのできる構成要素の部分集合を識別することが可能である。これを行うために、訓練サンプルはそれの予め定められた特徴に基づいて複数のサンプルグループ（または“クラス”）にグループ化され、この予め定められた特徴において、各サンプルのメンバーは共通の特徴を有しており、共通のグループ識別子を割当てられる。尤度関数は、線形結合（グループ化された訓練サンプルから発生されたデータを含む）を条件とする多項または二項ロジスティック回帰に基づいて公式化される。特徴は、訓練サンプルがグループ化される任意の所望の分類であってもよい。たとえば、組織サンプルを分類する特徴は、その組織が正常、悪性または良性であることでもよいし、あるいはその細胞が白血病細胞または健康な細胞であることでもよいし、もしくは訓練サンプルがある症状を有する患者または有しない患者のいずれの血液から得られたことでもよく、あるいは訓練サンプルが正常な細胞と比較されるいくつかのタイプの癌の１つからの細胞から得られたことであってもよい。
【００１３】
尤度関数は、
【数１６】
2005524124.tif 000003
【００１４】
の形態のロジスティック回帰に基づいていることが好ましく、ここで、
ｘiTβgは、構成要素加重βgを有する訓練サンプルｉからの入力データから発生された線形結合であり、
ｘiTは、Ｘのｉ番目の行に対する構成要素であり、βgはサンプルクラスｇに対する構成要素加重のセットであり、
訓練サンプルｉがクラスｇのメンバーである場合、ｅig＝１であり、そうでない場合にはｅig＝０であり、
Ｘは、ｐ個の構成要素を含むｎ個の訓練サンプルからのデータである。
【００１５】
別の実施形態において、尤度関数は順序付けられたカテゴリー的ロジスティック回帰に基づいている。この順序付けられたカテゴリー的ロジスティック回帰は、クラスが特定の順序にされている（たとえば、病気の重大度が増加または減少していくクラスのような順序付けられたクラスの）多項分布をモデル化する。順序付けられたカテゴリー的ロジスティック回帰に基づいて尤度関数を規定することにより、あるサンプルをあるクラスに分類することのできる構成要素の部分集合を識別することが可能であり、この場合そのクラスは複数の予め規定された順序付けられたクラスの１つである。各グループ識別子がある順序付けられたクラスのあるメンバーに対応している一連のグループ識別子を規定し、訓練サンプルをその訓練サンプルの予め定められた特徴に基づいてその順序付けられたクラスの１つにグループ化することにより、尤度関数は、線形結合（グループ化された訓練サンプルから発生されたデータを含む）を条件とする順序付けられたカテゴリー的ロジスティック回帰に基づいて公式化されることができる。
【００１６】
尤度関数は、
【数１７】
2005524124.tif 000004
【００１７】
の形態の順序付けられたカテゴリー的ロジスティック回帰に基づいていることが好ましく、ここで、
γikは、訓練サンプルｉがｋ以下である識別子を有するクラスに属している確率であり（ここで、順序付けられたクラスの合計はＧである）、
ｘiTβ*は、構成要素加重β*iを有する訓練サンプルｉからの入力データから発生された線形結合であり、
ｘiTは、Ｘのｉ番目の行に対する構成要素であり、
ｒijは、
【数１８】
2005524124.tif 000005
【００１８】
として定義されており、ここで、
観察ｉがクラスｊ中のものである場合はｃij＝１であり、そうでない場合はｃij＝０である。
【００１９】
本発明の別の実施形態において、尤度関数は一般化された線形モデルに基づいている。一般化された線形モデルは、分布の正規指数族に属する分布を有する特徴をモデル化することが好ましい。正規指数族分布の例には、正規分布、ガウス分布、ポアソン分布、ガンマ分布および逆ガンマ分布が含まれる。したがって、本発明の方法の別の実施形態では、予測されるべき特性をモデル化する一般化された線形モデルを規定することにより、分布の正規指数族の範囲内であるサンプルの予め規定された特性を予測することのできる構成要素の部分集合が識別される。一般化された線形モデルを使用して予測されることのできる特性の例には、たとえば、あるサンプルの重量、寸法、カウント、グループメンバーシップその他のディメンションまたは量あるいは特性のような指定された分布を示すあるサンプルの任意の量が含まれる。
【００２０】
一般化された線形モデルは、
【数１９】
2005524124.tif 000006
【００２１】
の形態のものであることが好ましく、ここで、
ｙ＝（ｙ1，・・・，ｙn）Tであり、ｙiはｉ番目のサンプルに関して測定された特性であり、
ａi（φ）＝φ／ｗiであり、ｗiは既知の加重の固定されたセットであり、φは単一の
スケールパラメータであり、
関数ｂ（．）およびｃ（．）はNelderおよびWedderburn氏等による文献（1972）において定義されたものであり、
Ｅ｛ｙi｝＝ｂ’（θi）
Ｖａｒ｛ｙ｝＝ｂ’’（θi）ａi（φ）＝τi2ａi（φ）であることが好ましい。
【００２２】
各観察はコバリエイトｘiと線形予測子ηi＝ｘiTβとのセットを有していることが好ましい。ｉ番目の観察の平均とそれの線形予測子との間の関係は、リンク関数ηi＝ｇ（μi）＝ｇ（ｂ’（θi））によって与えられることが好ましい。このリンクの逆数はｈによって表され、それは：
Ｅ｛ｙi｝＝ｂ’（θi）＝ｈ（ηi）
であることが好ましい。
【００２３】
別の実施形態において、本発明の方法は、ある事象がデータを得た時点では発生していないと仮定してその事象までの時間の確率を推定することが好ましいハザードモデルに基づいた尤度関数を使用することにより、あるサンプルに対して事象までの時間を予測するために使用されてもよい。１実施形態において、尤度関数は、コックスの比例ハザードモデル、パラメトリック生存モデルおよび加速された故障時間モデルからなるグループから選択されたモデルに基づいている。コックスの比例ハザードモデルは、ハザード関数の形式に関して制限的仮定を行わずに、ある事象までの時間が構成要素と構成要素加重のセットに関してモデル化されることを可能にする。加速された故障モデルは、構成要素の測定値が時間スケールに乗法的に作用すると仮定されているので、個体が時間軸に沿って進行する速度に影響を与える生存時間からなるデータに対する一般的なモデルであることが好ましい。したがって、加速された生存モデルは、たとえば、病気の進行速度に関して解釈されることができる。パラメトリック生存モデルは、ある事象までの時間（たとえば、生存時間）に対する分布関数が既知の分布によってモデル化されているか、あるいは特定されたパラメトリック公式を有するものである。一般に使用されている生存分布には、ワイブル指数および極値分布がある。
【００２４】
あるサンプルに対するある事象までの時間を予測することのできる構成要素の部分集合は、コックスの比例ハザードモデル、パラメトリック生存モデルまたは加速された生存時間モデルに基づいて尤度を規定することにより識別され、この尤度は、複数のサンプルに対してそのサンプルが得られた時間から事象の時間までの経過した時間を測定することを含んでいる。
【００２５】
ある事象までの時間を予測する尤度関数は、以下の形態のものであることが好ましい：
【数２０】
2005524124.tif 000007
【００２６】
コックスの比例ハザードモデルに基づいた尤度関数は、
【数２１】
2005524124.tif 000008
【００２７】
の形態のものであることが好ましく、ここで、
ＺはＸの行の再配列であるマトリックスであることが好ましく、ここでＺの行の順序付けは生存時間の順序付けにより誘発された順序付けに対応し、ｄは生存時間を順序付けるために必要とされる同じ置換により検閲削除インデックスを順序付けた結果である。また、ＺjはマトリックスＺのｊ番目の行であり、ｄjはｄのｊ番目のエレメントであり、ここで、
【数２２】
2005524124.tif 000009
【００２８】
であり、Ｒj＝｛ｉ：ｉ＝ｊ，ｊ＋１，・・・，Ｎ｝＝ｊ番目の順序付けられた事象時間ｔ(j)におけるリスクセットである。
【００２９】
パラメトリック生存モデルに基づいた対数尤度関数は、以下の形態のものであることが好ましい：
【数２３】
2005524124.tif 000010
【００３０】
ｉ番目のサンプルが検閲削除されていない場合はｃi＝１、ｉ番目のサンプルが検閲削除されていない場合はｃi＝０である。尤度関数のこの形態は、ワイブル指数および極値分布によって共用される。関数λ（．）およびΛ（．）は、AitkinおよびClayton氏等による文献（1980）において定義されたものである。
【００３１】
任意の規定されたモデルに対して、典型的にベイズの統計的モデルを使用して構成要素加重が推定され（KotzおよびJohnson氏等による文献、1983）、その場合、尤度関数と事前分布を組合せる構成要素加重の事後分布が公式化される。構成要素加重は、各訓練サンプルに対して発生されたデータを与えられた加重の事後分布を最大化することによって推定される。したがって、最大化されるべき目標関数は、上述した特徴に対するモデルと加重に対する事前分布とに基づいた尤度関数からなる。
【００３２】
事前分布は、
【数２４】
2005524124.tif 000011
【００３３】
の形態のものであることが好ましく、
ここで、ｖはハイパーパラメータのｐ×１ベクトルであり、また、ｐ（β|ｖ2）はＮ（０，ｄｉａｇ｛ｖ2｝）であり、ｐ（ｖ2）はｖ2に対するあるハイパープライア分布である。このハイパープライア分布（この方法の全ての実施形態に対して同じであることが好ましい）は、種々の表記規定を使用して表現されてもよく、好ましい実施形態の詳細な説明（以下、参照）においては、以下の表記規定は単に特定の好ましい実施形態に対する便宜性のために採用されたに過ぎない。
【００３４】
ここで使用されているように、確率分布に対する尤度関数が多項または二項ロジスティック回帰に基づいているとき、事前分布に対する表記は：
【数２５】
2005524124.tif 000012
【００３５】
また、ｐ（βg|τg2）はＮ（０，ｄｉａｇ｛τg2｝）であり、ｐ（τg2）はτg2に対するあるハイパープライア分布である。
【００３６】
ここで使用されているように、確率分布に対する尤度関数が順序付けられたカテゴリー的ロジスティック回帰に基づいているとき、事前分布に対する表記は：
【数２６】
2005524124.tif 000013
【００３７】
ここで、β1，β2，・・・，βnは構成要素加重であり、ｐ（βi|τi）はＮ（０，τi2）であり、ｐ（τi）はτiに対するあるハイパープライア分布である。
【００３８】
ここで使用されているように、分布に対する尤度関数が一般化された線形モデルに基づいているとき、事前分布に対する表記は：
【数２７】
2005524124.tif 000014
【００３９】
ここで、ｖはハイパーパラメータのｐ×１ベクトルであり、また、ｐ（β|ｖ2）はＮ（０，ｄｉａｇ｛ｖ2｝）であり、ｐ（ｖ2）はｖ2に対するある事前分布である。
【００４０】
ここで使用されているように、分布に対する尤度関数がハザードモデルに基づいているとき、事前分布に対する表記は：
【数２８】
2005524124.tif 000015
【００４１】
ここで、ｐ（β*|ｖ2）はＮ（０，ｄｉａｇ｛ｖ2｝）であり、ｐ（ｖ2）はｖ2に対するあるハイパープライア分布である。
【００４２】
事前分布は、可能な限り常に０加重が好ましいことを保証するハイパープライアを含んでいる。  
ハイパープライアは、ジェフェリーのハイパープライアである（KotzおよびJohnson氏等の文献、1983）。
【００４３】
上述したように、事後分布を生成するために事前分布と尤度関数が組合せられる。事後分布は以下の形態のものであることが好ましい： 
【数２９】
2005524124.tif 000016
【００４４】
事後分布における構成要素加重は、事後分布の確率密度が最大化されるように反復的な処理手続きで推定されることが好ましい。反復的な処理手続き中、予め定められたしきい値より小さい値を有する構成要素加重は、好ましくはこれらの構成要素をゼロに設定することによって推定される。この結果、対応した構成要素が除去される。
【００４５】
反復的な処理手続きは、ＥＭアルゴリズムであることが好ましい。ＥＭアルゴリズムは、事後分布の確率密度を最大化する所定の構成要素加重に収斂する構成要素加重推定値のシーケンスを処理する。ＥＭアルゴリズムは、Ｅまたは期待ステップおよびＭまたは最大化ステップとして知られている２つのステップから構成されている。Ｅステップにおいて、観察されたデータおよび現在のパラメータ値を条件とする対数−事後関数期待値が決定される。Ｍステップにおいて、尤度を増加させる更新された構成要素加重推定値を求めるために、期待された対数−事後関数が最大化される。２つのステップは、ＥステップおよびＭステップの収斂が達成されるまで、換言すると、期待値および対数−事後関数の最大化された値が収斂するまで交互に実施される。
【００４６】
本発明の方法は、測定値が得られることのできる任意のシステムに対して、好ましくは、非常に大量のデータが発生されるシステムに対して適用可能であると考えられる。本発明の方法が適用可能なシステムの例には、たとえば、クレジットリスク査定システム、保険システム、マーケティングシステムまたは企業記録システム、電子システム、物理的システム、天体物理学的システムおよび機械的システムが含まれる。たとえば、金融システムにおいては、サンプルは特定の株であることができ、構成要素は、企業収益、従業員数、株主の数等の、株価に影響を与える可能性のある任意の数のファクタに対してなされた測定値であることができる。
【００４７】
本発明の方法は、とくに、生物学的システムの解析における使用に適している。本発明の方法は、構成要素に対する測定可能な値を生成し、その構成要素が特有のラベルを付けられる任意の生物学的システムからのサンプルを分類するために、構成要素の部分集合を識別するために使用されることができる。換言すると、構成要素は、１つの構成要素からのデータと別の構成要素からのデータとの区別を可能にするようにラベル付けされるか、あるいは編成される。たとえば、構成要素は、たとえば、各構成要素からのデータと別のものとを空間的な位置で区別することをを可能にするアレイで空間的に編成されてもよいし、あるいは各構成要素は、それと関連付けられた、識別信号またはタグのようなある特有の識別物を有していてもよい。たとえば、構成要素は個々のキャリアに結び付けられてもよく、各キャリアは量子ドット［たとえば、Rosenthal氏による文献（2001，Nature Biotech 19: 621-622）；Han氏他による文献（2001，Nature Biotechnology 19: 631-635）を参照］、蛍光マーカ［たとえば、Fu氏他による文献（1999，Nature Biotech 17: 1109-1111）を参照］、バーコード付きタグ［たとえば、LockhartおよびTrulson氏等による文献（2001，Nature Biotechnology 19: 1122-1123）を参照］のような検出可能な識別シグナチャを有している。
【００４８】
とくに好ましい実施形態において、生物学的システムはバイオテクノロジーアレイである。バイオテクノロジーアレイ［それらの例はSchena氏他による文献（1995，Science 270: 467-470）；Lockhart氏他による文献（1996，Nature Biotechnology 14: 1649）；米国特許第 5,569,5880号明細書に記載されている］の例には、オリゴヌクレオチドアレイ、ＤＮＡアレイ、ＤＮＡマイクロアレイ、ＲＮＡアレイ、ＲＮＡマイクロアレイ、ＤＮＡマイクロチップ、ＲＮＡマイクロチップ、蛋白質アレイ、蛋白質マイクロチップ、抗体アレイ、化学物質アレイ、炭水化物アレイ、プロテオミクスアレイ、脂質アレイが含まれる。別の実施形態において、生物学的システムは、たとえば、ＤＮＡまたはＲＮＡ電気泳動ゲル、蛋白質またはプロテオミクス電気泳動ゲル、Ｂｉａｃｏｒｅ解析、アミノ酸解析、ＡＤＭＥＴｏｘスクリーニング［たとえば、High-throughput ADMETox estimation： In Vitro and In Silico approaches(2002)，Ference Darvas and Gyorgy Dorman(Eds)，Biotechniques Pressを参照］のような生体分子相互作用解析、蛋白質電気泳動ゲルおよびプロテオミクス電気泳動ゲルを含むグループから選択されてもよい。
【００４９】
構成要素は、システムの任意の測定可能な構成要素であることができる。生物学的システムの場合、構成要素は、たとえば、遺伝子またはその一部、ＤＮＡシーケンス、ＲＮＡシーケンス、ペプチド、蛋白質、炭水化物分子、脂質またはその混合物、生理学的構成要素、解剖学的構成要素、疫学的構成要素または化学的構成要素であってもよい。
【００５０】
訓練サンプルは、サンプルの特徴が知られているシステムから得られた任意のデータであることができる。たとえば、訓練サンプルは発生されたデータであってもよい。たとえば、訓練サンプルは、生物学的システムに適用されたサンプルから得られた任意のデータであることができる。たとえば、生物学的システムがＤＮＡマイクロアレイであるとき、訓練サンプルは、既知の特徴を有する細胞から抽出されたＲＮＡとのアレイのハイブリダイゼーションに後続してそのアレイから得られたデータ、または細胞から抽出されたＲＮＡから合成されたｃＤＮＡであることができ、あるいは生物学的システムがプロテオミクス電気泳動ゲルである場合には、訓練サンプルは、そのシステムに適用された蛋白質または細胞抽出から発生されてもよい。
【００５１】
本発明者は、本発明の方法がテスト治療に反応して混合した結果を示した被験者からのテストデータを再推定または推定するときに１実施形態において使用されてもよいと考えている。したがって、第２の特徴において、本発明は、被験者を複数の予め規定されたグループの１つに分類することのできる被験者の構成要素の部分集合を識別する方法を提供し、各グループは、
（ａ）複数の被験者をテスト治療に露出し、その治療に対する反応に基づいて被験者を反応グループにグループ化し、
（ｂ）被験者の構成要素を測定し、
（ｃ）統計的解析方法を使用して被験者を反応グループに分類することのできる構成要素の部分集合を識別するステップを含むテスト治療に対する反応によって規定される。  
統計的解析方法は、本発明の第１の特徴による方法であることが好ましい。
【００５２】
構成要素の部分集合が識別されると、その部分集合を使用して、テスト治療に反応する可能性の高いものおよび反応しないもののようなグループに被験者を分類することができる。このようにして、本発明の方法は、人口のほんの一部分に対して有効である可能性がある治療が識別されることを可能にし、また、テスト治療に対して反応するであろう人口のその一部分の識別を可能にする。
【００５３】
第３の特徴において、本発明は、各反応グループが複数の被験者をテスト治療に露出し、被験者をその治療に対する反応に基づいて反応グループにグループ化することにより形成された複数の予め規定された反応グループの１つに被験者を分類することのできる被験者の構成要素の部分集合を識別する装置を提供し、この装置は、
（ａ）被験者の測定された構成要素を受取る手段と、
（ｂ）統計的解析方法を使用して被験者を反応グループに分類することのできる構成要素の部分集合を識別する手段とを備えている。  
統計的解析方法は、本発明の第１または第２の特徴による方法であることが好ましい。  
第４の特徴において、本発明は、被験者をテスト化合物による治療に反応した者あるいは反応しなかった者として分類することのできる被験者の構成要素の部分集合を識別する方法を提供し、この方法は、
（ａ）複数の被験者を化合物に露出し、その化合物に対する各被験者の反応に基づいて被験者を反応グループにグループ化し、
（ｂ）被験者の構成要素を測定し、
（ｃ）統計的解析方法を使用して被験者を反応グループに分類することのできる構成要素の部分集合を識別するステップを含んでいる。  
統計的解析方法は、第１の特徴による方法であることが好ましい。
【００５４】
第５の特徴において、本発明は、各反応グループが複数の被験者を化合物に露出し、被験者をその化合物に対する反応に基づいて反応グループにグループ化することにより形成された複数の予め規定された反応グループの１つに被験者を分類することのできる被験者の構成要素の部分集合を識別する装置を提供し、この装置は、
（ｃ）被験者の測定された構成要素を受取る手段と、
（ｄ）統計的解析方法を使用して被験者を反応グループに分類することのできる構成要素の部分集合を識別する手段とを備えている。  
統計的解析方法は、本発明の第１または第２の特徴による方法であることが好ましい。  
本発明の第２乃至第５の特徴において測定された構成要素は、たとえば、その被験者の遺伝子または小さいヌクレオチド多型（ＳＮＰｓ）、蛋白質、抗体、炭水化物、脂質、あるいは任意の別の測定可能な構成要素であることができる。
【００５５】
とくに好ましい実施形態において、化合物は、薬剤化合物、または薬剤化合物と薬剤的に許容可能なキャリアとを含む合成物である。
【００５６】
本発明の識別方法は、適切なコンピュータソフトウェアおよびハードウェアによって実施されてもよい。
【００５７】
第６の特徴において、本発明は、テストサンプルの特徴を予測することのできるシステムから発生されたデータからのシステムの構成要素の部分集合を、そのシステムからの複数のサンプルから識別する装置を提供し、この装置は、
（ａ）各訓練サンプルが１つの既知の特徴を有している複数の訓練サンプルから発生されたデータから各構成要素に対する値が導かれる、構成要素と構成要素加重との線形結合を発生する手段と、
（ｂ）線形結合を条件とし、また、２クラス反応に対する二項分布を、線形結合とこの反応の期待値とをリンクするプロビット関数と組合せたものではない、ある特徴の確率分布に対するモデルを規定する手段と、
（ｃ）ゼロに近い高い確率密度を有するハイパープライアを含む線形結合の構成要素加重に対して事前分布を構成する手段と、
（ｄ）この事前分布とモデルとを組合せて事後分布を生成する手段と、
（ｅ）この事後分布を最大にする構成要素加重を有する構成要素の部分集合を識別する手段とを具備している。  
この装置は、適切にプログラムされた計算装置を備えていてもよい。
【００５８】
第７の特徴によると、本発明は、計算装置にロードされたときにこの計算装置を制御して、本発明の第１の特徴による方法を実施するように構成されているコンピュータプログラムを提供する。
【００５９】
コンピュータプログラムは、上述した本発明の第１または第２の特徴の好ましいアルゴリズムおよび方法ステップの任意のものを実施することができる。
【００６０】
本発明の第８の特徴によると、本発明の第４の特徴によるコンピュータプログラムを提供するコンピュータ読出し可能な媒体が提供される。
【００６１】
本発明の第９の特徴によると、システムからのサンプルをテストしてそのサンプルの特徴を識別する方法が提供され、この方法は、特徴の診断に役立つ構成要素の部分集合をテストするステップを含んでおり、その構成要素の部分集合は本発明の第１または第２の特徴による方法によって決定されている。  
システムは生物学的システムであることが好ましい。
【００６２】
本発明の第１０の特徴によると、サンプルの特徴を決定するためにシステムからのサンプルをテストする装置が提供され、この装置は本発明の第１または第２の特徴による方法にしたがって識別された構成要素をテストする手段を備えている。
【００６３】
本発明の第１１の特徴によると、本発明は、計算装置上で実行されたときに、システムからのテストサンプルの特徴を予測することのできるそのシステムからの構成要素を識別する方法でこの計算装置を制御するように構成されており、構成要素と構成要素加重との線形結合は、各訓練サンプルが１つの既知の特徴を有する複数の訓練サンプルから発生されたデータから発生され、ゼロに近い高い確率密度を有するハイパープライアを含む構成要素加重に対する事前分布と、線形結合を条件とするモデルとを組合せることによって事後分布が発生され、それによって事後分布を最大にする構成要素加重を推定し、前記モデルは２クラス反応に対する二項分布を、線形結合とこの反応の期待値とをリンクするプロビット関数と組合せたものではないコンピュータプログラムを提供する。
【００６４】
本発明の特徴が計算装置によって実施される場合、たとえば、ＰＣまたは主フレームあるいはネットワーク計算インフラストラクチャ等の、任意の適切なコンピュータハードウェアが使用されてもよいことが認識されるであろう。
【００６５】
第１２の特徴において、本発明は、生物学的システムからのテストサンプルの特徴を予測することのできる生物学的システムの構成要素の部分集合を識別する方法を提供し、この方法は、
（ａ）構成要素と構成要素加重との線形結合を発生し、この線形結合において、各訓練サンプルが１つの既知の特徴を有している複数の訓練サンプルから発生されたデータから各構成要素に対する値が決定され、
（ｂ）ある特徴の確率分布に対するモデルを規定し、このモデルは線形結合を条件とし、
（ｃ）ゼロに近い高い確率密度を有するハイパープライアを含む線形結合の構成要素加重に対して事前分布を構成し、
（ｄ）この事前分布とモデルとを組合せて事後分布を生成し、
この事後分布を最大にする構成要素加重を有する構成要素の部分集合を識別するステップを含んでいる。
【発明を実施するための最良の形態】
【００６６】
図１は、本発明の１実施形態の予測の成功に関する置換テストの結果を示している。クラスラベルは２００回ランダムに入れ替えられ、各置換に対して解析が繰返された。このヒストグラムは置換の下での予測の成功分布を示している。ｘ軸は正しく分類されたサンプルの数を示し、ｙ軸は頻度を示している。
【００６７】
図２は、本発明の１実施形態の予測の成功に関する置換テストの結果を示している。クラスラベルは２００回ランダムに入れ替えられ、各置換に対して解析が繰返された。このヒストグラムはクラスラベルの置換の下での予測の成功分布を示している。ｘ軸は合計サンプルの割合を示し、ｙ軸（ラムダ）は正しく分類されるケースの割合を示している。
【００６８】
図３は、本発明の方法の１実施形態において使用される一般化された線形モデルに対する曲線のプロットを示している。この方法により選択された５つの構成要素がこのモデル中で使用されたとき、あてはめられた曲線（実線）が生成され、真の曲線（点線）は破線として示されており、その５つの構成要素に基づく２００回の観察（ｘ軸）からのデータ（ｎｆ，ｙ軸）が円として示されている。
【００６９】
図４は、本発明の方法の１実施形態を使用して識別された単一の遺伝子に対するあてはめられた確率のプロットを示している。ｘ軸は遺伝子インデックスを示し、ｙ軸は特定の順序にされたクラスに属するサンプルの確率を示している。ラインは以下のようにクラスを表している：破線＝クラス１、実線＝クラス２、点線＝クラス３、点線と破線を組合せた線＝クラス４。
【００７０】
図５は、本発明によるシステムを構成するために使用されるパーソナルコンピュータの概略図である。
【００７１】
本発明は、特定の訓練サンプルが特定の特徴を有しているか否かを識別するために使用されることのできる、好ましくは最少数の、構成要素を識別する。この最少数の構成要素はその特徴の“診断に役立ち”、あるいは異なった特徴を有するサンプルを区別することを可能にする。本質的に、システムから発生されたデータの全てから、本発明の方法は、特定の特徴をテストするために使用されることのできる最少数の構成要素を識別することを可能にする。この方法によってこれらの構成要素が識別されると、その構成要素は新しいサンプルを評価するために将来使用されることができる。本発明の方法は統計的方法を使用して、特徴を正しく予測するのに必要とされない構成要素を除去する。
【００７２】
本発明者は、訓練サンプルから発生されたデータの構成要素の線形結合の構成要素加重が、訓練サンプルの特徴を正しく予測するのに必要とされない構成要素を消去する方法で推定可能であることを発見した。その結果、訓練サンプルの特徴を正しく予測することのできる構成要素の部分集合が識別される。このようにして、本発明の方法は、ある特徴を正しく予測することのできる比較的少数の構成要素を大量のデータから識別することを可能にする。
【００７３】
また、本発明の方法には、使用する必要のあるコンピュータメモリが、構成要素に関する限界情報ではなく共有（ｊｏｉｎｔ）情報を使用する従来技術の方法よりも少ないという利点がある。したがって、本発明の方法は、たとえば、ラップトップマシンのようなコンピュータ上で迅速に行われることができる。本発明の方法はまた、少ないメモリを使用することにより、たとえば、生物学的データの解析のために構成要素に関する共有情報（限界情報ではなく）を使用する従来技術の方法よりも迅速に行われることを可能にする。
【００７４】
以下、多クラスロジスティック回帰モデルに関する第１の実施形態を説明する。  
Ａ．多クラスロジスティック回帰モデル
この実施形態の方法は、訓練サンプルを予め規定されたグループに分類することのできる構成要素の部分集合を識別するために訓練サンプルを使用する。続いて、たとえば、臨床テストのようなテストに対して構成要素の部分集合の情報を使用して、サンプルに病気クラスのようなグループに分類することができる。たとえば、ＤＮＡマイクロアレイの構成要素の部分集合は、臨床サンプルを、たとえば、健康または病気のような臨床的に適切なクラスにグループ化するために使用されることができる。
【００７５】
このようにして、本発明は、ある特定の訓練サンプルがある特定のグループに属しているか否かを識別するために使用されることのできる、好ましくは最少数の構成要素を識別する。最少数の構成要素はそのグループの“診断に役立ち”、あるいはグループの区別を可能にする。本質的に、システムから発生されたデータの全てから、本発明の方法は、特定のグループをテストするために使用されることのできる最少数の構成要素を識別することを可能にする。この方法によってこれらの構成要素が識別されると、その構成要素は新しいサンプルをグループに分類するために将来使用されることができる。本発明の方法は統計的方法を使用して、サンプルが属しているグループを正しく識別するのに必要とされない構成要素を除去することが好ましい。
【００７６】
サンプルは、予め定められた分類に基づいてサンプルグループ（または“クラス”）にグループ化される。この分類は、訓練サンプルがグループ化される任意の所望の分類であってよい。この分類は、たとえば、訓練サンプルが得られたのは白血病細胞または健康な細胞のどちらからかでもよく、あるいは訓練サンプルが得られたのはある症状を有する患者または有しない患者のどちらの血液からかでもよく、もしくは訓練サンプルが得られたのは正常な細胞と比較されるいくつかのタイプの癌の１つからの細胞からであってもよい。
【００７７】
１実施形態において、入力データは、ｎ個の訓練サンプルおよびｐ個の構成要素を有するｎ×ｐデータマトリックスＸ＝（ｘij）に編成される。典型的に、ｐはｎよりはるかに大きくなる。
【００７８】
別の実施形態において、データマトリックスＸは、線形予測子の代りの予測子としてＸの平滑な関数を得るためにｎ×ｎカーネルマトリックスＫによって置換されてもよい。カーネルマトリックスＫの一例はｋij＝ｅｘｐ（−０．５＊（ｘi−ｘj）t（ｘi−ｘj）／σ2）であり、ここでｘの添字はマトリックスＸ中の行の数を示している。理想的には、これらの平滑な関数の希薄な（ｓｐａｒｓｅ）表現を提供するＫの列の部分集合が選択される。カーネルマトリックスＫの別の例は以下の表２に示されている。（表３は少しでも必要か。）
各サンプルクラス（グループ）にはあるクラスラベルｙiが関連付けられており、ここでｙi＝ｋ，ｋ∈｛１，・・・，Ｇ｝であり、これは、訓練サンプルがＧのサンプルクラスのいずれに属しているかを示す。エレメントｙiを有するｎ×ｌベクトルをｙと書く。そのベクトル：
【数３０】
2005524124.tif 000017
【００７９】
を仮定すると、インジケータ変数が規定されることができる：
ｙi＝ｇならば、ｅig＝１
その他場合、ｅig＝０                          （１Ａ）
【００８０】
１実施形態において、構成要素加重は、ベイズの統計的モデル（KotzJohnson氏等による1983年の文献を参照されたい）を使用して推定される。この加重は、各訓練サンプルから発生されたデータを与えられた加重の事後分布を最大化することによって推定されることが好ましい。その結果、最大化されるべき目標関数は２つの部分から構成されたものとなる。第１の部分は尤度関数であり、第２の部分は、可能な限り常にゼロ加重が好ましいことを保証する加重に対する事前分布である。好ましい実施形態において、尤度関数は多クラスロジスティックモデルから導出される。尤度関数は、確率から計算されることが好ましい：
【数３１】
2005524124.tif 000018
【００８１】
ここで、
ｐigは、入力データＸiを有する訓練サンプルがサンプルクラスｇ中に存在する確率であり、
ｘiTβgは、構成要素加重βgを有する訓練サンプルｉからの入力データから発生された線形結合であり、
ｘiTはＸのｉ番目の行に対する構成要素であり、βgはサンプルクラスｇに対する構成要素加重のセットである。
【００８２】
上述したように、構成要素加重は典型的に、構成要素加重の大部分がゼロであるアプリオリ仮定を考慮した方式で推定される。
【００８３】
１実施形態において、式（２Ａ）中の構成要素加重βgは、値の大部分がゼロであるが、そのサンプルが依然として正確に分類されることのできる方式で推定される。
【００８４】
１実施形態において、パラメータβ1，・・・，βG-1に対して特定されたプライアは以下の形態のものである：
【数３２】
2005524124.tif 000019
【００８５】
また、ｐ（βg|τg2）はＮ（０，ｄｉａｇ｛τg2｝）であり、
【数３３】
2005524124.tif 000020
【００８６】
はジェフェリーのハイパープライアである（KotzおよびJohnson氏等の文献、1983）。
【００８７】
１実施形態において、尤度関数は式（８Ａ）中の形の
【数３４】
2005524124.tif 000021
【００８８】
であり、ｙを与えられたβおよび
【数３５】
2005524124.tif 000022
【００８９】
の事後分布は、
【数３６】
2005524124.tif 000023
【００９０】
１実施形態において、第１の導関数は以下の式から決定される：
【数３７】
2005524124.tif 000024
【００９１】
ここで、
【数３８】
2005524124.tif 000025
【００９２】
は、サンプルクラスｇのメンバーシップおよびクラスｇの確率をそれぞれ示すベクトルである。
【００９３】
１実施形態において、第２の導関数は以下のアルゴリズムから決定される：
【数３９】
2005524124.tif 000026
【００９４】
式６および式７は、次のように導出されてもよい：
（ａ）式（１Ａ），（２Ａ）および（３Ａ）を使用することにより、データの尤度関数を以下のように書くことができる：
【数４０】
2005524124.tif 000027
【００９５】
（ｂ）式（８Ａ）の対数をとり、全てのｉに対して
【数４１】
2005524124.tif 000028
【００９６】
であることを使用することによって、
【数４２】
2005524124.tif 000029
【００９７】
が得られる。
【００９８】
（ｃ）式（９Ａ）をβgに関して微分することにより、次式が得られる：
【数４３】
2005524124.tif 000030
【００９９】
ここで、
【数４４】
2005524124.tif 000031
【０１００】
は、サンプルクラスｇのメンバーシップおよびクラスｇの確率をそれぞれ示すベクトルである。
【０１０１】
（ｄ）式（９Ａ）の第２の微分はエレメント：
【数４５】
2005524124.tif 000032
【０１０２】
を有しており、ここで、
ｈ＝ｇならば、δhg＝１
そうでない場合、δhg＝０である。
【０１０３】
尤度関数の事後分布を最大化する構成要素加重は、ＥステップおよびＭステップを含むＥＭアルゴリズムを使用して特定されてもよい。
【０１０４】
典型的に、ＥＭアルゴリズムは：
（ａ）関数：
【数４６】
2005524124.tif 000033
【０１０５】
を使用して構成要素加重の事後分布の条件付き期待値を計算することによってＥステップを行い、
（ｂ）γの関数としてＱを最大化するために反復的な処理手続きを適用することによりＭステップを行い、このとき、
【数４７】
2005524124.tif 000034
【０１０６】
この場合、α’はステップ長であり、このとき０≦α’≦１である；
βg＝Ｐgγg；
ここで、Ｐgはゼロと１との、ＰTgβgがβgのゼロでないエレメントを選択するようなマトリックスである；
γ＝（γg，ｇ＝１，・・・，Ｇ−１）
【０１０７】
式（１２Ａ）は以下のように導出されることができる：
パラメータ推定値のセット：
【数４８】
2005524124.tif 000035
【０１０８】
と観察されたデータｙとを与えられた式（５Ａ）の条件付き期待値を計算する。  
【数４９】
2005524124.tif 000036
【０１０９】
β（および
【数５０】
2005524124.tif 000037
【０１１０】
）の構成要素がゼロに設定される場合、すなわち、ｇ＝１，・・・，Ｇ−１、βg＝Ｐgγgおよび
【数５１】
2005524124.tif 000038
【０１１１】
に対する場合を検討する。ここで、Ｐgはゼロと１との、ＰTgβgがβgのゼロでないエレメントを選択するようなマトリックスである。以下において、γ＝（γg，ｇ＝１，・・・，Ｇ−１）と書く。γgは実際にβgの構成要素の部分集合であることに注意されたい。それらを使用して、表記を可能な限り単純なものにしておく。
【０１１２】
γを含まない項を無視し、式（４Ａ），（５Ａ），（９Ａ）を使用することにより、次式が得られる：
【数５２】
2005524124.tif 000039
【０１１３】
条件付き期待値は、式（４Ａ）を与えられた第１の原理から推定されることができることに注意されたい。
【０１１４】
反復的な処理手続きは以下のように得られることができる：
式（１３Ａ）において必要とされる導関数を得るために、最初に、式（８Ａ）、（９Ａ）および（１０Ａ）から以下の式を得ることに注意されたい：
【数５３】
2005524124.tif 000040
【０１１５】
ここで、
ｇ＝ｈの場合、   δgh＝１
それ以外の場合、δgh＝０
および
ＸgT＝ＰgTＸT  ， ｇ＝１，・・・，Ｇ−１          （１７Ａ）
【０１１６】
好ましい実施形態において、反復的な処理手続きは、式（１３Ａ）において式（１６Ａ）のブロック対角線上の成分（ｄｉａｇｏｎａｌｓ）だけを使用することによって簡単化されることができる。これによって、ｇ＝１，・・・，Ｇ−１に対して、以下の関係が得られる：
【数５４】
2005524124.tif 000041
【０１１７】
式（１８Ａ）を再構成することにより、次式が得られる：
【数５５】
2005524124.tif 000042
【０１１８】
Ｙgの列の数に対してｐ（ｇ）を書くと、式（１９Ａ）は、非常に大きいものとなる可能性のあるｐ（ｇ）×ｐ（ｇ）マトリックスの逆数を必要とする。これは、ｐ（ｇ）＞ｎに対しては、以下のように表すことによりｎ×ｎマトリックスに減少されることができる：
【数５６】
2005524124.tif 000043
【０１１９】
ここで、Ｚg＝Δ1/2ggＹgである。ｐ（ｇ）＜ｎのときは、式（１９Ａ）が使用され、ｐ（ｇ）≧ｎであるときには、式（１９Ａ）に式（２０Ａ）を代入したものが使用されることが好ましい。
【０１２０】
好ましい実施形態において、ＥＭアルゴリズムは以下のように行われる：
１．ｎ＝０、Ｐg＝Ｉと設定し、
【数５７】
2005524124.tif 000044
【０１２１】
に対する初期値を選択する。これは、ｘiに関するｌｏｇ（ｐig／ｐiG）のリッジ（ｒｉｄｇｅ）回帰により行われ、ここでｐigは、全ての確率が合計すると１になるという制約を条件として、グループｇにおいて観察されたものに対して１に近くなり、その他の場合は、わずかな量＞０となるように選択される。
【０１２２】
２．Ｅステップを行う。すなわち、
【数５８】
2005524124.tif 000045
【０１２３】
を推定する。
【０１２４】
３．ｔ＝０を設定する。ｇ＝１，・・・，Ｇ−１に対して、以下を計算する：
（ａ）ｐ（ｇ）≧ｎであるとき、式（１９Ａ）に式（２０Ａ）を代入したものを使用して、δtg＝γt+1g−γtg
（ｂ）δt＝（δtg，ｇ＝１，・・・，Ｇ−１）と書く。式（１２Ａ）をαtの関数として最大化する（あるいは単に増加させる）
【数５９】
2005524124.tif 000046
【０１２５】
中のαtの値を見つけるためにラインサーチを行う。  
（ｃ）
【数６０】
2005524124.tif 000047
【０１２６】
およびｔ＝ｔ＋１を設定する。  
収斂するまでステップ（ａ）乃至（ｃ）を繰り返す。  
これによって、たとえば、現在のＱ関数をγの関数として最大化するγ*n+1が生成される。  
ｇ＝１，・・・，Ｇ−１に対して、
【数６１】
2005524124.tif 000048
【０１２７】
を決定する。
【０１２８】
ここで、ε＜＜１、たとえば１０-5である。ｉ∈Ｓgおよび
【数６２】
2005524124.tif 000049
【０１２９】
に対してΒig＝０となるようにＰgを定める。  
このステップにより、小さい係数を有する変数がモデルから除去される。
【０１３０】
４．ｎ＝ｎ＋１を設定し、収斂するまで２を行う。
【０１３１】
以下、順序付けられたカテゴリー的ロジスティック回帰に関する第２の実施形態を説明する。  
Ｂ．順序付けられたカテゴリー的モデル
この実施形態の方法は訓練サンプルを使用して、あるテストサンプルが特定のクラスに属しているか否かを決定するために使用されることのできる構成要素の部分集合を識別することができる。たとえば、マイクロアレイ解析を使用して組織生検サンプルを評価するために遺伝子を識別するために、正常な組織、良性の組織、局所化した腫瘍および転移した腫瘍組織のような、前に病気の重大度が増加するまたは減少するクラスに順序付けられた組織からの一連のサンプルからのマイクロアレイデータは、訓練サンプルに関連した病気の重大度を示すことのできる構成要素の部分集合を識別するために訓練サンプルとして使用される。その後、構成要素の部分集合は、以前分類されていなかったテストサンプルが正常、良性、局所化した腫瘍、または転移した腫瘍組織のいずれとして分類されるかを決定するために後続的に使用されることができる。したがって、構成要素の部分集合は、テストサンプルが順序付けられたクラスセット内の特定のクラスに属するかどうかの診断に役立つ。構成要素の部分集合がひとたび識別されると、将来の診断処理手続きにおいて、あるサンプルが順序付けられたクラスのどれに属しているのかを決定するのに、その構成要素の部分集合がテストされるだけでよいことが認識される。
【０１３２】
本発明の方法はとくに、非常に大量のデータの解析に適している。典型的に、テストサンプルから得られた大量のデータセットは非常に多種多様であり、また、訓練サンプルから得られたものとは著しく異なっていることが多い。本発明の方法は訓練サンプルから発生された非常に大量のデータからの構成要素の部分集合を識別することが可能であり、また、この方法によって識別された構成要素の部分集合は、その後、テストサンプルから発生されたデータが同じクラスに属している訓練サンプルから発生されたデータと比較して大幅に異なっている場合でも、そのテストサンプルを分類するために使用されることができる。したがって、本発明の方法は、データが品質の低いものであっても、および、または順序付けられた同じクラスのサンプル間に大きいばらつきがある場合であってさえも、サンプルを正しく分類する可能性が高い構成要素の部分集合を識別することができる。
【０１３３】
最少数の構成要素は、その特定の順序付けられたクラスに対する“予測的な”ものである。本質的に、システムから発生された全てのデータから、本発明の方法は訓練データを分類するために使用されることのできる最少数の構成要素を識別することができる。この方法によってこれらの構成要素が識別されると、その構成要素は、将来、テストサンプルを分類するために使用されることができる。本発明の方法は統計的方法を使用して、サンプルを順序付けられたクラスのメンバーであるクラスに正しく分類するのに必要でない構成要素を除去することが好ましい。
【０１３４】
以下において、Ｎ個のサンプルが存在し、ｙ、ｚおよびμのようなベクトルはｉ＝１，・・・，Ｎに対して構成要素ｙi、ｚiおよびμiを有している。ベクトル乗算および除算が構成要素に関して定義され、ｄｉａｇ｛・｝は、対角線上の成分がアーギュメントに等しい対角行列を表す。ユークリッドのノルムを表すために｜｜・｜｜もまた使用される。
【０１３５】
Ｎ個の観察ｙiが存在していることが好ましく、ここでｙiは整数値１，・・・，Ｇである。これらの値は、たとえば、病気の重大度のようなある方式で順序付けられたクラスを表す。Ｎ行およびｐ列を有するマトリックスＸに配列されたあるセットのコバリエイト（たとえば、遺伝子発現値等の変数）が各観察と関連付けられており、ここでＮはサンプルであり、ｐは構成要素である。表記ｘiTはＸのｉ番目の行を表す。個体（サンプル）ｉは、πik＝πk（ｘi）で与えられるクラスｋに属する確率を有している。
【０１３６】
累積確率を規定する：
【数６３】
2005524124.tif 000050
【０１３７】
γikは観察ｉがｋ以下のインデックスを有するクラスに属する確率であることに注意されたい。Ｃを、
観察ｉがクラスｊ中のものである場合、ｃij＝１
そうでない場合、ｃij＝０
により与えられるエレメントｃijを有するｎ×ｐマトリックスとし、Ｒを
【数６４】
2005524124.tif 000051
【０１３８】
により与えられるエレメントｒijを有するｎ×Ｐマトリックスとする。  
これらは、行内におけるＣの列の累積的な和である。
【０１３９】
独立した観察（サンプル）に関して、データの尤度は以下のように書かれることができる：
【数６５】
2005524124.tif 000052
【０１４０】
また、対数尤度（ｌｏｇ（Ｌ））ｌは以下のように書かれることができる：
【数６６】
2005524124.tif 000053
【０１４１】
ここにおいて、以下のように、ｋ＝２，・・・，Ｇに対して、連続率（ｃｏｎｔｉｎｕａｔｉｏｎ  ｒａｔｉｏ）モデルが採用されてもよい：
【数６７】
2005524124.tif 000054
【０１４２】
McCullaghおよびNelder氏等の文献(1989)、McCullagh氏の文献(1980)ならびにこの明細書における説明を参照されたい。
【数６８】
2005524124.tif 000055
【０１４３】
であることに注意されたい。
【０１４４】
尤度は、反応ベクトルｙおよびコバリエイトマトリックスＸを有する対数回帰尤度に等しい
【数６９】
2005524124.tif 000056
【０１４５】
ここで、ＩG-1は（Ｇ−１）×（Ｇ−１）識別子マトリックスであり、ｌG-1は１の（Ｇ−１）×ｌベクトルである。  
ここにおいて、ｖｅｃ｛  ｝はマトリックスをとり、ベクトル行×行を成す。
【０１４６】
典型的に、上述したように、構成要素加重は、構成要素加重の大部分がゼロであるというアプリオリ仮定を考慮したやり方で推定される。
【０１４７】
Figueiredo氏の文献(2001)にならって、冗長な変数（コバリエイト）を除去するために、ハイパーパラメータ（ｈｙｐｅｒｐａｒａｍｅｔｅｒｓ）のｐ×ｌベクトルを導入することにより、あるプライアがパラメータβ*に対して特定される。
【０１４８】
構成要素加重に対して特定されたプライアは、以下の形のものであることが好ましい：
【数７０】
2005524124.tif 000057
【０１４９】
ここで、KotzおよびJohnson氏等の文献(1983)に示されているように、ｐ（β*｜ｖ2）はＮ（０，ｄｉａｇ｛ｖ2｝）であり、
【数７１】
2005524124.tif 000058
【０１５０】
はジェフェリーのプライアである。θ＝（θ2，・・・θG）Tのエレメントは非有益プライアを有している。
【０１５１】
尤度関数を
【数７２】
2005524124.tif 000059
【０１５２】
と書くことにより、ベイズのフレームワークにおいて、ｙを与えられたβ*、θおよびｖの事後分布は、以下のようになる：
【数７３】
2005524124.tif 000060
【０１５３】
Ｖを失われたデータのベクトルとして処理することにより、式（６Ｂ）を最大化して局所極大であるβ*およびθの事後推定値を生成するために、ＥＭアルゴリズム［Dempster氏の文献(1977)］のような反復的アルゴリズムが使用可能であることが好ましい。上記のプライアは最大である事後推定値が希薄（ｓｐａｒｓｅ）になる傾向がある、すなわち、多数のパラメータが冗長である場合、β*の多くの構成要素がゼロになるようなものである。以下においてβT＝（θT，β*T）であり、ｄｉａｇ（）は対角行列を表すことが好ましい：
上記の順序付けられたカテゴリーに関して、
【数７４】
2005524124.tif 000061
【０１５４】
であることが証明されることができる。
【０１５５】
上述したように、事後分布を最大化する構成要素加重は、反復的な処理手続きを使用して決定されることができる。構成要素および構成要素加重の事後分布を最大化する反復的な処理手続きは、たとえば、Dempster氏の文献(1977)に記載されているようなＥＭアルゴリズムであることが好ましい。ＥＭアルゴリズムは以下のようにして行われることが好ましい：
１．ｎ＝０，Ｓ0＝｛１，２，・・・，ｐ｝，φ(0) および（たとえば）ε＝１０-5を設定する。レギュラリゼーション（regularization）パラメータκを１よりはるかに大きい値、たとえば、１００に設定する。これは、１／κ2を第２の微分マトリックスの第１のＧ−１対角エレメントに以下のＭステップで加算することに対応する。
【０１５６】
ｐ≦Ｎである場合、
【数７５】
2005524124.tif 000062
【０１５７】
によって初期値β*を計算し、
ｐ＞Ｎである場合、
【数７６】
2005524124.tif 000063
【０１５８】
によって初期値β*を計算する。ここでリッジパラメータλは、０＜λ≦１を満足させ、ζは小さく、ロジット（ｌｏｇｉｔ）リンク関数ｇがｙ＋ζにおいて明確に定義されるように選択される。
【０１５９】
２．
【数７７】
2005524124.tif 000064
【０１６０】
と定義し、Ｐnを、β(n)のゼロでないエレメントγ(n)が
γ(n)＝ＰnTβ(n) ，    β(n)＝Ｐnγ(n)
γ＝ＰnTβ  ，      β＝Ｐnγ
を満足させるような０と１とのマトリックスとする。
【０１６１】
ｉ≧Ｇである場合、ｗβi＝１
それ以外の場合、ｗβi＝０
であるようにｗβ＝（ｗβi，ｉ＝１，ｐ）と定義し、ｗγ＝Ｐnｗβとする。
【０１６２】
３．
【数７８】
2005524124.tif 000065
【０１６３】
を計算することによってＥステップを行う。ここで、ｌはｙの対数尤度関数である。
【０１６４】
β＝Ｐnγおよびβ(n)＝Ｐnγ(n)を使用することにより、式（１１Ｂ）は以下のように書かれることができる：
【数７９】
2005524124.tif 000066
【０１６５】
４．Ｍステップを行う。これは、以下のようにＮｅｗｔｏｎ  Ｒａｐｈｓｏｎ反復により行われることができる。γ0＝γ(n) を設定し、ｒ＝０，１，２，・・・に対してγr+1＝γr＋αrδrであり、ここでαrは、
Ｑ（γr+1｜γ(n)）＞Ｑ（γr｜γ(n)）
を保証するようにラインサーチアルゴリズムによって選択される。
【０１６６】
ｐ≦Ｎに対しては、以下の式を使用する：
【数８０】
2005524124.tif 000067
【０１６７】
ｐ＞Ｎに対しては、以下の式を使用する：
【数８１】
2005524124.tif 000068
【０１６８】
Ｖrおよびｚrは前のように定義されている。
【０１６９】
γ*を、ある収斂基準が満足されたとき、たとえば、以下のようなときのγrの値とする：
｜｜γr−γr+1｜｜＜ε（たとえば、１０-5）
５．
【数８２】
2005524124.tif 000069
【０１７０】
と定義する。ここで、ε1は、たとえば、ｌｅ−５等の小さい定数である。ｎ＝ｎ＋１を設定する。
【０１７１】
６．収斂をチェックする。｜｜γ*−γ(n)｜｜＜ε2である場合、中止し、ここでε2は適切に小さく、その他の場合、上記のステップ２に進む。
【０１７２】
確率を再生する。
【０１７３】
パラメータβの推定値が得られると、ｉ＝１，・・・，Ｎおよびｋ＝２，・・・，Ｇに対して、
【数８３】
2005524124.tif 000070
【０１７４】
を計算する。
【０１７５】
確率を得るために、以下の反復が使用されることが好ましい：
【数８４】
2005524124.tif 000071
【０１７６】
また、ｉ＝１，・・・，Ｎに対して、この確率の合計は１になる。
【０１７７】
１実施形態において、行ｘiTを有するコバリエイトマトリックスＸは、ｉｊ番目のエレメントｋijを有し、あるカーネル関数κに対してｋij＝κ（ｘi−ｘj）であるマトリックスＫによって置換されることができる。このマトリックスはまた１のベクトルにより増加されることができる。ある例のカーネルが以下の表１に与えられている［Evgeniou氏他による文献(1999)参照］。
【０１７８】
カーネル関数                          κ（ｘ−ｙ）に対する公式
ガウスのラジアル基礎関数              ｅｘｐ（−｜｜ｘ−ｙ｜｜2／ａ），ａ＞０
逆マルチクアドリック(multiquadric)      （｜｜ｘ−ｙ｜｜2＋ｃ2）-1/2
マルチクアドリック                  （｜｜ｘ−ｙ｜｜2＋ｃ2）1/2
薄板スプライン                    ｜｜ｘ−ｙ｜｜2n+1
｜｜ｘ−ｙ｜｜2nｌｎ（｜｜ｘ−ｙ｜｜）
マルチレーヤ知覚                  適切なθに対して、ｔａｎｈ（ｘ’ｙ−θ）
度数ｄの多項式                    （１＋ｘ’ｙ）d
Ｂスプライン                      Ｂ2n+1（ｘ−ｙ）
三角多項式                 ｓｉｎ（（ｄ＋１／２）（ｘ−ｙ））/ｓｉｎ（（ｘ−
ｙ）／２）
表１：カーネル関数の例
【０１７９】
表１において、最後の２つのカーネルは１次元、すなわち、Ｘが列を１つしか有しない場合に対して１次元であることが好ましい。多変量バージョンはこれらのカーネル関数の積から得られることができる。Ｂ2n+1の定義はDe Boor氏による文献（1978）に記載されている。カーネル関数を使用した結果、推定された確率はコバリエイトＸの平滑な（線形の変換とは対照的に）関数となる。このようなモデルにより、データに対する実質的に良好な適合が得られる。
【０１８０】
以下、一般化された線形モデルに関する第３の実施形態を説明する。  
Ｃ．一般化された線形モデル
この実施形態の方法は、あるサンプルの特徴を予測することのできる構成要素の部分集合を識別するために訓練サンプルを使用する。後続的に、構成要素の部分集合の知識は、重要な特徴の未知の値を予測するための、たとえば、臨床テスト等の、テストのために使用されることができる。たとえば、ＤＮＡマイクロアレイの構成要素の部分集合は、たとえば、血液グルコースレベル、白血球細胞数、腫瘍の寸法、腫瘍の成長速度または生存時間等の、臨床的に関連した特徴を予測するために使用されることができる。
【０１８１】
このようにして、本発明は、特定のサンプルに対する特徴を予測するために使用されることのできる最少数の構成要素を識別することが好ましい。最少数の構成要素が、その特徴に対する“予測的な”ものである。本質的に、システムから発生された全てのデータから、本発明の方法は、特定の特徴を予測するために使用されることのできる最少数の構成要素を識別することができる。この方法によってこれらの構成要素が識別されると、その構成要素は、将来、新しいサンプルに対する特徴を予測するために使用されることができる。本発明の方法は統計的方法を使用して、サンプルに対する特徴を正しく予測するのに必要でない構成要素を除去することが好ましい。
【０１８２】
本発明は、訓練サンプルから発生されたデータの構成要素の線形結合の構成要素加重が、ある訓練サンプルに対する特徴を予測するのに必要とされない構成要素を除去する方式で推定されることができることを見出している。その結果、訓練セット中のサンプルに対する特徴を正しく予測することのできる構成要素の部分集合が識別される。したがって、本発明の方法は、たとえば、重要な量等の、訓練サンプルに対する特徴を正しく予測することのできる比較的少数の構成要素を大量のデータから識別することを可能にする。
【０１８３】
特徴は任意の重要な特徴であることができる。１実施形態において、特徴は量または測定値である。別の実施形態において、それらはあるグループのインデックス番号であることができ、ここにおいて、サンプルは予め定められた分類に基づいて２つのサンプルグループ（または“クラス”）にグループ化されている。この分類は、訓練サンプルがグループ化される任意の所望の分類であることができる。たとえば、分類は、訓練サンプルが白血病細胞からのものか、あるいは健康な細胞からのものかであってもよく、もしくは訓練サンプルが得られたのはある症状を有する患者または有しない患者かでもよく、あるいは訓練サンプルが得られたのは正常な細胞と比較されるいくつかのタイプの癌の１つからの細胞からであってもよい。別の実施形態において、特徴は、特定の患者が少なくとも所定の日数生存していることを示す検定された生存時間であってもよい。別の実施形態において、量は、たとえば、血圧等の、測定可能なサンプルの連続的に変化する任意の特徴であってもよい。
【０１８４】
１実施形態において、データは量ｙiであってもよく、ここでｉ∈｛１，・・・，Ｎ｝。エレメントｙiを有するＮ×１ベクトルをｙと書く。構成要素加重のｐ×１パラメータベクトルβ（それらの多くはゼロであると期待される）およびパラメータφのｑ×１ベクトル（ゼロでないと期待される）が規定される。ｑはゼロであることができる（すなわち、ゼロでないと期待されたパラメータのセットは空であることができる）ことに注意されたい。
【０１８５】
１実施形態において、入力データはＮ個のテスト訓練サンプルとｐ個の構成要素とを有するＮ×ｐデータマトリックスＸ＝（ｘij）に編成される。典型的に、ｐはＮよりはるかに大きい。
【０１８６】
別の実施形態において、データマトリックスＸは、線形予測子の代りの予測子としてＸの平滑な関数を得るためにＮ×ＮカーネルマトリックスＫによって置換されることができる。カーネルマトリックスＫの一例はｋij＝ｅｘｐ（−０．５＊（ｘi−ｘj）t（ｘi−ｘj）／σ2）であり、ここでｘの添字はマトリックスＸ中の行番号を示す。理想的には、これらの平滑な関数のスパース表現を提供するＫの構成要素の部分集合が選択される。
【０１８７】
典型的に、上述したように、構成要素加重は、構成要素加重の大部分がゼロであるというアプリオリ仮定を考慮した方式で推定される。
【０１８８】
１実施形態において、構成要素加重に対して特定されたプライアは以下の形態のものである：
【数８５】
2005524124.tif 000072
【０１８９】
また、ｐ（β|ｖ2）はＮ（０，ｄｉａｇ｛ｖ2｝）であり、
【数８６】
2005524124.tif 000073
【０１９０】
はジェフェリーのプライアである（KotzおよびJohnson氏等の文献、1983）。φに対する非有益プライアが特定されることが好ましい。
【０１９１】
尤度関数は、データの分布に基づいてそのデータに適合するモデルを規定する。尤度関数は、一般化された線形モデルから得られることが好ましい。たとえば、尤度関数：
【数８７】
2005524124.tif 000074
【０１９２】
は、たとえば、ＮｅｌｄｅｒおよびＷｅｄｄｅｒｂｕｒｎ氏等の文献(1972)に記載されているもののような、一般化された線形モデル（ＧＬＭ）に対して適した形態であることができる。尤度関数は、以下の形態のものであることが好ましい：
【数８８】
2005524124.tif 000075
【０１９３】
ここで、ｙ＝（ｙ1，・・・，ｙn）Tおよびａi（φ）＝φ／ｗiであり、ｗiは既知の加重の固定されたセットであり、φは単一のスケールパラメータである。
【０１９４】
尤度関数は、以下のように特定されることが好ましい：
【数８９】
2005524124.tif 000076
【０１９５】
各観察はコバリエイトのセットｘiおよび線形予測子ηi＝ｘiTβを有している。ｉ番目の観察の平均とそれの線形予測子との間の関係は、リンク関数ηi＝ｇ（μi）＝ｇ（ｂ’（θi））によって与えられる。このリンクの逆数はｈによって表される。すなわち、
μi＝＝ｂ’（θi）＝ｈ（ηi）
である。
【０１９６】
スケールパラメータに加えて、一般化された線形モデルは４つの構成要素によって特定されることができる：
・尤度または（スケールされた）逸脱関数
・リンク関数
・リンク関数の微分
・分散関数。
【０１９７】
一般化された線形モデルのいくつかの共通した例は、以下の表２に与えられている。
【表１】
2005524124.tif 000077
【０１９８】
別の実施形態において、尤度関数は多クラスロジスティックモデルから得られる。
【０１９９】
別の実施形態において、リンク関数および分散関数だけが規定される擬似尤度モデルが特定される。いくつかの例において、このような仕様は結果的に上記の表中のモデルを生じさせる。別の例では、分散は特定されない。
【０２００】
１実施形態において、ｙを与えられたβφおよびｖの事後分布は、以下の式を使用して推定される：
【数９０】
2005524124.tif 000078
【０２０１】
１実施形態において、ｖは失われたデータのベクトルとして処理されてもよく、反復的な処理手続きは局所極大であるβの事後推定値を生成するように式（２Ｃ）を最大化するために使用されてもよい。式（５Ｃ）のプライアは、最大である事後推定値が希薄になるようなもの、すなわち、非常に多くのパラメータが冗長である場合にβの多くの構成要素がゼロになるようなものである。
【０２０２】
上述したように、事後分布を最大化する構成要素加重は、反復的処理手続きを使用して決定されることができる。構成要素および構成要素加重の事後分布を最大化する反復的処理手続きは、たとえば、Dempster氏他の文献(1977)に記載されているようなＥＭアルゴリズムである。
【０２０３】
１実施形態において、ＥＭアルゴリズムは以下のステップを含んでいる：
（ｃ）ｎ＝０，Ｓ0＝｛１，２，・・・，ｐ｝を設定することによりアルゴリズムを初期化し、φ(0) ，β*を初期化し、たとえばε＝１０-5のようなεに対する値を適用する；
（ｄ）ｉεＳnのときは、βi(n) ＝βi*，
それ以外のときは、βi(n) ＝０                （５Ｃ）
と定義し、ｐnを、β(n)のゼロでないエレメントγ(n)が
γ(n)＝ＰnTβ(n) ，    β(n)＝Ｐnγ(n)
γ＝ＰnTβ  ，      β＝Ｐnγ
を満足させるような０と１とのマトリックスとする。  
（ｅ）以下の関数を使用して構成要素加重の事後分布の条件付き期待値を計算することにより推定（Ｅ）ステップを行う：
【数９１】
2005524124.tif 000079
【０２０４】
ここで、ｌはｙの対数尤度関数である。
【０２０５】
β＝Ｐnγおよびβ(n)＝Ｐnγ(n)を使用することにより、この式（６Ｃ）は以下のように書かれることができる：
【数９２】
2005524124.tif 000080
【０２０６】
（ｆ）反復的な処理手続きを適用してＱをγの関数として最大化することにより、最大化（Ｍ）ステップを行う。このとき、γ0＝γ(n) であり、ｒ＝０，１，２，・・・に対してγr+1＝γr＋αrδrであり、ここでαrは、
Ｑ（γr+1｜γ(n)，φ(n)）＞Ｑ（γr｜γ(n)，φ(n)）
を保証するようにラインサーチアルゴリズムによって選択され、また、
【数９３】
2005524124.tif 000081
【０２０７】
（ｇ）γ*を、ある収斂基準が満足されたとき、たとえば、以下のようなときのγrの値とする：
｜｜γr−γr+1｜｜＜ε（たとえば、１０-5）；
（ｈ）
【数９４】
2005524124.tif 000082
【０２０８】
と定義する。ここで、ε1は、たとえば、ｌｅ−５等の小さい定数である。  
（ｉ）ｎ＝ｎ＋１を設定し、φ(n+1)＝φ(n)＋κn（φ*−φ(n)）を選択し、ここでφ*は、
【数９５】
2005524124.tif 000083
【０２０９】
を満足させ、κnは０＜κn≦１となるような減衰定数である；
（ｊ）収斂をチェックする。｜｜γ*−γ(n)｜｜＜ε2である場合、ε2は適切に小さく、停止し、その他の場合、上記のステップ（ｂ）に進む。
【０２１０】
別の実施形態において、最大化ステップにおけるステップ（ｄ）は、
【数９６】
2005524124.tif 000084
【０２１１】
をその期待値：
【数９７】
2005524124.tif 000085
【０２１２】
で置換することにより推定されることができる。これは、データのモデルが一般化された線形モデルであるときに好ましい。
【０２１３】
一般化された線形モデルについて、期待値：
【数９８】
2005524124.tif 000086
【０２１４】
は以下のように計算されることができる：
【数９９】
2005524124.tif 000087
【０２１５】
ここで、ＸはＮ×ｐマトリックスであり、ｉ番目の行がｘiTであり、
【数１００】
2005524124.tif 000088
【０２１６】
これは、以下のように書かれることができる：
【数１０１】
2005524124.tif 000089
【０２１７】
ＥＭアルゴリズムは、以下のステップを含んでいることが好ましい：
（ａ）ｎ＝０，Ｓ0＝｛１，２，・・・，ｐ｝，φ(0) を設定することによりアルゴリズムを初期化し、たとえばε＝１０-5のようなεに対する値を適用し、
ｐ≦Ｎである場合、
【数１０２】
2005524124.tif 000090
【０２１８】
によって初期値β*を計算し、
ｐ＞Ｎである場合、
【数１０３】
2005524124.tif 000091
【０２１９】
によって初期値β*を計算する。ここでリッジパラメータλは、０＜λ≦１を満足させ、ζは小さく、リンク関数ｇがｙ＋ζにおいて明確に定義されるように選択される。
【０２２０】
（ｂ）ｉεＳnのときは、βi(n) ＝βi*，
それ以外のときは、βi(n) ＝０
と定義し、ｐnを、β(n)のゼロでないエレメントγ(n)が
γ(n)＝ＰnTβ(n) ，    β(n)＝Ｐnγ(n)
γ＝ＰnTβ  ，      β＝Ｐnγ
を満足させるような０と１とのマトリックスとする。  
（ｃ）以下の関数を使用して構成要素加重の事後分布の条件付き期待値を計算することにより推定（Ｅ）ステップを行う：
【数１０４】
2005524124.tif 000092
【０２２１】
ここで、ｌはｙの対数尤度関数である。
【０２２２】
β＝Ｐnγおよびβ(n)＝Ｐnγ(n)を使用することにより、この式（１６Ｃ）は以下のように書かれることができる：
【数１０５】
2005524124.tif 000093
【０２２３】
（ｄ）たとえば、Ｎｅｗｔｏｎ  Ｒａｐｈｓｏｎ反復等の、反復的な処理手続きを適用してＱをγの関数として最大化することにより、最大化（Ｍ）ステップを行う。このとき、γ0＝γ(n) であり、ｒ＝０，１，２，・・・に対してγr+1＝γr＋αrδrであり、ここでαrは、
Ｑ（γr+1｜γ(n)，φ(n)）＞Ｑ（γr｜γ(n)，φ(n)）
を確実にするようにラインサーチアルゴリズムによって選択され、また、
ｐ≦Ｎである場合、
【数１０６】
2005524124.tif 000094
【０２２４】
添字ｒはこれらの量がμ＝ｈ（ＸＰnγr）で推定されることを示す。
【０２２５】
ｐ＞Ｎである場合、
【数１０７】
2005524124.tif 000095
【０２２６】
Ｖrおよびｚrは前のように定義されている。
【０２２７】
γ*を、ある収斂基準が満足されたとき、たとえば、以下のようなときのγrの値とする：
｜｜γr−γr+1｜｜＜ε（たとえば、１０-5）
１）
【数１０８】
2005524124.tif 000096
【０２２８】
と定義する。ここで、ε1は、たとえば、ｌｅ−５等の小さい定数である。
【０２２９】
ｎ＝ｎ＋１を設定し、φn+1＝φn＋κn（φ*−φn）を選択し、ここでφ*は、
【数１０９】
2005524124.tif 000097
【０２３０】
を満足させ、κnは０＜κn≦１となるような制動定数である。いくつかの場合において、スケールパラメータは知られており、あるいはこの式は明示的に解かれることができるため、φに対する更新方程式が得られることに注意されたい。
【０２３１】
上記の実施形態は、Wedderburn氏による文献(1974)ならびにMcCullaghおよびNelder氏等による文献(1983)に記載されている擬似尤度方法を含むように拡張されることができる。このようなの実施形態において、詳細に上述された同じ反復的な処理手続きは適切であるが、しかし尤度が上記および、たとえば、McCullaghおよびNelder氏等による文献(1983)中の表８．１に示されている擬似尤度によって置換される。１実施形態において、スケールパラメータφに対する修正された更新方法が存在する。これらのモデルを規定するには、分散関数τ2、リンク関数ｇ、およびリンク関数の導関数：
【数１１０】
2005524124.tif 000098
【０２３２】
の仕様が必要である。
【０２３３】
これらがひとたび規定されると、上記のアルゴリズムは適用されることができる。１実施形態において、擬似尤度モデルに対して、上記のアルゴリズムのステップ５は、スケールパラメータが以下を計算することにより更新されるように修正される：
【数１１１】
2005524124.tif 000099
【０２３４】
ここで、μおよびτはβ*＝Ｐnγ*で推定される。この更新は、モデル中のパラメータの数ｓがＮより小さいときに行われることが好ましい。Ｎ−ｓの除数は、ｓがＮよりはるかに小さいときに使用されることができる。
【０２３５】
別の実施形態において、一般化された線形モデルおよび擬似尤度モデルの両者に対して、行ｘiTを有するコバリエイトマトリックスＸは、ｉｊ番目のエレメントｋijを有し、あるカーネル関数κに対してｋij＝κ（ｘi−ｘj）であるマトリックスＫによって置換されることができる。このマトリックスはまた１のベクトルにより増加されることができる。ある例のカーネルが以下の表３に与えられている［Evgeniou氏他による文献(1999)参照］。
【０２３６】
カーネル関数                          κ（ｘ−ｙ）に対する公式
ガウスのラジアル基礎関数              ｅｘｐ（−｜｜ｘ−ｙ｜｜2／ａ），ａ＞０
逆マルチクアドリック(multiquadric)      （｜｜ｘ−ｙ｜｜2＋ｃ2）-1/2
マルチクアドリック                  （｜｜ｘ−ｙ｜｜2＋ｃ2）1/2
薄板スプライン                    ｜｜ｘ−ｙ｜｜2n+1
｜｜ｘ−ｙ｜｜2nｌｎ（｜｜ｘ−ｙ｜｜）
マルチレーヤ知覚                  適切なθに対して、ｔａｎｈ（ｘ’ｙ−θ）
度数ｄの多項式                    （１＋ｘ’ｙ）d
Ｂスプライン                      Ｂ2n+1（ｘ−ｙ）
三角多項式                 ｓｉｎ（（ｄ＋１／２）（ｘ−ｙ））/ｓｉｎ（（ｘ−
ｙ）／２）
表３：カーネル関数の例
【０２３７】
表３において、最後の２つのカーネルは１次元である、すなわち、Ｘが列を１つしか有しない場合に対して１次元である。多変量バージョンはこれらのカーネル関数の積から得られることができる。Ｂ2n+1の定義はDe Boor氏による文献（1978）に記載されている。一般化された線形モデルまたは擬似尤度モデルのいずれかにおいてカーネル関数を使用した結果、推定された確率はコバリエイトＸの平滑な（線形の変換とは対照的に）関数である平均値が得られる。このようなモデルにより、データに対する実質的に良好な適合が得られる。
【０２３８】
以下、比例ハザードモデルに関する第４の実施形態を説明する。  
Ｄ．比例ハザードモデル
この実施形態の方法は、規定された事象（たとえば、死、回復等）がある時間期間内に発生する確率に影響を与えることのできる構成要素の部分集合を識別するために訓練サンプルを使用することができる。訓練サンプルはシステムと、その訓練サンプルが得られたときからその事象が発生したときまで測定された時間とから得られる。統計的方法を使用して、複数の訓練サンプルから得られたデータにより時間を事象に関連付けることにより、事象までの時間の分布を予想することのできる構成要素の部分集合が識別されることができる。後続的に、その構成要素の部分集合の情報は、たとえば、死までの時間または病気の再発までの時間の統計的特徴を予測するために、たとえば、臨床テスト等のテストに対して使用されることができる。たとえば、システムの構成要素の部分集合からのデータは、ＤＮＡマイクロアレイから得られることができる。このデータは、たとえば、期待されるまたは中間的な患者の生存時間のような臨床的に関連した事象を予測するために、あるいはある症状の始まりまたは病気の再発を予測するために使用されることができる。
【０２３９】
このようにして、本発明は、あるシステムのある事象までの時間の分布を予測するために使用されることのできる最少数の構成要素を識別することが好ましい。この最少数の構成要素が、そのある事象に対する時間の“予測的な”ものである。本質的に、本発明の方法は、システムから発生された全てのデータから、ある事象までの時間を予測するために使用されることのできる最少数の構成要素を識別することができる。この方法によってこれらの構成要素が識別されると、その構成要素は、将来、あるシステムのある事象までの時間の統計的特徴を新しいサンプルから予測するために使用されることができる。本発明の方法は統計的方法を使用して、あるシステムのある事象までの時間を正しく予測するのに必要でない構成要素を除去することが好ましい。
【０２４０】
この明細書において使用されているように、“ある事象までの時間”とは、本発明の方法が適用されるサンプルを獲得したときからある事象の時間までの時間の尺度のことである。ある事象はどのような観察可能な事象であってもよい。システムが生物学的システムである場合、事象は、たとえば、あるシステムの故障までの時間、死までの時間、特定の症状の始まり、病気の始まりまたは再発、表現型または遺伝子型の変化、生化学的特徴の変化、組織体または組織の形態学的特徴の変化、行動の変化であることができる。
【０２４１】
ある事象までの以前の時間からのサンプルがある事象までの特定の時間と関連付けられる。ある事象までの時間は、たとえば、サンプリングしたときから死までの時間が知られている患者から得られたデータ、換言すると“本当の”生存時間と、およびサンプルが最後に得られたときに生きていたことが唯一の情報である患者から得られたデータ、換言すると、特定の患者が少なくとも所定の日数生存していることを示す“検閲削除された”生存時間とから決定された時間であってもよい。
【０２４２】
１実施形態において、入力データはＮ×ｐマトリックスＸ＝（ｘij）に編成され、ここで、Ｎはテスト訓練サンプルであり、ｐは構成要素である。典型的に、ｐはＮよりはるかに大きい。
【０２４３】
たとえば、マイクロアレイ実験からのＮ×ｐデータマトリックスＸ＝（ｘij）を考慮する。ここで、Ｎは個人（またはサンプル）であり、各個人に対して同じｐの遺伝子が存在する。たとえば、生存時間等のある事象までの時間を示す変数ｙi（ｙi≧０）が各個人ｉ（ｉ＝１，２，・・・，Ｎ）と関連付けられることが好ましい。各個人に対して、その個人の生存時間が本当の生存時間であるか、あるいは検閲削除された生存時間であるかを示す変数が規定されることができることも好ましい。ｃiのような検閲削除インジケータは、
ｙiが検閲削除されていない場合、ｃi＝１、
ｙiが検閲削除された場合、ｃi＝０
を表す。
【０２４４】
生存時間ｙiを有するＮ×１ベクトルは、
【数１１２】
2005524124.tif 000100
【０２４５】
と書くことができ、検閲削除インジケータｃiを有するＮ×１ベクトルは、
【数１１３】
2005524124.tif 000101
【０２４６】
と書くことができる。
【０２４７】
典型的に、上述したように、構成要素加重は、その構成要素加重の大部分がゼロであるというアプリオリ仮定が何を考慮したかで推定される。
【０２４８】
構成要素加重に対して特定されたプライアは以下の形態のものである：
【数１１４】
2005524124.tif 000102
【０２４９】
ここで、β1，β2，・・・，βnは構成要素加重であり、Ｐ（βi｜τi）はＮ（０，τi2）であり、Ｐ（τi）α１／τi2はジェフェリーのプライアである（KotzおよびJohnson氏等の文献、1983）。
【０２５０】
尤度関数は、データの分布に基づいてデータに適合したモデルを規定する。尤度関数は以下の形態のものであることが好ましい：
【数１１５】
2005524124.tif 000103
【０２５１】
尤度関数によって規定されたモデルは、あるシステムのある事象までの時間を予測する任意のモデルであってもよい。  
１実施形態において、尤度関数によって規定されたモデルは、コックスの比例ハザードモデルである。コックスの比例ハザードモデルはコックス氏によって導入され(1972年)、好ましいことに、生存データに対する回帰モデルとして使用されることができる。コックスの比例ハザードモデル：
【数１１６】
2005524124.tif 000104
【０２５２】
は、構成要素と関連付けられた（説明的な）パラメータのベクトルである。本発明の方法は、データ：
【数１１７】
2005524124.tif 000105
【０２５３】
を与えられたコックスの比例ハザードモデルに対して、パラメータ：
【数１１８】
2005524124.tif 000106
【０２５４】
からの非常に少数のものの選択（および除去）を行うことが好ましい。
【０２５５】
コックスの比例ハザードモデルの適用は、同じ生存時間に対して種々のデータがシステムから得られる状況において、換言すると、同点の生存時間が生じた場合に問題となる可能性が高い。同点の生存時間は予備処理ステップを実施され、それによって特有の生存時間にされることができる。提案された予備処理により、コックスの比例ハザードモデルが後続的に適用されたときに同点の生存時間に関する問題が回避されるため、必然的に続いて発生するアルゴリズムは簡単なものになる。
【０２５６】
生存時間の予備処理は、非常に小量のランダム雑音を追加することにより適用される。この処理手続きは、同点の時間のセットを選び、ゼロ平均と、分類された生存時間の間の最小の非ゼロ距離に比例した分散とを有する正規分布から引き出されたランダムな量を同点の時間のセット内の各同点時間に付加することであることが好ましい。このような予備処理により、生存時間の厳しい摂動を生ぜずに同点時間の除去が行われる。
【０２５７】
予備処理により、異なった生存時間が発生する。好ましいことに、これらの時間は、以下に表されるように大きさが増加していく順に順序付けられることができる：
【数１１９】
2005524124.tif 000107
【０２５８】
Ｘの行の再配列であるＮ×ｐマトリックスをＺで表し、ここで、Ｚの行の順序付けは、
【数１２０】
2005524124.tif 000108
【０２５９】
の順序付けにより誘発された順序付けに対応し、また、マトリックスＺのｊ番目の行はＺjで表される。ｄを、
【数１２１】
2005524124.tif 000109
【０２６０】
を順序付けるために必要とされる同じ置換によりｃを順序付けた結果とする。
【０２６１】
同点の生存時間に対する予備処理を考慮し、生存データ解析に関する標準テキスト（たとえば、CoxおよびOakes氏等の文献，1984）への参照を行った後、確率ハザードモデルに対する尤度関数を以下のように書くことができることが好ましい：
【数１２２】
2005524124.tif 000110
【０２６２】
尤度の対数（すなわち、ｌ＝ｌｏｇ（Ｌ））は以下のように書くことができることが好ましい：
【数１２３】
2005524124.tif 000111
【０２６３】
ここで、
ｊ＜ｉである場合、ζi,j＝０
ｊ≧ｉである場合、ζi,j＝１
このモデルは、生存分布のパラメトリック形態が特定されておらず、生存時間の通常の特性だけが使用されることが好ましい（リスクセットの決定において）ため、非パラメトリックであることに注意されたい。これは非パラメトリックなケースであるため、
【数１２４】
2005524124.tif 000112
【０２６４】
は必要とされない（すなわち、ｑ＝０）。
【０２６５】
本発明の方法の別の実施形態において、尤度関数により規定されるモデルはパラメトリック生存モデルである。パラメトリック生存モデルでは、
【数１２５】
2005524124.tif 000113
【０２６６】
は構成要素に関連付けられた（説明的な）パラメータのベクトルであり、
【数１２６】
2005524124.tif 000114
【０２６７】
は生存密度関数の関数形態に関連付けられたパラメータのベクトルであることが好ましい。好ましくは、本発明の方法は、データ：
【数１２７】
2005524124.tif 000115
【０２６８】
を与えられたパラメトリック生存モデルに対して、パラメータ：
【数１２８】
2005524124.tif 000116
【０２６９】
および
【数１２９】
2005524124.tif 000117
【０２７０】
の推定値からの非常に少数のものの選択（および除去）を行うことが好ましい。
【０２７１】
パラメトリック生存モデルを適用するとき、生存時間は予備処理を必要とせず、
【数１３０】
2005524124.tif 000118
【０２７２】
として表される。パラメトリック生存モデルは次のように適用される：
生存時間のパラメトリック密度関数を、
【数１３１】
2005524124.tif 000119
【０２７３】
で表し、その生存関数を、
【数１３２】
2005524124.tif 000120
【０２７４】
で表し、ここで、
【数１３３】
2005524124.tif 000121
【０２７５】
は密度関数のパラメトリック形態に関連したパラメータであり、
【数１３４】
2005524124.tif 000122
【０２７６】
は上記に規定されたものである。ハザード関数は以下のように規定される：
【数１３５】
2005524124.tif 000123
【０２７７】
対数尤度関数の一般的な公式は、検閲削除されたデータを考慮して、以下のように表されることが好ましい：
【数１３６】
2005524124.tif 000124
【０２７８】
パラメトリック回帰生存モデルを介した生存時間データの解析に関する標準的なテキストへの参照により、使用されることのできる生存時間分布の集りが示される。使用されることのできる生存分布には、たとえば、ワイブル、指数または極値分布が含まれる。
【０２７９】
ハザード関数が
【数１３７】
2005524124.tif 000125
【０２８０】
として書かれることができる場合、
【数１３８】
2005524124.tif 000126
【０２８１】
ワイブル、指数または極値分布は、このようにして表された形態で書かれることのできる密度およびハザード関数を有している。
【０２８２】
詳細な適用は部分的に、AitkenおよびClayton氏等の文献(1980)に記載されているアルゴリズムに依存するが、しかしながら、それによってユーザは任意のパラメトリック基礎ハザード関数を指定することが可能になる。
【０２８３】
AitkinおよびClayton氏等の文献(1980)によると、パラメトリック生存モデルをモデル化する好ましい尤度関数は以下のように表される：
【数１３９】
2005524124.tif 000127
【０２８４】
AitkinおよびClayton氏等の文献(1980)には、式（５Ｄ）の結果として、ｃiは平均μiを有するポアソン変量として処理されてもよく、また、式（１１Ｄ）中の最後の項は、
【数１４０】
2005524124.tif 000128
【０２８５】
に依存しない（もっとも、それは
【数１４１】
2005524124.tif 000129
【０２８６】
に依存する）ことが示されている。
【数１４２】
2005524124.tif 000130
【０２８７】
の事後分布は以下のように表されることが好ましい：
【数１４３】
2005524124.tif 000131
【０２８８】
１実施形態において、
【数１４４】
2005524124.tif 000132
【０２８９】
は失われたデータのベクトルとして処理され、式（６Ｄ）を最大化して、
【数１４５】
2005524124.tif 000133
【０２９０】
の事後推定値を生成するために反復的な処理手続きが使用されることができる。式（１Ｄ）のプライアは、極大である事後推定値が希薄になる、すなわち、多数のパラメータが冗長である場合に、
【数１４６】
2005524124.tif 000134
【０２９１】
の多くの構成要素がゼロになるようなものである。
【数１４７】
2005524124.tif 000135
【０２９２】
の多くの構成要素がゼロであるプライア期待値が存在するため、推定は、推定されたβiの大部分がゼロであり、残りのゼロでない推定値が生存時間の適切な説明を提供するように行われる。
【０２９３】
マイクロアレイデータの状況においては、この動作の結果、その事象時間に対する適切な説明を提供する非常に少数の遺伝子のセットが識別されることとなる。
【０２９４】
上述したように、事後分布を最大化する構成要素加重は、反復的な処理手続きを使用して決定されることができる。構成要素および構成要素加重の事後分布を最大化するための反復的な処理手続きは、たとえば、Dempster氏他による文献(1977)に記載されているようなＥＭアルゴリズムであることが好ましい。
【０２９５】
１実施形態において、ＥＭアルゴリズムは以下のステップを含んでいる：
１．ｎ＝０，Ｓ0＝｛１，２，・・・，ｐ｝を設定することによりアルゴリズムを初期化し、
【数１４８】
2005524124.tif 000136
【０２９６】
を初期化する
２．    ｉεＳnのときは、βi(n) ＝βi*，
それ以外のときは、βi(n) ＝０
と定義し、Ｐnを、
【数１４９】
2005524124.tif 000137
【０２９７】
ような０と１とのマトリックスとする。  
３．構成要素加重の事後分布の期待値を計算することにより推定ステップを行う。これは以下の関数を使用して行われる：
【数１５０】
2005524124.tif 000138
【０２９８】
４．最大化ステップを行う。これは、以下のようにＮｅｗｔｏｎ  Ｒａｐｈｓｏｎ反復を使用して行われることができる：
【数１５１】
2005524124.tif 000139
【０２９９】
を設定し、ｒ＝０，１，２，・・・に対して
【数１５２】
2005524124.tif 000140
【０３００】
であり、ここでαrは、以下を保証するようにラインサーチアルゴリズムによって選択される：
【数１５３】
2005524124.tif 000141
【０３０１】
５．
【数１５４】
2005524124.tif 000142
【０３０２】
と定義する。ここでε1は、たとえば、１０-5等の小さい定数である。ｎ＝ｎ＋１を設定し、
【数１５５】
2005524124.tif 000143
【０３０３】
を選択する。ここで、
【数１５６】
2005524124.tif 000144
【０３０４】
また、κnは、０＜κn＜１であるような減衰定数である。
【０３０５】
６．収斂をチェックする。
【数１５７】
2005524124.tif 000145
【０３０６】
である場合、中止し、ここでε2は適切に小さく、その他の場合、上記のステップ２に進む。
【０３０７】
別の実施形態において、最大化ステップの中のステップ（４）は、
【数１５８】
2005524124.tif 000146
【０３０８】
をその期待値：
【数１５９】
2005524124.tif 000147
【０３０９】
で置換することにより推定されることができる。
【０３１０】
１実施形態において、ＥＭアルゴリズムは、モデルがコックスの比例ハザードモデルであるときの事後分布を最大化するために適用される。
【０３１１】
モデルがコックスの比例ハザードモデルであるときのＥＭアルゴリズムの適用の説明を助けるために、“動的加重”およびこれらの加重に基づいたマトリックスを規定することが好ましい。加重は以下のように表される：
【数１６０】
2005524124.tif 000148
【０３１２】
これらの加重に基づいたマトリックスは以下のように表される：
【数１６１】
2005524124.tif 000149
【０３１３】
加重のマトリックスに関して、ｌの第１および第２の微分を以下のように書いてもよい：
【数１６２】
2005524124.tif 000150
【０３１４】
ここで、Ｋ＝Ｗ**−Δ（Ｗ*）である。したがって、ＥＭアルゴリズムのステップ（２）の一部（式７Ｄ）（式１０Ｄも参照）として示された変換マトリックスＰnから、
【数１６３】
2005524124.tif 000151
【０３１５】
となることに注意しなければならない。
【０３１６】
モデルがコックスの比例ハザードモデルであるとき、ＥＭアルゴリズムのＥステップおよびＭステップは、以下のように行われることが好ましい：
１．ｎ＝０，Ｓ0＝｛１，２，・・・，ｐ｝を設定する。ｖを、構成要素Ｖiを有するベクトルとする。ここで、たとえば０．００１のような、ある小さいεに対して、
ｃi＝１である場合、Ｖi＝１−ε，
ｃi＝０である場合、Ｖi＝ε
である。ｆはｌｏｇ（ｖ／ｔ）であると定義する。
【０３１７】
ｐ≦Ｎである場合、
【数１６４】
2005524124.tif 000152
【０３１８】
によって初期値：
【数１６５】
2005524124.tif 000153
【０３１９】
を計算し、
ｐ＞Ｎである場合、
【数１６６】
2005524124.tif 000154
【０３２０】
によって初期値：
【数１６７】
2005524124.tif 000155
【０３２１】
を計算する。ここで、リッジパラメータλは、０＜λ≦１を満足させる。
【０３２２】
２．    ｉ∈Ｓnの場合、βi(n) ＝βi*，
それ以外の場合、βi(n) ＝０
と定義する。  
Ｐnを、
【数１６８】
2005524124.tif 000156
【０３２３】
ような０と１のマトリックスとする。
【０３２４】
３．以下を計算することによりＥステップを行う： 
【数１６９】
2005524124.tif 000157
【０３２５】
β＝Ｐnγおよびβ(n)＝Ｐnγ(n)を使用することにより、以下の式が得られる：
【数１７０】
2005524124.tif 000158
【０３２６】
４．Ｍステップを行う。これは、以下のようにＮｅｗｔｏｎ  Ｒａｐｈｓｏｎ反復により行われることができる。
【数１７１】
2005524124.tif 000159
【０３２７】
を設定し、ｒ＝０，１，２，・・・に対して、
【数１７２】
2005524124.tif 000160
【０３２８】
ここで、αrは、
【数１７３】
2005524124.tif 000161
【０３２９】
を保証するようにラインサーチアルゴリズムによって選択される。  
ｐ≦Ｎに対しては、以下の式を使用する：
【数１７４】
2005524124.tif 000162
【０３３０】
ｐ＞Ｎに対しては、以下の式を使用する：
【数１７５】
2005524124.tif 000163
【０３３１】
γ*を、ある収斂基準が満足されたとき、たとえば、
｜｜γr−γr+1｜｜＜ε（たとえば、１０-5）
であるときのγrの値とする。
【０３３２】
５．
【数１７６】
2005524124.tif 000164
【０３３３】
と定義する。ここで、ε1は、たとえば、１０-5等の小さい定数である。このステップにより、非常に小さい係数を有する変数が除去される。
【０３３４】
６．収斂をチェックする。
【数１７７】
2005524124.tif 000165
【０３３５】
である場合、中止し、ここでε2は適切に小さく、その他の場合には、ｎ＝ｎ＋１を設定し、上記のステップ２に進んで、収斂が発生するまで処理手続きを繰返す。
【０３３６】
別の実施形態において、ＥＭアルゴリズムは、モデルがパラメトリック生存モデルであるときの事後分布を最大化するために適用される。
【０３３７】
パラメトリック生存モデルにＥＭアルゴリズムを適用したとき、式（５Ｄ）の結果として、ｃiは平均μiを有するポアソン変量として処理されてもよく、また、式（５Ｄ）中の最後の項はβに依存しない（もっとも、それはφに依存する）。
【数１７８】
2005524124.tif 000166
【０３３８】
であり、したがって、ポアソンのような平均に対する対数・線形モデルで問題を表すことができることに注意されたい。対数尤度関数の反復的な最大化が行われ、
【数１７９】
2005524124.tif 000167
【０３３９】
の初期推定値が与えられた場合に
【数１８０】
2005524124.tif 000168
【０３４０】
の推定値が得られることが好ましい。それ故、
【数１８１】
2005524124.tif 000169
【０３４１】
のこれらの推定値が与えられた場合には、
【数１８２】
2005524124.tif 000170
【０３４２】
の更新された推定値が得られる。この処理手続きは、収斂が発生するまで続けられる。
【０３４３】
上述した事後分を適用することにより、
【数１８３】
2005524124.tif 000171
【０３４４】
であることが認識される。  
結果的に、式（１１Ｄ）および（１２Ｄ）から、
【数１８４】
2005524124.tif 000172
【０３４５】
パラメトリック生存モデルに関連した式（１２Ｄ）の別の形は以下のとおりである：
【数１８５】
2005524124.tif 000173
【０３４６】
ＥＭアルゴリズムの各Ｍステップの後に
【数１８６】
2005524124.tif 000174
【０３４７】
について解く（以下の、ステップ５を参照）ために、
【数１８７】
2005524124.tif 000175
【０３４８】
とし、ここで、０＜κn≦１に対して
【数１８８】
2005524124.tif 000176
【０３４９】
βは前のＭステップから得られた値に固定されていることが好ましい。
【０３５０】
パラメトリック生存モデルおよびマイクロアレイデータの状況においてパラメータ選択のためのＥＭアルゴリズムを提供することができる。このＥＭアルゴリズムは以下のようなものであることが好ましい：
１．
【数１８９】
2005524124.tif 000177
【０３５１】
を設定する。ｖを、構成要素Ｖiを有するベクトルとする。ここで、たとえば０．００１のような、ある小さいεに対して、
ｃi＝１である場合、Ｖi＝１−ε，
ｃi＝０である場合、Ｖi＝ε
である。ｆはｌｏｇ（ｖ／Λ（ｙ，φ））であると定義する。  
ｐ≦Ｎである場合、
【数１９０】
2005524124.tif 000178
【０３５２】
によって初期値：
【数１９１】
2005524124.tif 000179
【０３５３】
を計算し、
ｐ＞Ｎである場合、
【数１９２】
2005524124.tif 000180
【０３５４】
によって初期値：
【数１９３】
2005524124.tif 000181
【０３５５】
を計算する。ここで、リッジパラメータλは、０＜λ≦１を満足させる。
【０３５６】
２．    ｉ∈Ｓnの場合、βi(n) ＝βi*，
それ以外の場合、βi(n) ＝０
と定義する。  
Ｐnを、
【数１９４】
2005524124.tif 000182
【０３５７】
ような０と１のマトリックスとする。
【０３５８】
３．以下を計算することによりＥステップを行う： 
【数１９５】
2005524124.tif 000183
【０３５９】
β＝Ｐnγおよびβ(n)＝Ｐnγ(n)を使用することにより、以下の式が得られる：
【数１９６】
2005524124.tif 000184
【０３６０】
４．Ｍステップを行う。これは、以下のようにＮｅｗｔｏｎ  Ｒａｐｈｓｏｎ反復により行われることができる。
【数１９７】
2005524124.tif 000185
【０３６１】
を設定し、ｒ＝０，１，２，・・・に対して、
【数１９８】
2005524124.tif 000186
【０３６２】
ここで、αrは、
【数１９９】
2005524124.tif 000187
【０３６３】
を保証するようにラインサーチアルゴリズムによって選択される。
【０３６４】
ｐ≦Ｎに対しては、以下の式を使用する：
【数２００】
2005524124.tif 000188
【０３６５】
ｐ＞Ｎに対しては、以下の式を使用する：
【数２０１】
2005524124.tif 000189
【０３６６】
γ*を、ある収斂基準が満足されたとき、たとえば、
｜｜γr−γr+1｜｜＜ε（たとえば、１０-5）
であるときのγrの値とする。
【０３６７】
５．
【数２０２】
2005524124.tif 000190
【０３６８】
と定義する。ここで、ε1は、たとえば、１０-5等の小さい定数である。ｎ＝ｎ＋１を設定し、
【数２０３】
2005524124.tif 000191
【０３６９】
を選択し、ここで、
【数２０４】
2005524124.tif 000192
【０３７０】
κnは０＜κn＜１となるような減衰定数である。
【０３７１】
６．収斂をチェックする。
【数２０５】
2005524124.tif 000193
【０３７２】
である場合、中止し、ここでε2は適切に小さく、その他の場合には、ステップ２に進む。
【０３７３】
別の実施形態において、生存時間はワイブル生存密度関数によって表される。ワイブルの場合、
【数２０６】
2005524124.tif 000194
【０３７４】
は１次元であり、
【数２０７】
2005524124.tif 000195
【０３７５】
であることが好ましい。
【０３７６】
αの更新された値を提供するために、各Ｍステップの後に
【数２０８】
2005524124.tif 000196
【０３７７】
が解かれることが好ましい。コックスの比例ハザードモデルに対して適用されたステップに後続して、αが推定されることができ、また、生存時間がワイブル分布に従ったものである場合にその生存時間を適切に説明することのできる
【数２０９】
2005524124.tif 000197
【０３７８】
からのパラメータの非常に少数の部分集合が選択されることができる。
【０３７９】
本発明の特徴および利点は、以下の例の説明から明らかになるであろう。
【０３８０】
例
例１：対数回帰モデルを使用する前立腺癌に対する２グループ分類
組織を前立腺および非前立腺グループに分類することのできる遺伝子の部分集合を識別するために、Luo氏他の文献(2001)で報告され、解析されたマイクロアレイデータセットは、二項ロジスティック回帰がモデルとして使用される本発明の方法を使用して解析された。このデータセットは６５００の人の遺伝子に関するマイクロアレイを含んでいる。研究対象には、前立腺癌を有することが知られている１６人の被験者と、良性の前立腺細胞過形成である９人の被験者が含まれている。しかしながら、表現を単に簡潔にするために、５０の遺伝子が解析のために選択された。この５０の全ての遺伝子の遺伝子発現率（行）および２５人の患者（列）は表４に示されている。
【０３８１】
以下、この方法を適用した結果を示す。モデルはＧ＝２クラスを有しており、５０の遺伝子の全てをそのモデル中の前立腺変数（構成要素または基礎関数）として開始された。２１回の反復後（以下を参照）、アルゴリズムは、完全な分類を与える２つの遺伝子（表５の番号３６および４７）を発見した。この結果がデータセット中で利用可能な多数の遺伝子（変数）による人為結果だったのか否かを決定するために、クラスラベルがランダムに置換され、後続的にアルゴリズムが適用される置換テストがランされた。これは２００回繰返された。図１は、正しく分類された場合の数のヒストグラムを示している。実際のデータセットに対する１００％の正確さはこの置換分布の末尾にあり、０．０１５のｐ値を有する。これは、この結果が偶然によるものではないことを示唆している。
【０３８２】
以下に交換されていないデータに対する反復の詳細を示す：
【表２】
2005524124.tif 000198
【表３】
2005524124.tif 000199
【表４】
2005524124.tif 000200
【表５】
2005524124.tif 000201
【表６】
2005524124.tif 000202
【表７】
2005524124.tif 000203
【表８】
2005524124.tif 000204
【表９】
2005524124.tif 000205
【表１０】
2005524124.tif 000206
【表１１】
2005524124.tif 000207
【表１２】
2005524124.tif 000208
【表１３】
2005524124.tif 000209
【表１４】
2005524124.tif 000210
【表１５】
2005524124.tif 000211
【表１６】
2005524124.tif 000212
【０３８３】
例２：大きいデータセットおよび二項ロジスティック回帰モデルを使用する２グループ分類
組織を異なった臨床タイプのリンパ腫に分類することのできる遺伝子の部分集合を識別するために、Alizadeh,A.A.氏他の文献［Distinct types of diffuse large B-cell lymphoma identified by gene expression profiling,Nature 403:503-511(2000)］で報告され、解析されたデータセットは、二項ロジスティック回帰がモデルとして使用される本発明の方法を使用して解析された。
【０３８４】
データセットには、ｎ＝４０２６の遺伝子およびｎ＝４２のサンプルが存在している。以下において、ＤＬＢＣＬとは“びまん性大Ｂ細胞リンパ腫”のことである。サンプルは、ＧＣ  Ｂ様ＤＬＢＣＬ（２１個のサンプル）および活性化されたＢ様ＤＬＢＣＬ（２１個のサンプル）の２つの疾病タイプに分類されている。このセットを使用して、異なった疾病タイプの診断に役立つ遺伝子を迅速に発見するための上記の方法の使用法を説明する。
【０３８５】
以下において、この方法を適用した結果が与えられている。モデルはＧ＝２クラスを有しており、全ての遺伝子をそのモデル中の潜在的変数（基礎関数）として開始された。２０回の反復後、アルゴリズムは、以下の誤分類（表５）を与える番号１２８１および１３１２（GENE3332XおよびGENE3258X）の２つの遺伝子を発見し、全体的な分類成功率は９８％であった。この例は、ラップトップマシンで約２０秒で稼動した。
【０３８６】
表  ５
予測されたクラス１            予測されたクラス２
真のクラス１              ２０                                １
真のクラス２                ０                              ２１
【０３８７】
結果がデータセット中で利用可能な多数の遺伝子（変数）による人為結果だったのか否かを決定するために、クラスラベルがランダムに置換され、後続的にアルゴリズムが適用される置換テストがランされた。これは１０００回繰返された。図２は、正しく分類された場合の割合（ラムダ）のヒストグラムを示している。実際のデータセットに対する９７．６％の正確さはこの置換分布の末尾にあり、０．０１３のｐ値を有する。これらの観察は、この結果が偶然によるものではないことを示唆している。
【０３８８】
例３：多グループ分類
サンプルを多数のクラスの１つに分類することのできる遺伝子を識別するために、Yeoh氏他の文献［Cancer Cell v1 : 133-143(2002)］で報告され、解析されたデータセットは、尤度が多項ロジスティック回帰に基づいて使用される本発明の方法を使用して解析された。Yeoh氏他の文献に記載されているものと同じ予備処理が適用されている。これは、以下から構成されている：
・８つのアレイ：ＢＣＲ．ＡＢＬ．Ｒ４，ＭＬＬ．Ｒ５，Ｎｏｒｍａｌ．Ｒ４，Ｔ．ＡＬＬ．Ｒ７，Ｔ．ＡＬＬ．Ｒ８，Ｈｙｐｅｒｄｉｐ．５０．２Ｍ．３，Ｈｙｐｏｄｉｐ．２Ｍ．３，およびＨｙｐｏｄｉｐ．２Ｍ．２をドロップする
・各アレイの平均反応値を２５００に設定する
・しきい値化する；１００より小さい４５０００個の値に設定された４５０００個を越える値が１に設定される
・存在が０．０１未満の遺伝子が除去される；これは合計１６０７個の遺伝子である
・最大値と最小値との間の差が１００未満であった遺伝子（１６０４個）が除去される
【０３８９】
予備処理後、ｎ＝１１００５個の遺伝子およびｎ＝２４８個のサンプルが存在している。サンプルは６つの疾病タイプに分類されている：
１．ＢＣＲ−ＡＢＬ；
２．Ｅ２Ａ−ＰＢＸ１；
３．Ｈｙｐｅｒｄｉｐ＞５０；
４．ＭＬＬ；
５．Ｔ−ＡＬＬおよび
６．ＴＥＬ−ＡＭＬ１
【０３９０】
このセットは、異なった疾病タイプの診断に役立つ遺伝子を迅速に発見するための方法の使用法を示すために使用された。この方法を適用した結果を以下に示す。モデルはＧ＝６のクラスを有しており、全ての遺伝子をそのモデル中の潜在的変数（基礎関数）として開始された。２０回の反復の後、アルゴリズムは、以下の１０個の遺伝子がクラスを分けることを発見した：
【表１７】
2005524124.tif 000213
【０３９１】
１５重の相互確認（ｃｒｏｓｓ  ｖａｌｉｄａｔｉｏｎ）から、分類成功率が９４％である以下の誤分類表（表６）が得られた：
【表１８】
2005524124.tif 000214
【０３９２】
交換テスト（クラスラベルを交換する）は、０．９４％の相互確認されたエラーレートが非常に重要であること（ｐ＝０．００）を示した。
【０３９３】
例４：一般化された線形モデルを使用する標準回帰
この例は、この方法が一般化された線形モデルフレームワークでどのようにして実施されることができるかを示す。この例は、２００の観察および４１の変数（基礎関数）に関する標準回帰問題である。真の曲線はエラー（または雑音）を伴なって観察され、変数の一部だけに依存することが知られている。反応は連続し、正規分布される。これらのデータは、一般化された線形モデル変数選択のために本発明のアルゴリズムを使用して解析される。
【０３９４】
これは、以下を有する一般化された線形モデルである：
リンク関数：ｇ（μ）＝μ
リンク関数の偏微分：
【数２１０】
2005524124.tif 000215
【０３９５】
分散関数：τ2＝１
スケールパラメータφ＝σ2
逸脱（尤度関数）：
【数２１１】
2005524124.tif 000216
【０３９６】
σ2に対する更新公式は
【数２１２】
2005524124.tif 000217
【０３９７】
であり、ここでμi*はアルゴリズムのステップ５においてβ*で推定された平均である。
【０３９８】
以下にアルゴリズムの出力を示す。  
ＥＭ反復：１    期待事後分布：-55.45434    基礎関数  ４１    シグマの２乗  0.5607509
ＥＭ反復：２    期待事後分布：-43.96193    基礎関数  ４１    シグマの２乗  0.5773566
ＥＭ反復：３    期待事後分布：-48.87198    基礎関数  ３９    シグマの２乗  0.5943395
ＥＭ反復：４    期待事後分布：-52.79632    基礎関数  ３１    シグマの２乗  0.6072137
ＥＭ反復：５    期待事後分布：-55.18578    基礎関数  ２８    シグマの２乗  0.6161707
ＥＭ反復：６    期待事後分布：-56.5303     基礎関数  ２３    シグマの２乗  0.6224545
ＥＭ反復：７    期待事後分布：-57.47589    基礎関数  １７    シグマの２乗  0.626674
ＥＭ反復：８    期待事後分布：-58.0566     基礎関数  １５    シグマの２乗  0.6293923
ＥＭ反復：９    期待事後分布：-58.41912    基礎関数  １３    シグマの２乗  0.6315789
ＥＭ反復：１０  期待事後分布：-58.6923     基礎関数  １１    シグマの２乗  0.633089
ＥＭ反復：１１  期待事後分布：-58.88766    基礎関数  １０    シグマの２乗  0.6343793
ＥＭ反復：１２  期待事後分布：-59.05261    基礎関数  １０    シグマの２乗  0.635997
ＥＭ反復：１３  期待事後分布：-59.24126    基礎関数  ９    シグマの２乗  0.6381456
ＥＭ反復：１４  期待事後分布：-59.47668    基礎関数  ９    シグマの２乗  0.640962
ＥＭ反復：１５  期待事後分布：-59.7677     基礎関数  ９    シグマの２乗  0.6443392
ＥＭ反復：１６  期待事後分布：-60.10277    基礎関数  ９    シグマの２乗  0.6477088
ＥＭ反復：１７  期待事後分布：-60.44193    基礎関数  ９    シグマの２乗  0.6508144
ＥＭ反復：１８  期待事後分布：-60.7684     基礎関数  ９    シグマの２乗  0.6539145
ＥＭ反復：１９  期待事後分布：-61.09251    基礎関数  ９    シグマの２乗  0.6565873
ＥＭ反復：２０  期待事後分布：-61.38427    基礎関数  ８    シグマの２乗  0.6589498
ＥＭ反復：２１  期待事後分布：-61.65061    基礎関数  ８    シグマの２乗  0.6615976
ＥＭ反復：２２  期待事後分布：-61.92217    基礎関数  ８    シグマの２乗  0.664281
ＥＭ反復：２３  期待事後分布：-62.17683    基礎関数  ７    シグマの２乗  0.6663748
ＥＭ反復：２４  期待事後分布：-62.37402    基礎関数  ７    シグマの２乗  0.6679655
ＥＭ反復：２５  期待事後分布：-62.51645    基礎関数  ７    シグマの２乗  0.6689011
ＥＭ反復：２６  期待事後分布：-62.59567    基礎関数  ６    シグマの２乗  0.6689011
ＥＭ反復：２７  期待事後分布：-62.6151     基礎関数  ６    シグマの２乗  0.6690962
ＥＭ反復：２８  期待事後分布：-62.61717    基礎関数  ６    シグマの２乗  0.6691031
ＥＭ反復：２９  期待事後分布：-62.61739    基礎関数  ５    シグマの２乗  0.6691035
【０３９９】
アルゴリズムは、４１の基礎ベクトル（変数）の中の５つを含むモデルにより収斂する。図３には、アルゴリズムにより選択された５つの変数（基礎関数）を有するモデルに対するあてはめられた曲線（実線）、真の曲線（破線）、および雑音の多いデータのプロットが示されており、ここでｙ変数はｎｆで示されている。
【０４００】
例５：一般化された線形モデルを使用する小さい線形回帰の例
この例は例４に類似しているが、しかし簡潔にするために、少数の変数（１０）が使用される。これによって、全データセットが一覧表にされることが可能になる（表７参照）。従属変数は最初の４つの変数だけの関数であり、残りの変数は雑音である。
【０４０１】
データは、アイデンティティリンク、一定の分散および正常な反応を有する一般化された線形モデルとして解析された。１２回の反復後、アルゴリズムは、予測情報を有することが知られている４つの変数を含み、６つの雑音変数を全て廃棄した解に収斂した。
【表１９】
2005524124.tif 000218
【表２０】
2005524124.tif 000219
【表２１】
2005524124.tif 000220
【表２２】
2005524124.tif 000221
【表２３】
2005524124.tif 000222
【表２４】
2005524124.tif 000223
【表２５】
2005524124.tif 000224
【表２６】
2005524124.tif 000225
【表２７】
2005524124.tif 000226
【表２８】
2005524124.tif 000227
【表２９】
2005524124.tif 000228
【表３０】
2005524124.tif 000229
【０４０２】
例６：リンパ種生存解析
この例は、http://llmpp.nih.gov/lymphoma/data.shtml からの実際の生存データを使用する。
【０４０３】
参考文献は、Alizadeh AA氏他の文献(2000)［Distinct types of diffuse large B-cell lymphoma identified by gene expression profiling. Nature 403(6769):503-11］である。
【０４０４】
データは４０２６個の遺伝子および４０個のサンプル（個人）に対するデータからなるマイクロアレイデータであり、生存時間および検閲削除インジケータが各サンプルに対して利用可能である。結果は、アルゴリズムを使用し、コックスの比例ハザードモデルを実施して解析された。
【０４０５】
アルゴリズムは３つの遺伝子（遺伝子：３７９７Ｘ，３３０２Ｘ，３５６Ｘ）を生存時間と関連付けられたものとして選択したことに注意されたい。
【０４０６】
例７：減少されたリンパ腫生存解析
文書による証明を完全なものにするために、ランダムに選択された４７個の遺伝子と全データセットの解析で意味のあるものとして識別された３個の遺伝子とを含む５０個の遺伝子が、Alizadeh氏他の文献に記載されている遺伝子の部分集合に基づいて選択された例もまた示されている。データは以下の表９に示されている。
【０４０７】
データは、コックスの比例ハザード生存モデルを含むアルゴリズムの別形式を使用して解析された。２２回の反復後、全セットに対する解からの２つの遺伝子を含む５つの遺伝子が選択された。全ての結果（反復履歴を含む）を以下に示す：
【表３１】
2005524124.tif 000230
【表３２】
2005524124.tif 000231
【表３３】
2005524124.tif 000232
【表３４】
2005524124.tif 000233
【表３５】
2005524124.tif 000234
【表３６】
2005524124.tif 000235
【表３７】
2005524124.tif 000236
【表３８】
2005524124.tif 000237
【０４０８】
例８：パラメトリックハザードによる生存解析
データは、
【数２１３】
2005524124.tif 000238
【０４０９】
からの1694w.datである。これは、黒色腫の生存に関するデータである。ｎ＝２５５の個人が存在し、その中の１００人は生存時間を検閲削除されている。各個人は４つのコバリエイト、すなわち、治療、シックネス（ｔｈｉｃｋｎｅｓｓ）、年齢および性別を有している。方法を示すために、４０００個のダミーの遺伝子がこのデータセットに追加され、その結果４００４列と２５５行を有するデータマトリックスが得られた。故意に、４０００個の“遺伝子”は生存時間と関連付けられていない。そのほとんどが情報を持たない４００４個の可能性のある予測子から重要な変数を識別することは、アルゴリズム的に難題である。ハザード関数に対してパラメトリックワイブルモデルを使用して、データが解析された。
【０４１０】
アルゴリズムは変数：年齢に関してのみ選択を行った。擬似遺伝子変数は全て迅速に廃棄された。ワイブル形状パラメータは０．６８と推定された。
【０４１１】
例９：前立腺癌に対する順序付けられたカテゴリー的解析
この例は、Dhanasekaran氏他による文献(2001)からのものである。
【数２１４】
2005524124.tif 000239
【０４１２】
および
【数２１５】
2005524124.tif 000240
【０４１３】
における補足ファイルもまた参照されたい。
【０４１４】
９６０５個の遺伝子と共に１５個のサンプル（個人）が存在する。失われた値は、行平均＋列平均マイナス総（ｇｒａｎｄ）平均によって置換された。４つの順序付けられたカテゴリー（Ｇ＝４）、すなわち、
１．ＮＡＰ正常な
２．ＢＰＨ良性の
３．ＰＣＡ局所化した
４．ＭＥＴ転移した
カテゴリーが存在する。
【０４１５】
アルゴリズムは、１つの誤分類を除いて全ての個人を正しく分類することのできる１つの遺伝子（遺伝子番号６６１１、それらの取得ＩＤ  Ｒ３１６７９）を発見した。
【０４１６】
ＥＭアルゴリズムからの反復は以下のとおりである：
【表３９】
2005524124.tif 000241
【表４０】
2005524124.tif 000242
【表４１】
2005524124.tif 000243
【表４２】
2005524124.tif 000244
【表４３】
2005524124.tif 000245
【表４４】
2005524124.tif 000246
【表４５】
2005524124.tif 000247
【０４１７】
あてはめられた確率のプロットは以下の図６に示されている。ラインは以下のようにクラスを表している：破線＝クラス１、実線＝クラス２、点線＝クラス３、点線と破線を組合せた線＝クラス４。観察（インデックス）１乃至３はクラス２中にあり、観察４乃至７はクラス１中にあり、観察８乃至１１はクラス３中にあり、観察１２乃至１５はクラス４中にある。
【０４１８】
例１０：前立腺癌−選択された遺伝子に対する順序付けられたカテゴリー的解析
この例は、データセットが５０個の選択された遺伝子に減少されていることを除いて、例９のものと同一である。これらの遺伝子の１つは例９において重要であると認められた遺伝子であり、その他のものはランダムに選択された。この例の目的は、完全に表にされたデータセット（表１０）に基づいて例証を行うことである。
【０４１９】
失われた値は、行平均＋列平均マイナス総（ｇｒａｎｄ）平均によって置換された。４つの順序付けられたカテゴリー（Ｇ＝４）、すなわち、
１．ＮＡＰ正常な
２．ＢＰＨ良性の
３．ＰＣＡ局所化した
４．ＭＥＴ転移した
カテゴリーが存在する。
【０４２０】
アルゴリズムは、例９の遺伝子６６１１（取得  Ｒ３１６７９）に等価である１つの予測遺伝子（表１０の遺伝子１）を発見した。予測の成功は、当然ながら、例９のものと同じであった（それは同じ単一の遺伝子に基づいていたので）。
【表４６】
2005524124.tif 000248
【表４７】
2005524124.tif 000249
【０４２１】
例１１：方法を使用するための装置
図５を参照とすると、本発明の実施形態による方法を実施するのに適したパーソナルコンピュータ20が示されている。コンピュータ20は、ハードディスクデータ記憶装置21上に記憶されたソフトウェアタプログラムの命令の下で動作する。コンピュータ20はさらに、プロセッサ22、メモリ23、表示スクリーン24、プリンタ25、ならびにマウス26およびキーボード27のような入力装置を備えている。コンピュータはインターネット28へのネットワーク接続27のような通信手段、またはデータのダウンロードまたは収集および共用を容易にするデータ収集手段28を有していてもよい。
【０４２２】
データ収集手段はデータをシステムから収集し、あるいはダウンロードする。コンピュータはソフトウェアに埋込まれた操作手段を備えており、この操作手段はマウス26およびキーボード27と通信して、ユーザが本発明の実施形態による方法をデータに関して実施することを可能にする。システムは、本発明の実施形態による方法を実施するソフトウェアに埋込まれた手段と、およびグラフィックを生成する手段とを備えている。この方法が実施された後、出力は表示スクリーン24上にグラフィック的に示され、および、またはプリンタ25でプリントされることができる。
【０４２３】
上記の例では、生物学的システムに関して本発明の構成を説明してきた。前に述べたように、本発明は、サンプルの特徴が予測されることを必要とする任意の“システム”に適用されることができる。システムの例には、化学的システム、農業システム、気象システム、たとえば、クレジットリスク査定システムを含む金融システム、保険システム、マーケティングシステム、または企業記録システム、電子システム、物理的システム、天文物理学的システムおよび機械的システムが含まれる。
【０４２４】
当業者に明らかであるように、修正および変形は本発明の技術的範囲内であると考えられる。
【０４２５】
参照文献
【表４８】
2005524124.tif 000250
【表４９】
2005524124.tif 000251
【図面の簡単な説明】
【０４２６】
【図１】本発明の１実施形態の予測の成功に関する置換テストの結果を示すヒストグラム。
【図２】本発明の１実施形態の予測の成功に関する置換テストの結果を示すヒストグラム。
【図３】本発明の方法の１実施形態において使用される一般化された線形モデルに対する曲線のプロット。
【図４】本発明の方法の１実施形態を使用して識別された単一の遺伝子に対するあてはめられた確率のプロット。
【図５】本発明によるシステムを実施するために使用されるパーソナルコンピュータの概略図。
【図１】
2005524124.tif 000267
【図２】
2005524124.tif 000268
【図３】
2005524124.tif 000269
【図４】
2005524124.tif 000270
【図５】
2005524124.tif 000271
【国際調査報告】
2005524124.tif 000272
2005524124.tif 000273