(57)【要約】【課題】  隣接カテゴリーロジットモデルのパラメータ推定において、サンプル数が少ないときでも精度良く安定な推定法、システム及びプログラムを提供する。【解決手段】  パラメータ推定部２３は、入力装置１から与えられたデータに対して罰則付対数尤度関数ｌｐ(θ)を最大にするパラメータ［外２４］を推定し、推定されたパラメータ［外２４］を事後確率計算部２２へ送る。事後確率計算部は入力装置から与えられたデータに対して、パラメータ推定部から与えられたパラメータ［外２４］を用いてｙｉの事後確率を計算し、計算された事後確率を出力装置３へ送る。［外２４］2006024040.tif 000057【選択図】      図３
2006024040.tif 000002
【特許請求の範囲】
【請求項１】
入力データを入力するための入力装置と、該入力装置に入力された前記入力データを用いてパラメータを推定するデータ処理装置とを備えた隣接カテゴリーロジットモデル解析システムを用いたパラメータの推定方法において、
前記データ処理装置が罰則付最尤法によりパラメータを推定することを特徴とする隣接カテゴリーロジットモデルのパラメータ推定方法。
【請求項２】
請求項１に記載のパラメータ推定方法において、
順序カテゴリーデータとしてＹ∈｛０，・・・，Ｋ｝、共変量としてｘ＝(１，ｘ１，・・・，ｘｐ)Ｔを与えたときの応答変数Ｙの値が０となる確率が数１で定義され、
【数１】
2006024040.tif 000048
共変量ｘが与えられた状況下での応答変数Ｙの条件付確率が数２で定義され、
【数２】
2006024040.tif 000049
前記罰則付最尤法が、数３により表される尤度関数を用いてパラメータθを推定する方法である、
【数３】
2006024040.tif 000050
ことを特徴とするパラメータ推定方法。
【請求項３】
入力データを入力するための入力装置と、該入力装置に入力された前記入力データを用いてパラメータを推定するデータ処理装置とを備えた隣接カテゴリーロジットモデル解析システムにおいて、
前記データ処理装置が罰則付最尤法により前記パラメータを推定するパラメータ推定部を備えていることを特徴とする隣接カテゴリーロジットモデル解析システム。
【請求項４】
請求項３に記載の隣接カテゴリーロジットモデル解析システムにおいて、
前記パラメータ推定部内に、
順序カテゴリーデータとしてＹ∈｛０，・・・，Ｋ｝、共変量としてｘ＝(１，ｘ１，・・・，ｘｐ)Ｔを与えたときの応答変数Ｙの値が０となる確率を数４で定義し、
【数４】
2006024040.tif 000051
かつ、共変量ｘが与えられた状況下での応答変数Ｙの条件付確率を数５で定義し、
【数５】
2006024040.tif 000052
当該パラメータ推定部が、前記罰則付最尤法として、数６により表される尤度関数を用いてパラメータθを推定する、
【数６】
2006024040.tif 000053
ことを特徴とする隣接カテゴリーロジットモデル解析システム。
【請求項５】
請求項３又は４に記載された隣接カテゴリーロジットモデル解析システムにおいて、
前記パラメータ推定部が推定したパラメータを利用して、前記入力データについて、事後確率計算を行う事後確率計算部をさらに備えていることを特徴とする隣接カテゴリーロジットモデル解析システム。
【請求項６】
請求項３又は４に記載された隣接カテゴリーロジットモデル解析システムにおいて、
前記パラメータ推定部が推定したパラメータを外部出力する出力装置をさらに備えていることを特徴とする隣接カテゴリーロジットモデル解析システム。
【請求項７】
請求項５に記載された隣接カテゴリーロジットモデル解析システムにおいて、
前記パラメータ推定部が推定したパラメータ及び前記事後確率計算部が計算した事後確率を外部出力する出力装置をさらに備えていることを特徴とする隣接カテゴリーロジットモデル解析システム。
【請求項８】
コンピュータに隣接カテゴリーロジットモデルについてパラメータを推定させるプログラムにおいて、
罰則付最尤法により前記パラメータを推定するステップをコンピュータに実行させるプログラム。
【請求項９】
請求項８に記載のプログラムにおいて、
順序カテゴリーデータとしてＹ∈｛０，・・・，Ｋ｝、共変量としてｘ＝(１，ｘ１，・・・，ｘｐ)Ｔを与えたときの応答変数Ｙの値が０となる確率を数７で定義し、
【数７】
2006024040.tif 000054
共変量ｘが与えられた状況下での応答変数Ｙの条件付確率を数８で定義し、
【数８】
2006024040.tif 000055
数９により表される罰則付尤度関数を用いてパラメータθを推定する、
【数９】
2006024040.tif 000056
ステップをコンピュータに実行させることを特徴とするプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、隣接カテゴリーロジットモデルにおけるパラメータ推定方法、システム及びプログラムに関し、特に、パラメータを適切に推定し、サンプル数が少ないときでも安定した推定値を与えることができる、隣接カテゴリーロジットモデルにおけるパラメータの推定方法、システム及びプログラムに関する。
【背景技術】
【０００２】
順序カテゴリーデータ（ordered categorical data）は、順序をもったカテゴリーデータであって、医学や生物学あるいは社会学などにおいて非常に良く扱われるデータである。順序カテゴリーデータには、例えば、薬の臨床試験における症状を示す「悪化」、「不変」、「軽度改良」、「中等度改良」及び「著明改良」や、副作用を示す「グレード（Grade）０」、「グレード１」、「グレード２」、「グレード３」、「グレード４」及び「グレード５」などがある。これらは、まとめて多値反応データ（multiple response data）と呼ばれることもある。これに対して、血液型である「Ａ型」、「Ｂ型」、「ＡＢ型」及び「Ｏ型」はデータ間に順序関係がないので、順序カテゴリーデータではない。
【０００３】
順序カテゴリーデータに対する解析手法として、複数のカテゴリーを統合したり興味あるカテゴリーだけを抽出したりして、２カテゴリーのデータに作り直し、ロジスティック回帰分析が適用される。例えば、「グレード０」、「グレード１」及び「グレード２」を新しく[カテゴリー０]とし、「グレード３」、「グレード４」及び「グレード５」を新しく[カテゴリー１]とする。しかし、複数のカテゴリーを２カテゴリーに統合するので、本来データが持っている順序の情報が失われてしまう。
【０００４】
一方、順序カテゴリーデータを２カテゴリーに統合することなく解析できる統計モデルとして、隣接カテゴリーロジット（adjacent categories logit：ＡＣＬ）モデルがある（例えば、非特許文献１参照。）。ＡＣＬモデルは、順序カテゴリーデータに対する統計モデルとして医学や生物学あるいは社会学などにおいて一般的に広く使用されている。
【０００５】
以下、ＡＣＬモデルについて説明する。
【０００６】
順序カテゴリーデータとしてＹ∈｛０，・・・，Ｋ｝、共変量としてｘ＝(１，ｘ１，・・・，ｘｐ)Ｔを考える。このとき、応答変数Ｙの値が０となる確率を数１０で定義する。
【数１０】
2006024040.tif 000003
【０００７】
ここで、θ＝(γＴ，βＴ)Ｔ、γ＝(γ１，・・・，γＫ−１)Ｔおよびβ＝(β０，・・・，βｐ)Ｔとする。
【０００８】
いま、共変量ｘが与えられた状況下での応答変数Ｙの条件付確率を数１１で定義する。
【数１１】
2006024040.tif 000004
【０００９】
このＡＣＬモデルのパラメータθの推定法としては、例えば最尤法（以下、Ordinary likelihood method：ＯＬ法という）が用いられる。以下、図１及び図２を参照してＯＬ法を用いたＡＣＬモデルのパラメータ推定手順、及びデータ解析ついて説明する。
【００１０】
図１に従来の隣接カテゴリーロジットモデル解析システムを示す。図示の解析システムは、入力装置１、データ処理装置２、及び出力装置３を備えている。また、データ処理装置２は、通常の最尤法を用いるパラメータ推定部２１と事後確率計算部２２とを備えている。このデータ解析システムは、図２のフローチャートに示すように動作する。
【００１１】
はじめに、入力装置１からデータ（ｙ，ｘ）を入力する（ステップＢ１）。次に、入力装置１に入力されたデータ（ｙ，ｘ）に対して、データ処理装置２の通常の最尤法を用いるパラメータ推定部２１において、数１２で示す対数尤度関数ｌ０(θ)を最大化するパラメータ［外１］を求める（ステップＢ２）。
【数１２】
2006024040.tif 000005
【００１２】
［外１］
2006024040.tif 000006
【００１３】
次に、事後確率計算部２２において、推定されたパラメータ［外２］及び入力装置１により入力されたデータ（ｙｉ，ｘｉ）を用いて数１及び数２よりｙｉに対する事後確率Ａ(ｙｉ｜ｘｉ，θ)を計算する（ステップＢ３）。
【００１４】
［外２］
2006024040.tif 000007
【００１５】
全てのデータｙｉの事後確率を計算したかどうかを判定し、全てのｙｉに対する事後確率を計算していれば終了し、計算していなければステップＢ３に戻る（ステップＢ４）。
【００１６】
出力装置３においては、回帰係数であるβやパラメータγを出力したり、各データ（ｙｉ，ｘｉ）に対する事後確率Ａ(ｙｉ｜ｘｉ，θ)を出力したりする（ステップＢ５）。回帰係数βより目的変数ｙに対する説明変数ｘ＝(１，ｘ１，・・・，ｘｐ)Ｔの寄与の度合いを判断することができる。
【００１７】
ここで、ＯＬ法を用いたパラメータ推定には、不適切な推定値［外３］を与えたり、サンプル数が少ない場合には収束が不安定となる等の問題点があることが良く知られている。そこで、より安定してパラメータを推定する方法として、条件付最尤法（Conditional maximum likelihood method、以下ＣＭＬ法という）を用いる推定法が提案されている（例えば、非特許文献２及び３参照。）。しかしながら、ＣＭＬ法を用いたパラメータ推定は、ＡＣＬモデルの２カテゴリーの場合とみなせるロジスティック回帰モデルにおいてサンプルサイズが大きいときに計算量が膨大となり、実行が困難であることが知られている（非特許文献４及び５参照。）。
【００１８】
［外３］
2006024040.tif 000008
【００１９】
【非特許文献１】Agresti, A. (1984) Analysis of Ordinal Categorical Data, New York: John Wiley.
【非特許文献２】Hirji, K.F (1992) Computing exact distribution for polytomous response data, Journal of American Statistical Society, 87, 487-492
【非特許文献３】Agresti, A. (1999) Modeling ordered categorical data: recent advantages and future challenges, Statistics in Medicine, 18, 2191-2207.
【非特許文献４】Bull, S. B., et al, (1997). Jackknife bias reduction for polychotomous logistic regression. Statistics in Medicine, 16, 545-560.
【非特許文献５】Mehta, C. R. et al, (2000) Efficient Monte Carlo methods for conditional logistic regression, Journal of the American Statistical Association, 95, 99-108.
【発明の開示】
【発明が解決しようとする課題】
【００２０】
上述したようにＯＬ法を用いたＡＣＬモデルのパラメータ推定は、精度が悪く、特にサンプル数が少ないときにパラメータの推定が不安定になるという問題点がある。
【００２１】
また、ＣＭＬ法を用いたＡＣＬモデルのパラメータ推定は、計算量が非常に膨大になる場合があるという問題点がある。
【００２２】
本発明の目的は、ＡＣＬモデルのパラメータθの推定において、従来の推定法よりも精度良く、またサンプル数が少ないときにも解の収束が安定するようなパラメータ推定方法、システム及びプログラムを提供することである。
【００２３】
また、本発明の他の目的は、ＣＭＬ法を用いたパラメータ推定よりも計算量が少ないパラメータ推定方法、システム及びプログラムを提供することである。
【課題を解決するための手段】
【００２４】
本発明は、第１の要旨によれば、入力データを入力するための入力装置と、該入力装置に入力された前記入力データを用いてパラメータを推定するデータ処理装置とを備えた隣接カテゴリーロジットモデル解析システムを用いたパラメータの推定方法において、前記データ処理装置が罰則付最尤法によりパラメータを推定することを特徴とする。
【００２５】
具体的には、上記パラメータ推定方法において、順序カテゴリーデータとしてＹ∈｛０，・・・，Ｋ｝、共変量としてｘ＝(１，ｘ１，・・・，ｘｐ)Ｔを与えたときの応答変数Ｙの値が０となる確率が数１３で定義され、
【数１３】
2006024040.tif 000009
【００２６】
共変量ｘが与えられた状況下での応答変数Ｙの条件付確率が数１４で定義され、
【数１４】
2006024040.tif 000010
【００２７】
前記罰則付最尤法が、数１５により表される尤度関数を用いてパラメータθを推定する方法である、
【数１５】
2006024040.tif 000011
【００２８】
ことを特徴とする。
【００２９】
また、本発明は、第２の要旨によれば、入力データを入力するための入力装置と、該入力装置に入力された前記入力データを用いてパラメータを推定するデータ処理装置とを備えた隣接カテゴリーロジットモデル解析システムにおいて、前記データ処理装置が罰則付最尤法により前記パラメータを推定するパラメータ推定部を備えていることを特徴とする。
【００３０】
具体的には、上記隣接カテゴリーロジットモデル解析システムにおいて、前記パラメータ推定部内に、順序カテゴリーデータとしてＹ∈｛０，・・・，Ｋ｝、共変量としてｘ＝(１，ｘ１，・・・，ｘｐ)Ｔを与えたときの応答変数Ｙの値が０となる確率を数１６で定義し、
【数１６】
2006024040.tif 000012
【００３１】
かつ、共変量ｘが与えられた状況下での応答変数Ｙの条件付確率を数１７で定義し、
【数１７】
2006024040.tif 000013
【００３２】
当該パラメータ推定部が、前記罰則付最尤法として、数１８により表される尤度関数を用いてパラメータθを推定する、
【数１８】
2006024040.tif 000014
【００３３】
ことを特徴とする。
【００３４】
さらに、本発明は、第３の要旨によれば、コンピュータに隣接カテゴリーロジットモデルについてパラメータを推定させるプログラムにおいて、罰則付最尤法により前記パラメータを推定するステップをコンピュータに実行させることを特徴とする。
【００３５】
具体的には、上記プログラムにおいて、順序カテゴリーデータとしてＹ∈｛０，・・・，Ｋ｝、共変量としてｘ＝(１，ｘ１，・・・，ｘｐ)Ｔを与えたときの応答変数Ｙの値が０となる確率を数１９で定義し、
【数１９】
2006024040.tif 000015
【００３６】
共変量ｘが与えられた状況下での応答変数Ｙの条件付確率を数２０で定義し、
【数２０】
2006024040.tif 000016
【００３７】
数２１により表される罰則付尤度関数を用いてパラメータθを推定する、
【数２１】
2006024040.tif 000017
【００３８】
ステップをコンピュータに実行させることを特徴とする。
【発明の効果】
【００３９】
本発明の第１の効果は、隣接カテゴリーロジットモデルにおけるパラメータ推定において、対数尤度関数の代わりに罰則付尤度関数を用いてパラメータθの推定を行うようにしたことで、パラメータθを精度良く、またサンプル数が少ないときにも安定して推定することができることである。
【００４０】
本発明の第２の効果は、隣接カテゴリーロジットモデルにおけるパラメータ推定において、罰則付最尤法を用いていることにより、通常の最尤法と同じ計算量でパラメータを推定することができること、つまり、条件付最尤法よりも少ない計算量でパラメータを推定することができることである。
【発明を実施するための最良の形態】
【００４１】
以下、本発明の実施の形態について図面を参照して詳細に説明する。
【００４２】
図３を参照すると、本発明の一実施の形態に係る隣接カテゴリーロジットモデル解析システムは、キーボード等の入力装置１と、プログラム制御により動作するデータ処理装置２と、ディスプレイ装置や印刷装置等の出力装置３とを含む。
【００４３】
データ処理装置２は、罰則付最尤法を用いるパラメータ推定部２３と、事後確率計算部２２とを備えている。データ処理装置２は、例えば、コンピュータプログラムとそれを実行するコンピュータによって実現される。コンピュータプログラムは、フレキシブルディスクやハードディスク等の磁気ディスク、ＣＤ−ＲＯＭやＭＯ等の光ディスク、あるいは半導体メモリなど、コンピュータ読み取り可能な記録媒体により外部から提供されてもよいし、コンピュータに内蔵する記憶装置に記憶されていてもよい。あるいは、コンピュータにネットワーク接続された外部記憶装置や他のコンピュータ内の記憶装置に記憶されていてもよい。
【００４４】
このシステムの解析対象である隣接カテゴリーロジットモデルは、以下のように定義される。
【００４５】
順序カテゴリーデータとしてＹ∈｛０，・・・，Ｋ｝、共変量としてｘ＝(１，ｘ１，・・・，ｘｐ)Ｔが与えられるとき、応答変数Ｙの値が０となる確率を数２２で定義する。
【数２２】
2006024040.tif 000018
【００４６】
ここで、θ＝(γＴ，βＴ)Ｔ、γ＝(γ１，・・・，γＫ−１)Ｔおよびβ＝(β０，・・・，βｐ)Ｔである。
【００４７】
いま、共変量ｘが与えられた状況下での応答変数Ｙの条件付確率を数２３で定義する。
【数２３】
2006024040.tif 000019
【００４８】
本実施の形態に係る解析システムは、数２４に示す罰則付尤度関数を用いて、ＡＣＬモデルに対するパラメータθの推定を行う。なお、このような罰則付尤度関数を用いた最尤法を、以下、ＰＬ法（Penalized Likelihood method）という。
【数２４】
2006024040.tif 000020
【００４９】
ここで、ξは罰則パラメータ（非負の定数）である。なお、ξ＝０のとき、数２４は数１２に等しい。なお、基準化した共変数を用いる場合、上記罰則項の代わりに数２４−１を使用すると、解の偏りは大きくなるが、収束性がさらによくなる。
【数２４−１】
2006024040.tif 000021
以下、図３の解析システムの動作を図４を参照して説明する。
【００５０】
まず、入力装置１から順序カテゴリーデータ及び共変量が入力されると（ステップＣ１）、罰則付最尤法を用いるパラメータ推定部２３は、罰則パラメータξを初期化する（ステップＣ２）。
【００５１】
次に、パラメータ推定部２３は、入力装置１から与えられたデータに対して、上記数２４で示す罰則付対数尤度関数ｌｐ(θ)を最大にするパラメータ［外４］を推定する（ステップＣ３）。推定されたパラメータ［外４］は、事後確率計算部２２及び出力装置３へ送られる。
【００５２】
［外４］
2006024040.tif 000022
事後確率計算部２２は、入力装置１から与えられたデータに対して、罰則付最尤法を用いるパラメータ推定部２３から与えられたパラメータ［外５］を用いて、ｙｉの事後確率を計算する（ステップＣ４）。計算された事後確率は出力装置３へ送られる。
【００５３】
［外５］
2006024040.tif 000023
全てのデータｙｉの事後確率を計算したかどうかを判定し、全てのｙｉに対する事後確率を計算していれば終了し、計算していなければステップＣ４に戻る（ステップＣ５）。
【００５４】
出力装置３においては、回帰係数であるβやパラメータγを出力したり、各データ（ｙｉ，ｘｉ）に対する事後確率Ａ(ｙｉ｜ｘｉ，θ)を出力したりする（ステップＣ６）。
【００５５】
以上のようにして得られた回帰係数βより、目的変数ｙに対する説明変数ｘ＝(１，ｘ１，・・・，ｘｐ)Ｔの寄与の度合いを判断することができる。ここで、入力データとして遺伝子データ（説明変数ｘ）と臨床データ（目的変数ｙ）とが与えられたとする。遺伝子データは、例えば、マイクロアレイなどによって得られる遺伝子発現データである。また、臨床データは、例えば、副作用や薬効などの順序カテゴリーデータである。この場合、出力装置３から出力される回帰係数の推定値［外６］より、副作用や薬効などの臨床データに対する各遺伝子の寄与の度合いを判断することができる。つまり、［外６］が正となる遺伝子に対しては、副作用や薬効などの効果を増加させる働きがあることが予想できる。また、［外６］の絶対値の大きさの大小により、各遺伝子の副作用や薬効などに対する寄与の度合いを比較することができる。
【００５６】
［外６］
2006024040.tif 000024
【実施例１】
【００５７】
次に、本発明の実施例を、シミュレーションの結果を参照して具体的に説明する。ここでは、ＯＬ法とＰＬ法とをそれぞれ用いたパラメータ推定の性能を比較する。なお、本実施例におけるシミュレーションにおいてはｘ〜Ｎ３(０，Ｉ３)とする。ただし、Ｉ３は数２５で表され、Ｎ３(０，Ｉ３)は３次元の標準正規分布とする。
【数２５】
2006024040.tif 000025
【００５８】
（シミュレーション１）
シミュレーション１においては、数２６で表される確率に従って順序カテゴリーデータｙ∈｛０，１，２｝が生成される。
【数２６】
2006024040.tif 000026
【００５９】
ここで、γ＝１，β＝(０，０，１，−３)Ｔとする。また、Ａ(ｙ｜ｘ，（γ，β））は数２７（＝数２３の特別な場合）で表されるものとする。本モデルはＡＣＬモデルに対応している。
【数２７】
2006024040.tif 000027
【００６０】
（シミュレーション２）
シミュレーション２においては、数２８に従って順序カテゴリーデータｙ∈｛０，１，２｝が生成される。
【数２８】
2006024040.tif 000028
【００６１】
シミュレーション１及び２のそれぞれにおいて、サンプル数（Ｎ）を２０，３０，５０，１００と変化させながら１０００回ずつ繰り返し解析を行う。ただし、各カテゴリーに属するサンプル数が少なくとも５つあるデータセットのみをシミュレーションに使用する。また、ＰＬ法における罰則パラメータξをξ＝１．０とする。
【００６２】
表１は不適切な推定値が得られたり、解が収束しなかった頻度を示した表である。
【表１】
2006024040.tif 000029
【００６３】
ここで不適切な推定値とは、例えば数２９で示されるように、推定されたパラメータの下限が−１０を下回り、又は上限が１０を超えることをいう。
【数２９】
2006024040.tif 000030
【００６４】
表１から、サンプル数及びシミュレーションタイプに関係なくＰＬ法がＯＬ法よりも安定に適切な推定値を与えることが分かる。特にシミュレーション２においては、サンプル数が少ない（ｎ＝２０，３０）場合においては、ＯＬ法では不適切な推定値を与える頻度が３０％を超えるにもかかわらず、ＰＬ法では１％程度の頻度である。
【００６５】
図５乃至図８にシミュレーション１（Ｎ＝５０）において推定された回帰係数［外７］及び［外８］の分布を示す。なお、β０＝０，β１＝０，β２＝１，β３＝−３，γ＝１の直線はそれぞれのパラメータの真の値を示している。図５乃至図８を見て分かるように、罰則の項を尤度関数に導入することによって、ＰＬ法によって推定されたパラメータ［外７］と［外８］の相関がＯＬ法によって推定されたパラメータ［外７］と［外８］の相関よりも小さくなっていることがわかる。即ち、ＰＬ法によって推定されたパラメータの方が、ＯＬ法によって推定されたパラメータよりも安定している。
【００６６】
［外７］
2006024040.tif 000031
【００６７】
［外８］
2006024040.tif 000032
【００６８】
図９乃至図１３にシミュレーション１（Ｎ＝２０，３０，５０，１００）において推定された回帰係数［外９］及び［外１０］の分布を箱型図で示す。図９乃至図１３から分かるように、サンプル数（Ｎ）が大きくなるにつれてＯＬ法、ＰＬ法ともに推定値の分布が安定することがわかる。また、サンプル数が少ない場合、ＯＬ法よりもＰＬ法の方が推定値の分布が安定していることが分かる。特に図１４のＮ＝２０，３０の場合を比較すると、サンプル数が少ない場合においてもＰＬ法の方がより高精度にβ３を推定できることが分かる。
【００６９】
［外９］
2006024040.tif 000033
【００７０】
［外１０］
2006024040.tif 000034
【図面の簡単な説明】
【００７１】
【図１】従来の隣接カテゴリーロジットモデル解析システムの構成を示すブロック図である。
【図２】図１の解析システムの動作を説明するための流れ図である。
【図３】本発明の一実施の形態に係る隣接カテゴリーロジットモデル解析システムの構成を示すブロック図である。
【図４】図３の解析システムの動作を説明するための流れ図である。
【図５】本発明の実施例（シミュレーション１：Ｎ＝５０）により求められた［外１１］と［外１２］の相関を表す相関図である。
【００７２】
［外１１］
2006024040.tif 000035
［外１２］
2006024040.tif 000036
【図６】本発明の実施例（シミュレーション１：Ｎ＝５０）により求められた［外１３］と［外１４］の相関を表す相関図である。
【００７３】
［外１３］
2006024040.tif 000037
［外１４］
2006024040.tif 000038
【図７】本発明の実施例（シミュレーション１：Ｎ＝５０）により求められた［外１５］と［外１６］の相関を表す相関図である。
【００７４】
［外１５］
2006024040.tif 000039
［外１６］
2006024040.tif 000040
【図８】本発明の実施例（シミュレーション１：Ｎ＝５０）により求められた［外１７］と［外１８］の相関を表す相関図である。
【００７５】
［外１７］
2006024040.tif 000041
［外１８］
2006024040.tif 000042
【図９】本発明の実施例（シミュレーション１：Ｎ＝２０，３０，５０，１００）により求められた［外１９］の分布を示す図である。
【００７６】
［外１９］
2006024040.tif 000043
【図１０】本発明の実施例（シミュレーション１：Ｎ＝２０，３０，５０，１００）により求められた［外２０］の分布を示す図である。
【００７７】
［外２０］
2006024040.tif 000044
【図１１】本発明の実施例（シミュレーション１：Ｎ＝２０，３０，５０，１００）により求められた［外２１］の分布を示す図である。
【００７８】
［外２１］
2006024040.tif 000045
【図１２】本発明の実施例（シミュレーション１：Ｎ＝２０，３０，５０，１００）により求められた［外２２］の分布を示す図である。
【００７９】
［外２２］
2006024040.tif 000046
【図１３】本発明の実施例（シミュレーション１：Ｎ＝２０，３０，５０，１００）により求められた［外２３］の分布を示す図である。
【００８０】
［外２３］
2006024040.tif 000047
【符号の説明】
【００８１】
１    入力装置
２    データ処理装置
３    出力装置
２１    通常の最尤法を用いるパラメータ推定部
２２    事後確率計算部
２３    罰則付最尤法を用いるパラメータ推定部
【図１】
2006024040.tif 000058
【図２】
2006024040.tif 000059
【図３】
2006024040.tif 000060
【図４】
2006024040.tif 000061
【図５】
2006024040.tif 000062
【図６】
2006024040.tif 000063
【図７】
2006024040.tif 000064
【図８】
2006024040.tif 000065
【図９】
2006024040.tif 000066
【図１０】
2006024040.tif 000067
【図１１】
2006024040.tif 000068
【図１２】
2006024040.tif 000069
【図１３】
2006024040.tif 000070