(57)【特許請求の範囲】
【請求項１】
時刻（ｔ−１）での真の市場価格及び時刻（ｔ−１）までのＭ（Ｍは、正の整数）個の真の市場価格の移動平均である中心価格によって定義されるポテンシャル項と、変動誤差と、の和によって、時刻（ｔ）での真の市場価格Ｐ（ｔ）を算出するＰＵＣＫモデルに、粒子フィルタ法を適用する時系列データの変化点検出方法であって、
確率密度関数生成部に、それぞれ異なる値を有するＰＵＣＫモデルの状態を表すパラメータを持つ粒子群より次の時刻のパラメータの確率密度関数を得させる第１のステップと、
適合度評価部に、複数の粒子のそれぞれについて、時刻（ｔ）で観測した市場価格に対する前記時刻（ｔ）での真の市場価格の適合度を評価させる第２のステップと、
再サンプリング部に、前記適合度に応じて、前記複数の粒子の中から粒子を再サンプリングさせる第３のステップと、をコンピュータに実行させ、
前記第３のステップでは、
乱数を生成し、
第１の所定の値と前記乱数とを比較し、
前記乱数が前記第１の所定の値よりも大きい場合には、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする確率密度分布を成分とする確率密度関数を粒子として生成し、
前記乱数が前記第１の所定の値よりも小さい場合には、一様分布を成分とする確率密度関数を粒子として生成する、
時系列データの変化点検出方法。
【請求項２】
時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする確率密度分布は、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする正規分布である、
請求項１に記載の時系列データの変化点検出方法。
【請求項３】
前記第１のステップでは、ＰＵＣＫモデルの状態を表すパラメータを持つ粒子群より次の時刻のＰＵＣＫモデルのパラメータが更新される、
請求項１又は２に記載の時系列データの変化点検出方法。
【請求項４】
前記第３のステップでは、
それぞれの粒子に関して乱数を生成し、
前記第１の所定の値と前記乱数とを比較し、
前記乱数が前記第１の所定の値よりも大きい場合には、時刻（ｔ）におけるＭを前提とする条件付き確率密度関数を粒子として生成し、
前記乱数が前記第１の所定の値よりも小さい場合には、一様分布を成分とする確率密度関数を粒子として生成する、
請求項１乃至３のいずれか一項に記載の時系列データの変化点検出方法。
【請求項５】
時刻（ｔ＋１）での真の市場価格は、前記ポテンシャル項と、時刻（ｔ）での真の市場価格の変動誤差と、を加算することにより算出されることを特徴とする、
請求項１乃至４のいずれか一項に記載の時系列データの変化点検出方法。
【請求項６】
前記時刻（ｔ）で観測した市場価格は、前記時刻（ｔ）での真の市場価格と時刻（ｔ）での市場価格の観測誤差との和であることを特徴とする、
請求項１乃至５のいずれか一項に記載の時系列データの変化点検出方法。
【請求項７】
前記時刻（ｔ）での真の市場価格の変動誤差及び前記時刻（ｔ）での市場価格の観測誤差は、正規分布に従う確率密度関数で与えられることを特徴とする、
請求項５又は６に記載の時系列データの変化点検出方法。
【請求項８】
請求項１乃至７のいずれか一項に記載した時系列データの変化点検出方法により算出した時刻（ｔ）における前記ポテンシャル項のポテンシャル係数、前記Ｍの値、及び時刻（ｔ）での真の市場価格の変動誤差と、に基づいて、前記時刻（ｔ）での真の市場価格Ｐ（ｔ）を時間発展させることにより、時刻（ｔ＋Ｎ）での真の市場価格Ｐ（ｔ＋Ｎ）を算出する、
未来の時系列データ値の確率密度分布の予測方法。
【請求項９】
前記粒子の総数をＮｐ、
前記粒子を特定する番号をｊ（ただし、ｊは１≦ｊ≦Ｎｐを満たす整数）、
番号ｊの粒子の時刻（ｔ）でのポテンシャル係数をｂ（ｊ）ｉ（ｔ）（ただし、ｉはポテンシャルの次数を示す２以上の整数）、
前記番号ｊの粒子の時刻（ｔ）での真の市場価格の変動誤差をｆ（ｊ）Ｐ（ｔ）として、
前記番号ｊの粒子の時刻（ｔ）での真の市場価格をＰ（ｊ）（ｔ）、
前記番号ｊの粒子の時刻（ｔ＋Ｎ）での真の市場価格をＰ（ｊ）（ｔ＋Ｎ）、
前記番号ｊの粒子の時刻（ｔ）での市場価格の中心価格をＰ（ｊ）Ｍ（ｔ）、
時刻（ｔ）におけるＭの値をＭ（ｔ）とした場合に、
前記番号ｊの粒子の前記時刻（ｔ＋Ｎ）での真の市場価格Ｐ（ｊ）（ｔ＋Ｎ）は、Ｋを２以上の整数として、
【数６】
0006004284.tif 000045
で表される、
請求項８に記載の未来の時系列データ値の確率密度分布の予測方法。
【請求項１０】
番号ｊの粒子の時刻（ｔ＋１）でのポテンシャル係数ｂ（ｊ）ｉ（ｔ＋１）は、前記番号ｊの粒子の時刻（ｔ）でのポテンシャル係数ｂ（ｊ）ｉ（ｔ）を時間発展させることにより算出される、
請求項９に記載の未来の時系列データ値の確率密度分布の予測方法。
【請求項１１】
時刻（ｔ−１）での真の市場価格及び時刻（ｔ−１）までのＭ（Ｍは、正の整数）個の真の市場価格の移動平均である中心価格によって定義されるポテンシャル項と、変動誤差と、の和によって時刻（ｔ）での真の市場価格Ｐ（ｔ）を算出するＰＵＣＫモデルに粒子フィルタ法を適用する時系列データの変化点検出をコンピュータに実行させる時系列データの変化点検出プログラムであって、
確率密度関数生成部に、それぞれ異なる値を有するＰＵＣＫモデルの状態を表すパラメータを持つ粒子群を生成することによりパラメータの確率密度関数を得させる第１の処理と、
適合度評価部に、複数の粒子のそれぞれについて、時刻（ｔ）で観測した市場価格に対する前記時刻（ｔ）での真の市場価格の適合度を評価させる第２の処理と、
再サンプリング部に、前記適合度に応じて、前記複数の粒子の中から粒子を再サンプリングさせる第３の処理と、をコンピュータに実行させ、
前記第３の処理では、
乱数を生成する処理と、
第１の所定の値と前記乱数とを比較する処理と、
前記乱数が前記第１の所定の値よりも大きい場合には、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする確率密度分布を成分とする確率密度関数を粒子として再生成する処理と、
前記乱数が前記第１の所定の値よりも小さい場合には、一様分布にしたがって粒子を生成する処理と、を実行させる、
時系列データの変化点検出プログラム。
【請求項１２】
時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする確率密度分布は、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする正規分布である、
請求項１１に記載の時系列データの変化点検出プログラム。
【請求項１３】
前記第３の処理では、前記適合度が大きい粒子ほど、多くサンプリングされるだけでなく、一定の割合で、実現可能性のある適合度が低い粒子を生成することを特徴とする、
請求項１１又は１２に記載の時系列データの変化点検出プログラム。
【請求項１４】
前記第１の処理では、
乱数を生成する処理と、
時刻（ｔ）におけるＰＵＣＫモデルのパラメータを前提とする条件付き確率密度関数を粒子として生成する処理と、
前記乱数が前記第１の所定の値よりも小さい場合には、一様分布を成分とする確率密度関数を粒子として生成する処理と、を実行させることを特徴とする、
請求項１１乃至１３のいずれか一項に記載の時系列データの変化点検出プログラム。
【請求項１５】
前記時刻（ｔ）での真の市場価格は、前記ポテンシャル項と、時刻（ｔ）での真の市場価格の変動誤差と、を加算することにより算出されることを特徴とする、
請求項１１乃至１４のいずれか一項に記載の時系列データの変化点検出プログラム。
【請求項１６】
前記時刻（ｔ）で観測した市場価格は、前記時刻（ｔ）での真の市場価格と時刻（ｔ）での市場価格の観測誤差との和であることを特徴とする、
請求項１１乃至１５のいずれか一項に記載の時系列データの変化点検出プログラム。
【請求項１７】
前記時刻（ｔ）での真の市場価格の変動誤差及び前記時刻（ｔ）での市場価格の観測誤差は、正規分布に従う確率密度関数で与えられることを特徴とする、
請求項１５又は１６に記載の時系列データの変化点検出プログラム。
【請求項１８】
請求項１１乃至１７のいずれか一項に記載した時系列データの変化点検出プログラムにより算出した時刻（ｔ）における前記ポテンシャル項のポテンシャル係数、前記Ｍの値、及び時刻（ｔ）での真の市場価格の変動誤差と、に基づいて、前記時刻（ｔ）での真の市場価格Ｐ（ｔ）を時間発展させることにより、時刻（ｔ＋Ｎ）での真の市場価格Ｐ（ｔ＋Ｎ）を算出することを特徴とする、
未来の時系列データ値の確率密度分布の予測プログラム。
【請求項１９】
前記粒子の総数をＮｐ、
前記粒子を特定する番号をｊ（ただし、ｊは１≦ｊ≦Ｎｐを満たす整数）、
番号ｊの粒子の時刻（ｔ）でのポテンシャル係数をｂ（ｊ）ｉ（ｔ）（ただし、ｉはポテンシャルの次数を示す２以上の整数）、
前記番号ｊの時刻（ｔ）での真の市場価格の変動誤差をｆ（ｊ）Ｐ（ｔ）として、
前記番号ｊの粒子の時刻（ｔ）での真の市場価格をＰ（ｊ）（ｔ）、
前記番号ｊの粒子の時刻（ｔ＋Ｎ）での真の市場価格をＰ（ｊ）（ｔ＋Ｎ）、
前記番号ｊの粒子の時刻（ｔ）での市場価格の中心価格をＰ（ｊ）Ｍ（ｔ）、
時刻（ｔ）におけるＭの値をＭ（ｔ）とした場合に、
前記番号ｊの粒子の前記時刻（ｔ＋Ｎ）での真の市場価格Ｐ（ｊ）（ｔ＋Ｎ）は、Ｋを２以上の整数として、
【数７】
0006004284.tif 000046
で表される、
請求項１８に記載の未来の時系列データ値の確率密度分布の予測プログラム。
【請求項２０】
番号ｊの粒子の時刻（ｔ＋１）でのポテンシャル係数ｂ（ｊ）ｉ（ｔ＋１）は、前記番号ｊの粒子の時刻（ｔ）でのポテンシャル係数ｂ（ｊ）ｉ（ｔ）を時間発展させることにより算出される、
請求項１９に記載の未来の時系列データ値の確率密度分布の予測プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は時系列データの変化点検出方法及びプログラム、未来の時系列データ値の確率密度分布予測方法及びプログラムに関する。
【背景技術】
【０００２】
今日、様々な分野において、各種の時系列データが用いられている。例えば、国際金融における外国為替の変動や、証券取引における株価などが時系列データとして用いられている。他にも、例えば半導体製造などの製造業において、歩留まりトレンドの管理などが行われる。【０００３】
こうした時系列データは、一般に、トレンドの変化点を検出し、その変化に応じた行動及び対策を取るために用いられる。例えば為替取引における円／ドル取引においては、値動きの動向、特に大きくトレンドが変化する点を検出し、売買方針を適宜修正することが求められる。こうした時系列データでは、過去の時系列データの移動平均を用いて長期的なトレンドを見積もることが行われる。身近な例として、円／ドル取引では、為替チャートに２５日間や１３週間などの移動平均線が表示される。また、こうした時系列データの長期トレンド分析においては、移動平均を中心とするランダムウォークを適用する手法が広く用いられている。【０００４】
ところが、為替取引などでは、値動きに連動した売り注文又は買い注文が集中するなどの事象が生じることで、短期的に急激な値動きを生じる場合が有る。この場合、単に過去のデータの移動平均を中心とするランダムウォークのみでは、こうした急峻なトレンドの変化を検出することは困難である。つまり、移動平均は、過去の複数のデータについての平均であるので、直近のデータが大きな変動を示したとしても、母集団に比べて直近のデータ数が少なければ、ランダムウォークを適用しても急峻なトレンド変化を捉えることができない。【０００５】
上述の問題を解決するため、発明者らは、時系列データの急峻な変化を検出するＰＵＣＫ（Potentials of Unbalanced Complex Kinetics）モデルを提案した（非特許文献１）。ＰＵＣＫモデルでは、金融等の分野における時系列データの変動を、時刻（ｔ）での真の市場価格（Noise-Reduced Price）Ｐ（ｔ）、時刻（ｔ）での市場価格の観測誤差（Observation Noise）ｆＹ（ｔ）、時刻（ｔ）で観測した市場価格（best-bid, best-ask,  Mid-quote Data等）Ｙ（ｔ）の３要素により、以下の式（１）で表される。
【数１】
0006004284.tif 000002
【０００６】
ティックとは、取引市場において、約定成立による取引対象物の価格の変動が生じるイベントのことをいう。よって、１ティック間とは、ある時点での約定タイミングと、それよりも１つ前の約定タイミングとの間の時間間隔を意味する。また、ティック数とは、ある期間内における約定成立回数、すなわち取引対象物の価格の変動が生じるイベントの発生回数を意味する。以降（ｔ）は、ティック数で数えた時刻または実時間を表す。ＰＵＣＫモデルでは、数ティックスケールの変位の逆相関性を除去した後の真の市場価格Ｐ（ｔ）の変動を取り扱う。１ティックの間での真の市場価格の変動は、以下の式（２）で表すことができる。
【数２】
0006004284.tif 000003
但し、式（２）の右辺第１項は、１ティックの間での真の市場価格の変動に作用するポテンシャルの寄与分を表し、右辺第２項のｆＰ（ｔ）は、時刻（ｔ）での真の市場価格を変動させる変動誤差を表す。ｂ（ｔ）は時刻（ｔ）でのポテンシャル係数、Ｍは価格変動の中心を見積もるために必要な直近のデータ数、ＰＭ（ｔ）は直近Ｍ個の真の市場価格を用いて見積もった価格変動の中心価格である。また、式（２）における時刻（ｔ）でのポテンシャル係数ｂ（ｔ）は、直近のデータ数Ｍに対して依存性を有し、｛１／（Ｍ−１）｝に比例して減少することが経験的に知られている。よって、式（２）では、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）の直近のデータ数Ｍへの依存性を避けるため、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）に｛１／（Ｍ−１）｝を乗じている。【０００７】
中心価格ＰＭ（ｔ）は、以下の式（３）で表される。
【数３】
0006004284.tif 000004
【０００８】
さらに、式（２）は、時々刻々変化する線形な中心力が働く場の中でのランダムウォークとみなすと、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）を２次以上のポテンシャルの係数と捉え、式（２）の右辺第１項が表すポテンシャルの寄与分にかかるポテンシャルを、ポテンシャルＷ（ｑ,ｔ）として、以下の式（４）で表わすことができる。但し、ｉはポテンシャルの次数を示す。
【数４】
0006004284.tif 000005
価格変動の平均速度の絶対値は、ポテンシャル関数の接線勾配の大きさに比例する。価格変動の期待される方向は、勾配を下る方向である。したがって、ポテンシャル関数のどこに今の市場価格があるのかを理解するためには、中心価格と真の市場価格の乖離を計算し、ポテンシャル関数のどこに市場価格があるのかを理解する必要がある。
さらに、式（４）を用いると、式（２）を以下の式（５）に書き換えることができる。
【数５】
0006004284.tif 000006
【０００９】
上述のように、式（４）及び（５）において、時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ）の値により、ポテンシャルＷ（ｑ,ｔ）の寄与度が決定される。例えば、時刻（ｔ）で各次のポテンシャル係数ｂｉ（ｔ）＝０であれば、Ｐ（ｔ＋１）−Ｐ（ｔ）＝ｆＰ（ｔ）となり、変動誤差の成分のみが価格差に寄与する。一方、時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）の絶対値が大きくなるにつれて、ポテンシャルＷ（ｑ,ｔ）の寄与が大きくなる。その結果、ランダムウォークモデルで計算される拡散係数よりも速い拡散、遅い拡散が時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）の関数として定量的に推定できることが知られている（以下、異常拡散と称する）。特に、２次ポテンシャルの係数ｂ２（ｔ）≦−２の場合には、市場価格の変動は、もはや過去のデータ移動平均から逸脱し、発散してしまうことが知られている。【００１０】
このように、ＰＵＣＫモデルでは、２次以上のポテンシャル項を考える事により、上述のような短期間の異常拡散を表現することが可能となる。これにより、これまでのランダムウォークのみを用いた手法では実現できない、異常拡散の検出を行うことができる。
【先行技術文献】
【非特許文献】
【００１１】
【非特許文献１】Misako Takayasu at. el., "Econophysics Approaches to Large-Scale Bussiness Data and Financial Crisis", Springer Japan, June 2010, pp.79-98.
【非特許文献２】Kenta Yamada, Hideki Takayasu, Takatoshi Ito and Misako Takayasu, "Solvable Stochastic Dealer Model for Financial Markets", Physical Review E,79, Issue 5, 051120 (2009).
【非特許文献３】中村和幸、樋口知之、「最近のベイズ理論の進展と応用［ＩＩ］  −逐次ベイズとデータ同化−」、２００９年、電子情報通信学会誌、Ｖｏｌ．９２、Ｎｏ．１２、ｐ１０６２〜１０６７
【発明の概要】
【発明が解決しようとする課題】
【００１２】
ところが、発明者らは、上述のＰＵＣＫモデルを用いる手法では、以下のような問題点が有ることを見出した。今日の金融市場などでは、コンピューターシステムを用いた高速の自動取引などが導入されており、価格を示す時系列データの変動を迅速かつ正確に検出することが求められる。しかしながら、通常のＰＵＣＫモデルでは、時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）を見積もるには、１０００ティック程度のデータを用い、式（２）の｛Ｐ（ｔ＋１）−Ｐ（ｔ）｝と｛Ｐ（ｔ）−ＰＭ（ｔ）｝の値をプロットし、最小二乗法などでその散布図の傾きを求めて時刻（ｔ）でのポテンシャル係数ｂ（ｔ）（以下、特に断らない限り、各次のポテンシャル係数ｂi（ｔ）（２≦ｉ≦Ｋ）の集合を、単にポテンシャル係数ｂｉ（ｔ）と称する）を計測する必要があった。そのため、データのサンプリング開始から１０００個程度の時系列データが収集されるまでは、過去情報が入れ替わらず、急峻なトレンド変化への追随性については、実用上、未だ不足していた。【００１３】
また、モデル上、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）の値は、直近のデータ数Ｍに依存しないはずであるが、観測から式（２）を満たす適切な直近のデータ数Ｍの範囲が時間ごとに変わる効果が無視できず、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）の推定誤差が大きくなるという問題点があった。そこで、ＰＵＣＫモデルに合う最適な直近のデータ数Ｍを逐次探すことが重要な課題であった。【００１４】
さらに、時刻（ｔ）における最適な真の価格Ｐ（ｔ）を推定するために、数ティック程度の最適移動平均の重みを、１日程度ごとの全データを用いて別途計算を更新しておく必要があった。【００１５】
本発明は上記の事情に鑑みて成されたものであり、本発明の目的は、時系列データの変化点検出方法及びプログラム、未来の時系列データ値の確率密度分布予測方法及びプログラムを提供することである。
【課題を解決するための手段】
【００１６】
本発明の第１の態様である時系列データの変化点検出方法は、時刻（ｔ−１）での真の市場価格及び時刻（ｔ−１）までのＭ（Ｍは、正の整数）個の真の市場価格の移動平均である中心価格によって定義されるポテンシャル項と変動誤差との和によって、時刻（ｔ）での真の市場価格Ｐ（ｔ）を算出するＰＵＣＫモデルに、粒子フィルタ法を適用する時系列データの変化点検出方法であって、それぞれ異なる値を有するＰＵＣＫモデルの状態を表すパラメータを持つ粒子群より次の時刻のパラメータの確率密度関数を得る第１のステップと、前記複数の粒子のそれぞれについて、時刻（ｔ）で観測した市場価格に対する前記時刻（ｔ）での真の市場価格の適合度を評価する第２のステップと、
前記適合度に応じて、前記複数の粒子の中から粒子を再サンプリングする第３のステップと、を備え、前記第３のステップでは、乱数を生成し、第１の所定の値と前記乱数とを比較し、前記乱数が前記第１の所定の値よりも大きい場合には、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする確率密度分布を成分とする確率密度関数を粒子として生成し、前記乱数が前記第１の所定の値よりも小さい場合には、一様分布を成分とする確率密度関数を粒子として生成するものである。これにより、時系列データである市場価格に急峻なトレンド変化が生じても、良好な追随性を与えるポテンシャル係数を算出し、ポテンシャル係数の変動を検出することにより変化点検出ができる。【００１７】
本発明の第２の態様である時系列データの変化点検出方法は、上述の時系列データの変化点検出方法において、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする確率密度分布は、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする正規分布であるものである。これにより、正規分布を用いて、具体的に確率密度分布を生成することができる。【００１８】
本発明の第３の態様である時系列データの変化点検出方法は、上述の時系列データの変化点検出方法において、前記時刻（ｔ＋１）での真の市場価格は、前記ＰＵＣＫモデルに基づいたポテンシャル項と、時刻（ｔ）での真の市場価格の変動誤差と、を加算することにより算出されることを特徴とするものである。これにより、粒子フィルタ法を用いた再サンプリングを好適に行うことができる。【００１９】
本発明の第４の態様である時系列データの変化点検出方法は、上述の時系列データの変化点検出方法において、上述の時系列データの変化点検出方法において、前記第３のステップでは、それぞれの粒子に関して乱数を生成し、前記第１の所定の値と前記乱数とを比較し、前記乱数が前記第１の所定の値よりも大きい場合には、時刻（ｔ）におけるＭを前提とする条件付き確率密度関数を粒子として生成し、前記乱数が前記第１の所定の値よりも小さい場合には、一様分布を成分とする確率密度関数を粒子として生成するものである。これにより、ポテンシャル係数のみならず、時刻（ｔ）から見て直近の時系列データ数を示すＭについても同時に時間発展を適用することができる。【００２０】
本発明の第５の態様である時系列データの変化点検出方法は、上述の時系列データの変化点検出方法において、前記時刻（ｔ＋１）での真の市場価格は、前記ポテンシャル項と、時刻（ｔ）での真の市場価格の変動誤差と、を加算することにより算出されることを特徴とするものである。これにより、時刻（ｔ）での真の市場価格を精度よく見積もることが可能となる。【００２１】
本発明の第６の態様である時系列データの変化点検出方法は、上述の時系列データの変化点検出方法において、前記時刻（ｔ）で観測した市場価格は、前記時刻（ｔ）での真の市場価格と時刻（ｔ）での市場価格の観測誤差との和であることを特徴とするものである。これにより、前記時刻（ｔ）で観測した市場価格に対する時刻（ｔ）での真の市場価格の適合度を、精度よく見積もることができる。【００２２】
本発明の第７の態様である時系列データの変化点検出方法は、上述の時系列データの変化点検出方法において、前記時刻（ｔ）での真の市場価格の変動誤差及び前記時刻（ｔ）での市場価格の観測誤差は、正規分布等に従う確率密度関数で与えられることを特徴とするものである。これにより、時刻（ｔ）での真の市場価格及び前記時刻（ｔ）で観測した市場価格に対する時刻（ｔ）での真の市場価格の適合度を、精度よく見積もることができる。【００２３】
本発明の第８の態様である未来の時系列データ値の確率密度分布の予測方法は、上述の時系列データの変化点検出方法により算出した時刻（ｔ）における前記ポテンシャル項のポテンシャル係数、前記Ｍの値、及び時刻（ｔ）での真の市場価格の変動誤差と、に基づいて、前記時刻（ｔ）での真の市場価格Ｐ（ｔ）を時間発展させることにより、時刻（ｔ＋Ｎ）での前記時刻（ｔ）での真の市場価格Ｐ（ｔ）を算出するものである。これにより、ＰＵＣＫモデルを用いて、未来の時系列データ値の確率密度分布を予測することができる。【００２４】
本発明の第９の態様である未来の時系列データ値の確率密度分布の予測方法は、上述の未来の時系列データ値の確率密度分布の予測方法において、前記粒子の総数をＮｐ、前記粒子を特定する番号をｊ（ただし、ｊは１≦ｊ≦Ｎｐを満たす整数）、番号ｊの時刻（ｔ）での粒子のポテンシャル係数をｂ（ｊ）ｉ（ｔ）（ただし、ｉはポテンシャルの次数を示す２以上の整数）、前記番号ｊの粒子の時刻（ｔ）での真の市場価格の変動誤差をｆ（ｊ）Ｐ（ｔ）として、前記番号ｊの粒子の時刻（ｔ）での真の市場価格をＰ（ｊ）（ｔ）、前記番号ｊの粒子の時刻（ｔ＋Ｎ）での真の市場価格をＰ（ｊ）（ｔ＋Ｎ）、前記番号ｊの粒子の時刻（ｔ）での市場価格の中心価格をＰ（ｊ）Ｍ（ｔ）、時刻（ｔ）におけるＭの値をＭ（ｔ）とした場合に、前記番号ｊの粒子の前記時刻（ｔ＋Ｎ）での真の市場価格Ｐ（ｊ）（ｔ＋Ｎ）は、Ｋを２以上の整数として、
【数６】
0006004284.tif 000007
で表されるものである。これにより、ＰＵＣＫモデルを用いて、未来の時系列データ値の確率密度分布を具体的に予測することができる。【００２５】
本発明の第１０の態様である未来の時系列データ値の確率密度分布の予測方法は、上述の未来の時系列データ値の確率密度分布の予測方法において、番号ｊの粒子の時刻（ｔ＋１）でのポテンシャル係数ｂ（ｊ）ｉ（ｔ＋１）は、前記番号ｊの粒子の時刻（ｔ）でのポテンシャル係数ｂ（ｊ）ｉ（ｔ）を時間発展させることにより算出されるものである。これにより、ＰＵＣＫモデルを用いて、未来の時系列データ値の確率密度分布を具体的に予測することができる。【００２６】
本発明の第１１の態様である時系列データの変化点検出プログラムは時刻（ｔ−１）での真の市場価格及び時刻（ｔ−１）までのＭ（Ｍは、正の整数）個の真の市場価格の移動平均である中心価格によって定義されるポテンシャル項と変動誤差との和によって時刻（ｔ）での真の市場価格Ｐ（ｔ）を算出するＰＵＣＫモデルに粒子フィルタ法を適用する時系列データの変化点検出をコンピュータに実行させる時系列データの変化点検出プログラムであって、それぞれ異なる値を有するＰＵＣＫモデルの状態を表すパラメータを持つ粒子群を生成することによりパラメータの確率密度関数を得る第１の処理と、適合度評価部に、前記複数の粒子のそれぞれについて、時刻（ｔ）で観測した市場価格に対する前記時刻（ｔ）での真の市場価格の適合度を評価する第２の処理と、サンプリング部に、前記適合度に応じて、前記複数の粒子の中から粒子を再サンプリングさせる第３の処理と、をコンピュータに実行させ、前記乱数が前記第１の所定の値よりも大きい場合には、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする確率密度分布を成分とする確率密度関数を粒子として再生成する処理と、前記乱数が前記第１の所定の値よりも小さい場合には、一様分布にしたがって粒子を生成する処理と、を実行させるものである。これにより、時系列データである市場価格に急峻なトレンド変化が生じても、良好な追随性を与えるポテンシャル係数を算出し、ポテンシャル係数の変動を検出することにより変化点検出ができる。【００２７】
本発明の第１２の態様である時系列データの変化点検出プログラムは、上述の変化点検出プログラムにおいて、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする確率密度分布は、時刻（ｔ）でのＰＵＣＫモデルの状態を表すパラメータを平均値とする正規分布であるものである。これにより、正規分布を用いて、具体的に確率密度関数を生成することができる。【００２８】
本発明の第１３の態様である時系列データの変化点検出プログラムは、上述の変化点検出プログラムにおいて、前記第３のステップでは、前記適合度が大きい粒子ほど、多くサンプリングされるだけでなく、一定の割合で、実現可能性のある適合度が低い粒子を生成すること、を特徴とするものである。これにより、粒子フィルタ法を用いた再サンプリングを好適に行うことができる。【００２９】
本発明の第１４の態様である時系列データの変化点検出プログラムは、上述の変化点検出プログラムにおいて、前記第１のステップでは、乱数を生成する処理と、時刻（ｔ）におけるＰＵＣＫモデルのパラメータを前提とする条件付き確率密度関数を粒子として生成する処理と、前記乱数が前記第３の所定の値よりも小さい場合には、一様分布を成分とする確率密度関数を粒子として生成する処理と、を実行させることを特徴とするものである。これにより、ポテンシャル係数のみならず、時刻（ｔ）から見て直近の時系列データ数を示すＭについても同時に時間発展を適用することができる。【００３０】
本発明の第１５の態様である時系列データの変化点検出プログラムは、上述の変化点検出プログラムにおいて、前記時刻（ｔ）での真の市場価格は、前記ポテンシャル項と、時刻（ｔ）での真の市場価格の変動誤差と、を加算することにより算出されることを特徴とするものである。これにより、前記時刻（ｔ）で観測した市場価格に対する時刻（ｔ）での真の市場価格の適合度を、精度よく見積もることができる。【００３１】
本発明の第１６の態様である時系列データの変化点検出プログラムは、上述の変化点検出プログラムにおいて、前記時刻（ｔ）で観測した市場価格は、前記時刻（ｔ）での真の市場価格と時刻（ｔ）での市場価格の観測誤差との和であることを特徴とするものである。これにより、前記時刻（ｔ）で観測した市場価格に対する時刻（ｔ）での真の市場価格の適合度を、精度よく見積もることができる。【００３２】
本発明の第１７の態様である時系列データの変化点検出プログラムは、上述の変化点検出プログラムにおいて、前記時刻（ｔ）での真の市場価格の変動誤差及び前記時刻（ｔ）での市場価格の観測誤差は、正規分布等に従う確率密度関数で与えられることを特徴とするものである。これにより、時刻（ｔ）での真の市場価格及び前記時刻（ｔ）で観測した市場価格に対する時刻（ｔ）での真の市場価格の適合度を、精度よく見積もることができる。【００３３】
本発明の第１８の態様である未来の時系列データ値の確率密度分布の予測プログラムは、上述の時系列データの変化点検出プログラムにより算出した時刻（ｔ）における前記ポテンシャル項のポテンシャル係数、前記Ｍの値、及び時刻（ｔ）での真の市場価格の変動誤差と、に基づいて、前記時刻（ｔ）での真の市場価格Ｐ（ｔ）を時間発展させることにより、時刻（ｔ＋Ｎ）での前記時刻（ｔ）での真の市場価格Ｐ（ｔ）を算出するものである。これにより、ＰＵＣＫモデルを用いて、未来の時系列データ値の確率密度分布を予測することができる。【００３４】
本発明の第１９の態様である未来の時系列データ値の確率密度分布の予測プログラムは、上述の未来の時系列データ値の確率密度分布の予測プログラムにおいて、前記粒子の総数をＮｐ、前記粒子を特定する番号をｊ（ただし、ｊは１≦ｊ≦Ｎｐを満たす整数）、番号ｊの粒子の時刻（ｔ）でのポテンシャル係数をｂ（ｊ）ｉ（ｔ）（ただし、ｉはポテンシャルの次数を示す２以上の整数）、前記番号ｊの粒子の時刻（ｔ）での真の市場価格の変動誤差をｆ（ｊ）Ｐ（ｔ）として、前記番号ｊの粒子の時刻（ｔ）での真の市場価格をＰ（ｊ）（ｔ）、前記番号ｊの粒子の時刻（ｔ＋Ｎ）での真の市場価格をＰ（ｊ）（ｔ＋Ｎ）、前記番号ｊの粒子の時刻（ｔ）での市場価格の中心価格をＰ（ｊ）Ｍ（ｔ）、時刻（ｔ）におけるＭの値をＭ（ｔ）とした場合に、前記番号ｊの粒子の前記時刻（ｔ＋Ｎ）での真の市場価格Ｐ（ｊ）（ｔ＋Ｎ）は、Ｋを２以上の整数として、
【数７】
0006004284.tif 000008
で表されるものである。これにより、ＰＵＣＫモデルを用いて、具体的に未来の時系列データ値の確率密度分布を予測することができる。【００３５】
本発明の第２０の態様である未来の時系列データ値の確率密度分布の予測プログラムは、上述の未来の時系列データ値の確率密度分布の予測プログラムにおいて、番号ｊの粒子の時刻（ｔ＋１）でのポテンシャル係数ｂ（ｊ）ｉ（ｔ＋１）は、前記番号ｊの粒子の時刻（ｔ）でのポテンシャル係数ｂ（ｊ）ｉ（ｔ）を時間発展させることにより算出されるものである。
【発明の効果】
【００３６】
本発明によれば、時系列データの急峻な変化を精度よく検出できる時系列データの変化点検出方法及びプログラム、未来の時系列データ値の確率密度分布予測方法及びプログラムすることができる。
【図面の簡単な説明】
【００３７】
【図１】実施の形態１で用いる粒子フィルタ法の手順を示すフローチャートである。
【図２】実施の形態１にかかる式（２４）で与えられる確率分布を示すグラフである。
【図３】実施の形態１にかかる式（２５）で与えられる確率分布を示すグラフである。
【図４】実施の形態１における時刻（ｔ）でのポテンシャル係数ｂ（ｔ）の決定方法を示すフローチャートである。
【図５】実施の形態１における時刻（ｔ）での直近のデータ数Ｍ（ｔ）の決定方法を示すフローチャートである。
【図６】実施の形態２にかかる時系列データの変化点検出プログラムが実行される時系列データの変化点検出装置２００の構成を模式的に示すブロック図である。
【図７】実施の形態３にかかる粒子フィルタ法を適用したＰＵＣＫモデルの時間発展方法を示すフローチャートである。
【図８】実施の形態３にかかる未来の時系列データがとる値の確率密度分布の予測方法により時刻（ｔ＋Ｎ）における市場価格の確率密度分布を求めた場合の確率密度分布の推移を示すグラフである。
【図９】実施の形態３にかかる未来の時系列データがとる値の確率密度分布の予測方法により時刻（ｔ＋Ｎ）における市場価格がとる値の確率密度分布を示すグラフである。
【図１０】実施の形態４にかかる未来の時系列データがとる値の確率密度分布の予測プログラムが実行される未来の時系列データがとる値の確率密度分布の予測装置４００の構成を模式的に示すブロック図である。
【図１１】実施の形態５にかかる粒子フィルタ法を適用したＰＵＣＫモデルの時間発展方法を示すフローチャートである。
【図１２】現在の市場価格を数理モデルに代入して、粒子番号ｊ毎に時間発展させ、その経過を粒子毎に表示したグラフである。
【図１３】実施の形態６にかかる未来の時系列データがとる値の確率密度分布の予測プログラムが実行される未来の時系列データがとる値の確率密度分布の予測装置６００の構成を模式的に示すブロック図である。
【発明を実施するための形態】
【００３８】
以下、図面を参照して本発明の実施の形態について説明する。各図面においては、同一要素には同一の符号が付されており、必要に応じて重複説明は省略される。
【００３９】
実施の形態１
本発明の実施の形態１にかかる時系列データの変化点検出方法について説明する。本実施の形態における時系列データの変化点検出方法は、前述のＰＵＣＫモデルに、いわゆる粒子フィルタ法（非特許文献３）を適用したものである。粒子フィルタ法を適用することにより、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）及び直近のデータ数Ｍ（ｔ）の値を、逐次予測することができる。
【００４０】
まず、粒子フィルタ法について、概要を説明する。図１は、実施の形態１で用いる粒子フィルタ法の手順を示すフローチャートである。逐次モンテカルロ法の一種である粒子フィルタ法では、状態ベクトルｘｔが、以下の式（６）で与えられる。つまり、時刻（ｔ）における状態ベクトルｘｔは、１ティック前の状態ベクトルｘｔ−１にパラメータθを与えることにより得られる。
【数８】
0006004284.tif 000009
ここで、〜は、ｘｔ−１とθをパラメータとし、時刻（ｔ）に依存する確率密度関数、Ｑｔ（・）に従う確率変数ｘｔを生成する操作である。
【００４１】
また、観測ベクトルｙｔは、以下の式（７）で与えられる。ここで、Ｒｔ（・）はｘｔ−１とθをパラメータとし、時刻（ｔ）に依存する確率密度関数である。
【数９】
0006004284.tif 000010
【００４２】
粒子フィルタ法においては、状態ベクトルｘｔの確率密度関数は、観測ベクトル｛ｙ１,ｙ２,・・・,ｙｔ｝により決定される条件付き確率密度関数ＣＰＤＦとして、以下の式（８）で表される。
【数１０】
0006004284.tif 000011
【００４３】
ここで、粒子フィルタ法で想定する粒子（状態ベクトル）がＮ（Ｎは正の整数）個であるとすると、Ｎ個の状態ベクトルの条件付き確率密度関数（ＰＤＦ）は、以下の式（９）で表される。
【数１１】
0006004284.tif 000012
ここで、ｘは、ｊ（１≦ｊ≦Ｎ）番目の粒子において、時刻ｋの観測ベクトルから推定された時刻ｌの状態ベクトルである。
【００４４】
粒子フィルタ法では、Ｎ個の粒子の条件付き確率密度関数（ＰＤＦ）について、まず、ｔ＝０における初期値を与える（図１のステップＳ１）。なお、初期値は、以下の式（１０）で表される。
【数１２】
0006004284.tif 000013
そして、粒子フィルタ法による計算を開始するにあたり、時刻（ｔ）をｔ＝１に設定する。
【００４５】
次いで、Ｎ個の粒子について、式（６）を用いて、時刻（ｔ−１）における確率密度関数を基に時刻（ｔ）における予測粒子の確率密度関数（ＰＤＦ）を生成する（図１のステップＳ２）。なお、ｊ番目の粒子の状態ベクトルは、以下の式（１１）で表される。
【数１３】
0006004284.tif 000014
ここで、ｑｔ（・）は、式（９）によって定義される確率密度関数である。
【００４６】
次いで、各粒子の重み係数Ｗｊ（以下の式（１２））を、尤度Ｌｊ（以下の式（１３））を用いて計算する。なお、尤度Ｌｊは、式（７）より得られる。換言すると、ここでは、実際の時系列データから得られる観測ベクトルｙｔと、計算により得られた状態ベクトルｘｔと、の適合度を評価する（図１のステップＳ３）。この適合度は、尤度を用いることで評価することができる。
【数１４】
0006004284.tif 000015
【数１５】
0006004284.tif 000016
【００４７】
その後、各粒子の重み係数に従って、粒子の再サンプリングを行う（図１のステップＳ４）。つまり、式（１１）で示される各粒子を、重み係数に応じた個数だけサンプリングし、合計Ｎ個の粒子をサンプリングする。この際、適合度が小さければ重み係数は小さくなるので、その粒子は残らない可能性が高くなる。一方、適合度が大きければ重み係数は大きくなり、その粒子は複製され増加することとなる。以上が粒子フィルタ法の概要である。
【００４８】
次いで、上記の手順を、時刻（ｔ）を１ステップずつ進めて、最後のステップ（ｔ＝ｔｅｎｄ）に至るまで繰り返し行う（図１のステップＳ５及びＳ６）。これにより、時間の進展に応じて、逐次粒子フィルタ法を適用することが可能である。
【００４９】
次に、本実施の形態において、ＰＵＣＫモデルに粒子フィルタ法を用いる手順について説明する。本実施の形態では、図１のステップＳ２における、確率密度関数の生成時のＰＵＣＫモデルの取り扱いに特徴を有する。以下、ステップＳ２における具体例について説明する。以下では、状態ベクトルを用いて、時刻（ｔ）での各次のポテンシャル係数をｂｉ（ｔ）、直近のデータ数をＭ（ｔ）等で表す。時刻（ｔ−１）での各次のポテンシャル係数ｂｉ（ｔ）は以下の式（１４）、直近のデータ数Ｍ（ｔ−１）は以下の式（１５）で表される。但し、θｂ及びθＭはノイズを生成するためのパラメータである。
【数１６】
0006004284.tif 000017
【数１７】
0006004284.tif 000018
【００５０】
式（１５）を前述の式（３）に適用することにより、以下の式（１６）を得ることができる。
【数１８】
0006004284.tif 000019
【００５１】
また、式（１４）及び（１５）を前述の式（２）に適用することにより、以下の式（１７）を得ることができる。
【数１９】
0006004284.tif 000020
【００５２】
この場合においても、前述の式（１）と同様に、式（１８）が成立する。
【数２０】
0006004284.tif 000021
【００５３】
なお、時刻（ｔ）での真の市場価格の観測誤差ｆＰ（ｔ）及び時刻（ｔ）での市場価格の観測誤差ｆＹ（ｔ）は、以下の式（１９）のように、正規分布等に従うものとする。
【数２１】
0006004284.tif 000022
以下、一般に、Ｎ（・）の部分は、正規分布のみならず、さらに高次のモーメント等で特徴付けられる確率分布でもよい。また、ｆＰ（ｔ）、ｆｙ（ｔ）は、異なる確率分布でもよい。
【００５４】
時刻（ｔ）での真の市場価格の観測誤差ｆＰ（ｔ）の標準偏差σＰは、状態ベクトルを用いて、以下の式（２０）で表すことができる。
【数２２】
0006004284.tif 000023
【００５５】
時刻（ｔ）での市場価格の観測誤差ｆＹ（ｔ）の標準偏差σＹは、状態ベクトルを用いて、以下の式（２１）で表すことができる。
【数２３】
0006004284.tif 000024
【００５６】
式（２０）及び（２１）より、式（１９）を以下の式（２２）に書き換えることができる。
【数２４】
0006004284.tif 000025
【００５７】
さて、上述のモデルの時間発展について考える。以下、変数Ａ（ｔ）は、平均値Ａ（ｔ−１）、分散σＡ２の正規分布等で与えられるものとする。つまり、変数Ａ（ｔ）は、以下の式（２３）で表される。
【数２５】
0006004284.tif 000026
【００５８】
変数Ａ（ｔ）は、時刻（ｔ）よりも前の観測値｛ｙ（１）,・・・,ｙ（ｔ−１）｝により決定される。このときのＡ（ｔ）をＡｔ｜ｔ−１と表記すると、式（２３）は式（２４）に書き換えることができる。
【数２６】
0006004284.tif 000027
図２は、式（２４）で与えられる確率分布を示すグラフである。図２に示すように、Ａｔ｜ｔ−１は、ＡｍｉｎからＡｍａｘの範囲で、平均値Ａｔ−１｜ｔ−１、標準偏差σＡの正規分布等となる。
【００５９】
しかし、式（２４）及び図２に示すＡｔ｜ｔ−１は、２つの問題点が存在する。１つは、例えば実際の金融市場などでは、市場価格の急峻な変動が存在するが、ランダムウォークに基づく上述のモデルでは、こうした急峻な変化に追随することができない。また、現実系では、パラメータがとりうる値の範囲が有限である。よって、本実施の形態では、変数Ａ（ｔ）を、切断正規分布等と一様値分布との混合分布として、以下の式（２５）で表すこととする。
【数２７】
0006004284.tif 000028
式（２５）の｛（１−Ａｍｕｔ）ｔｒｕｎｃＮ（Ａ（ｔ−１）,σＡ２；Ａｍｉｎ,Ａｍａｘ）｝は、平均値Ａ（ｔ−１）、分散σＡ２、区間（Ａｍｉｎ,Ａｍａｘ）の規格化された切断された正規分布等である。ｔｒｕｎｃＮは、平均や標準偏差等で特徴づけられる一般の分布を上限と下限で切断したものに拡張することができる。また、Ａｍｕｔは、市場価格の急峻な変化などの突発現象に対応するために、式（２５）で新たに導入した係数である。つまり、Ａ（ｔ）は、正規分布等と一様分布Ｕ（Ａｍｉｎ，Ａｍａｘ）の２つの成分が混合した形式で表される。Ａｍｕｔの値により、一様分布Ｕ（Ａｍｉｎ，Ａｍａｘ）の混合率が決定される。
【００６０】
図３は、式（２５）で与えられる確率分布を示すグラフである。図３に示すように、Ａｔ｜ｔ−１は、ＡｍｉｎからＡｍａｘの範囲で、平均値Ａｔ−１｜ｔ−１、標準偏差σＡの正規分布等と一定値Ａｍｕｔとが混合された分布として表現される。
【００６１】
続いて、ＰＵＣＫモデルの各パラメータに、式（２５）を適用する具体例を説明する。ＰＵＣＫモデルの各パラメータに式（２５）を適用する処理については、上述のステップＳ１において行われる。図４は、実施の形態１における時刻（ｔ）でのポテンシャル係数ｂ（ｔ）の決定方法を示すフローチャートである。まず、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）について説明する。このとき、式（２５）を式（２６）に書き換えることができる。
【数２８】
0006004284.tif 000029
【００６２】
式（２６）を用いて、粒子フィルタ法を適用することにより、時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）を見積もることが可能である。ここで、時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）における正規分布等と一様分布との混合比を決定する必要がある。そのために、ここでは、乱数ｕを生成する（図４のステップＳ１１）。乱数ｕは、０≦ｕ≦１を満たす任意の実数である。式（２５）におけるｂｍｕｔは、時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）の時間発展において、突発現象が起きる割合を示す係数である。本実施の形態では、ｂｉ，ｍｕｔと乱数ｕを比較する（図４のステップＳ１２）ことにより、正規分布等と一様分布との混合比を決定する。ｂｉ，ｍｕｔは、無作為に生成される乱数ｕと比較されるので、結局のところ、混合比は（１−ｂｉ，ｍｕｔ）：ｂｉ，ｍｕｔに収束することとなる。
【００６３】
すなわち、ｕ≧ｂｉ，ｍｕｔである場合には、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）は、正規分布等のみから成り立つものとする（図４のステップＳ１３）。よって、この場合の時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）は、以下の式（２７）で表される。
【数２９】
0006004284.tif 000030
【００６４】
一方、ｕ＜ｂｉ，ｍｕｔである場合には、時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）は、一様分布のみから成り立つものとする（図４のステップＳ１４）。よって、この場合の時刻（ｔ）での各次のポテンシャル係数ｂｉ（ｔ）は、以下の式（２８）で表される。
【数３０】
0006004284.tif 000031
【００６５】
本実施の形態においては、上述のように、時系列データの急峻に変化する突発現象を考慮するために一様分布を導入している。換言すれば、一様分布が無い場合には、粒子の世代が進むごとに、適合度の低い粒子については、淘汰されて消滅してしまう。このような状態で突発現象が生じた場合には、突発状態に適合する粒子が既に消滅しており、突発現象に追随できない。ところが、本実施の形態におけるように、突発現象を考慮した一様分布を導入することにより、適合度の低い粒子の消滅を防ぐことができる。これにより、突発現象に対して適合度の高い粒子を残存させることが可能となり、通常の粒子フィルタ法と比べて、突発現象に対する追随能力を高めることが可能となる。このように、一様分布を導入することにより、通常の正規分布やローレンツ分布を用いることでは実現できない、特有の効果を実現することが可能となる。
【００６６】
次いで、直近のデータ数Ｍ（ｔ）について説明する。直近のデータ数Ｍ（ｔ）についても、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）と同様に、突発現象を考慮するための一様分布を導入する。本実施の形態では、直近のデータ数Ｍ（ｔ）を、時刻（ｔ−１）における確率密関数と一様分布との混合分布として表す。この直近のデータ数Ｍ（ｔ）についての処理は、上述のステップＳ１において行われる。図５は、実施の形態１における時刻（ｔ）での直近のデータ数Ｍ（ｔ）の決定方法を示すフローチャートである。この場合、直近のデータ数Ｍ（ｔ）は、式（２９）で表される。
【数３１】
0006004284.tif 000032
【００６７】
ここで、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）の場合と同様に、乱数ｕを生成し（図５のステップＳ２１）、Ｍｍｕｔと乱数ｕを比較する（図５のステップＳ２２）ことにより、時刻（ｔ−１）における確率密度関数と一様分布との混合比を決定する。
【００６８】
すなわち、ｕ≧Ｍｍｕｔである場合には、直近のデータ数Ｍ（ｔ）は、時刻（ｔ−１）における確率密度関数のみから成り立つものとする（図５のステップＳ２３）。よって、この場合の直近のデータ数Ｍ（ｔ）は、以下の式（３０）で表される。
【数３２】
0006004284.tif 000033
ただし、ここでの、Ｕ（Ｍｍｉｎ，Ｍｍａｘ）は、下限をＭｍｉｎ、上限をＭｍａｘとする整数値をとるような一様な確率密度関数である。
【００６９】
一方、ｕ＜Ｍｍｕｔである場合には、直近のデータ数Ｍ（ｔ）は、一様分布のみから成り立つものとする（図５のステップＳ２４）。よって、この場合の直近のデータ数Ｍ（ｔ）は、式（３１）で表される。
【００７０】
以下に想定される確率密度関数Ｇ（・）の具体例を示す。Ｍｍｉｎ＜Ｍ（ｔ−１）＜Ｍｍａｘの場合には、確率密度関数Ｇ（・）は、以下の式（３１）で表される。
【数３３】
0006004284.tif 000034
ただし、０＜ｒ＜１とする。
また、境界条件においては、確率密度関数Ｇ（・）は、以下の式（３２）で表される。
【数３４】
0006004284.tif 000035
上記の範囲以外でのＭ（ｔ−１）に関しては、確率密度関数Ｇ（・）は、以下の式（３３）で表される。
【数３５】
0006004284.tif 000036
【００７１】
式（２２）の変動誤差ｆＰ（ｔ）の標準偏差σＰ（ｔ）、（観測誤差ｆＹ（ｔ）の標準偏差σＹ（ｔ）も以下同様の手順）、は、次式で与えられる。
【数３６】
0006004284.tif 000037
ここで、σｐｍｕｔは、０と１の間の数値であり、σｐｍinとσｐｍax はそれぞれσｐ（ｔ）のとりうる範囲の下限と上限を表す。また、γｐ２は分散を表す数値である。
【数３７】
0006004284.tif 000038
ここで、σYｍｕｔは、０と１の間の数値であり、σYｍinとσYｍax はそれぞれσY（ｔ）のとりうる範囲の下限と上限を表す。また、γY２は分散を表す数値である。
【００７２】
以上、説明したように、本実施の形態にかかる時系列データの変化点検出方法によれば、時系列データが安定的なトレンドにある場合には、安定的なトレンドへの適合度が高いパラメータ（ポテンシャル係数ｂ、直近のデータ数Ｍ、標準偏差σＰ及びσＹ）を示す粒子を増加させ、トレンドに好適に追随することが可能である。
【００７３】
さらに、本実施の形態にかかる時系列データの変化点検出方法によれば、安定的なトレンドから逸脱した、急峻な変動、すなわち突発現象が生じた場合でも、トレンドへの追随が可能である。すなわち、通常の粒子フィルタ法であれば、安定的なトレンドが継続すると、安定的なトレンドに適合する粒子のみが増加し、安定的なトレンドに適合しないものの突発現象には適合する粒子は淘汰され、消滅してしまう。
【００７４】
ところが、本実施の形態にかかる時系列データの変化点検出方法では、最新の時刻におけるパラメータの確率密度関数を算出するにあたり、各パラメータに対して、式（２８）及び式（３０）に示すように、一定の割合で一様分布を付与する。これにより、安定的なトレンドに対する適合度が悪いために消滅する粒子にも、図１のステップＳ２で確率密度関数を生成する場合に、一定の割合で一様分布が付与される。これにより、安定的なトレンドに適合しないものの突発現象には適合する粒子は、粒子フィルタ法による逐次計算が進んでも残存することができる。そして、ある時点で突発現象が発生した場合に、残存している粒子の内で突発現象への適合度の高い粒子を増加させることで、好適に突発現象のような動力学的現象に追随することが可能である。
【００７５】
そして、ポテンシャル係数ｂｉ（ｔ）及び直近のデータ数Ｍ（ｔ）を逐次算出することにより、時系列データの変動の状態を検出することが可能である。２次のポテンシャル係数ｂ２（ｔ）が正の値をとり、その他のｂｉ（ｔ）が０近傍であれば時系列の変動は安定的である。しかし、逐次解析の結果、２次のポテンシャル係数ｂ２（ｔ）が負の値をとる粒子の適合度が増加する場合には、時系列の変動は不安定な状態となる。つまり、２次のポテンシャル係数ｂ２（ｔ）の値の変動を監視することにより、時系列データの変動が安定的であるか、不安定であるか、の境界点を検出することが可能となる。また、３次のポテンシャル係数ｂ３（ｔ）が正の値をとるときには時系列に急激な下落トレンドが、ｂ３（ｔ）が負の値をとるときには時系列に急激な上昇トレンドが、発生する可能性が予想される。
【００７６】
加えて、粒子フィルタ法では、ポテンシャル係数ｂｉ（ｔ）及び直近のデータ数Ｍ（ｔ）について、複数の粒子として表現するが、時刻の進展に応じて、適合度が高い粒子は変化するので、適合度が高いポテンシャル係数ｂｉ（ｔ）及び直近のデータ数Ｍ（ｔ）の値も変化する。適合度の高い直近のデータ数は可変であるので、仮に時系列データの蓄積数が少ない場合でも、利用可能なデータを用いて、適合度の高い直近のデータ数Ｍを算出することが可能である。つまり、通常のＰＵＣＫモデルでは、ある程度のデータ数（１０００点程度）のデータ数が必要であるが、本方法では、より少ないデータ数でも利用可能である。さらに、既存手法では不可欠だった最適移動平均等のデータのノイズ除去処理の計算が不必要となった。従って、ＰＵＣＫモデルのパラメータ推定の既存手法と比べて、データの蓄積が少ない段階でもＰＵＣＫモデルを適用することができるという利点を有する。
【００７７】
さらにまた、ポテンシャル係数ｂｉ（ｔ）の値と適合度、すなわち粒子数を参照することにより、時系列データの変動環境を解析することが可能である。具体的には、円／ドル取引を例とすると、２次のポテンシャル係数ｂ２（ｔ）が正の値をとる粒子についての適合度が高い傾向が有る場合には、価格の変化が安定であり、ランダムウォークと比較し、短い時間スケールでは逆方向に引き戻すように価格が変動しやすい。この場合には、市場では、価格の上昇あるいは下降のトレンドが反転しやすく、結果として、大きな価格変動は生じにくい。他方、２次のポテンシャル係数ｂ２（ｔ）が負の値をとる粒子についての適合度が高い傾向が有る場合には、価格の変化は不安定であり、ランダムウォークと比較し、短い時間スケールでは、同じ方向の価格変化が持続的に発生しやすく、上昇あるいは下降のトレンドは増幅されやすい。この場合には、トレンドを追随する戦略を持つトレーダーが市場において数的に支配的であると推測できる（非特許文献２）。このように、本方法を、為替取引のような金融市場における価格などの時系列データに適用することにより、市場におけるトレーダーの集団的な行動特性を定量的に評価することが可能となる。
【００７８】
上述の市場におけるトレーダーの集団的な行動特性を定量的に評価できるという効果は、本実施の形態にかかる時系列データの変化点検出方法特有の効果である。本実施の形態にかかる時系列データの変化点検出方法では、複数の値の異なる時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ）及び直近のデータ数Ｍについて、それぞれ適合度が異なる粒子として表現することができる。これにより、時刻（ｔ）において適合度が高く、時系列データの変動に支配的な影響を与えると考えられる、値の異なる時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ）及び直近のデータ数Ｍ（ｔ）を、複数個同時に評価することができる。これにより、異なる行動特性を有するトレーダーの、行動の態様（時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ））や行動を決定するため参照すると考えられるデータ数（直近のデータ数Ｍ（ｔ））を、同時に評価することができる。これに対し、前述のような粒子フィルタ法を適用しない通常のＰＵＣＫモデルでは、直近のデータ数Ｍを仮定し、仮定した直近のデータ数Ｍに対する時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ）を計算することしかできなかった。従って、原理的に、ポテンシャル係数ｂｉ（ｔ）及び中心を見積もるために必要な直近のデータ数Ｍ（ｔ）を同時に評価することはできなかった。つまり、時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ）及び直近のデータ数Ｍ（ｔ）等のＰＵＣＫモデルのパラメータの確率密度関数等を推定することはできなかった。つまり、ＰＵＣＫモデルによる複数個の時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ）及び直近のデータ数Ｍの同時評価は、本実施の形態にかかる時系列データの変化点検出方法を用いることで初めて実現することが可能となる。
【００７９】
なお、本実施の形態にかかる時系列データの変化点検出方法は、上述のステップＳ１〜Ｓ６を表現するアルゴリズムが記述されたプログラムとして提供することが可能である。このプログラムをコンピュータなどのハードウェアに実行させることにより、上述の時系列データの変化点検出などの、本実施の形態にかかる時系列データの変化点検出方法により得られるのと同様の効果をえることができる。例えば、当該プログラムを実行し、時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ）等を表示装置に逐次表示し、ポテンシャル係数の経時的な変化をリアルタイムで観察することが可能である。
【００８０】
実施の形態２
次に、実施の形態２にかかる時系列データの変化点検出プログラムについて説明する。図６は、実施の形態２にかかる時系列データの変化点検出プログラムが実行される時系列データの変化点検出装置２００の構成を模式的に示すブロック図である。時系列データの変化点検出装置２００は、実施の形態１にかかる時系列データの変化点検出方法のステップＳ１〜Ｓ６を表現するアルゴリズムが記述されたプログラムが実行される装置である。図６に示すように、時系列データの変化点検出装置２００は、記憶部１、演算部２、表示部３及びバス４を有する。
【００８１】
記憶部１は、ハードディスク、ＤＲＡＭ、ＳＲＡＭ、フラッシュメモリなどの記憶装置により構成され、過去の時系列データ、図１のステップＳ１で設定する初期値の情報、及び実施の形態１にかかる時系列データの変化点検出方法のステップＳ１〜Ｓ６を表現するアルゴリズムが記述されたプログラム２０が記憶されている。なお、記憶部に記憶されている過去の時系列データ、図１のステップＳ１で設定する初期値の情報及びプログラム２０は、外部から適宜書き換えることが可能である。また、時系列データについては、時間の進展に応じて、新しいデータ要素を逐次追加することが可能である。
【００８２】
演算部２は、バス４を介して、記憶部１から時系列データ、図１のステップＳ１で設定する初期値の情報、及びプログラム２０を読み込む。演算部２は、読み込んだプログラム２０を実行し、時系列データの変化点検出動作を行う。プログラム２０は、初期設定部２１、ＰＤＦ生成部２２、適合度評価部２３、再サンプリング部２４及びカウント部２５を有する。
【００８３】
具体的には、初期設定部２１では、図１のステップＳ１に対応する処理を行う。すなわち、Ｎ個の粒子の条件付き確率密度関数（ＰＤＦ）について、まず、式（１０）で示されるｔ＝０における初期値を与える。そして、粒子フィルタ法による計算を開始するにあたり、時刻（ｔ）をｔ＝１に設定する。
【００８４】
ＰＤＦ生成部２２では、図１のステップＳ２に対応する処理を行う。すなわち、Ｎ個の粒子について、式（６）を用いて、時刻（ｔ−１）における確率密度関数（ＰＤＦ）を基に、時刻（ｔ）における予測粒子の確率密度関数（ＰＤＦ）を生成する（式（１１））。
【００８５】
適合度評価部２３では、図１のステップＳ３に対応する処理を行う。すなわち、各粒子の重み係数Ｗｉ（式（１２））を、尤度Ｌｉ（式（１３））を用いて計算する。なお、尤度Ｌｉは、式（７）より得られる。つまり、実際の時系列データから得られる観測ベクトルｙｔと、計算により得られた状態ベクトルｘｔと、の適合度を評価する。
【００８６】
再サンプリング部２４では、図１のステップＳ４に対応する処理を行う。すなわち、各粒子の重み係数に従って、粒子の再サンプリングを行う。つまり、式（１１）で示される各粒子を、重み係数に応じた個数だけサンプリングし、合計Ｎ個の粒子をサンプリングする。この際、適合度が小さければ重み係数は小さくなるので、その粒子は残らない可能性が高くなる。一方、適合度が大きければ重み係数は大きくなり、その粒子は複製され増加することとなる。
【００８７】
演算部２は、再サンプリング部２４によって得られた重み係数（すなわち粒子数）及び各粒子で表される時刻（ｔ）でのポテンシャル係数ｂｉ（ｔ）や直近のデータ数Ｍ（ｔ）の情報を、バス４を介して、表示部３に出力する。
【００８８】
カウント部２５では、図１のステップＳ５及びＳ６に対応する処理を行う。すなわち、  最新の処理の時刻が最後のステップ（ｔ＝ｔｅｎｄ）であるかを判定する（図１のステップＳ５での処理）。そして、最新の処理の時刻が最後のステップ（ｔ＝ｔｅｎｄ）でなければ、時刻（ｔ）を１ステップ進めて、初期設定部２１に帰還させる。これにより、初期設定部２１、ＰＤＦ生成部２２、適合度評価部２３、再サンプリング部２４及びカウント部２５の間で、ループ処理が実現される。一方、最新の処理の時刻が最後のステップ（ｔ＝ｔｅｎｄ）であれば処理を終了し、処理完了報告を、バス４を介して、表示部３に出力する。
【００８９】
表示部３は、演算部２から出力される重み係数（すなわち粒子数）及び各粒子で表される時刻（ｔ）でのポテンシャル係数ｂ（ｔ）や直近のデータ数Ｍ（ｔ）の情報、処理完了報告を、視認可能な状態で、例えば液晶表示画面上に表示する。この際、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）等を表示装置に逐次表示し、ポテンシャル係数の経時的な変化をリアルタイムで観察することが可能である。
【００９０】
実施の形態３
次に、実施の形態３にかかる未来の時系列データがとる値の確率密度分布の予測方法について説明する。実施の形態１にかかる時系列データの変化点検出方法で説明したように、ＰＵＣＫモデル及び粒子フィルタ法を用いることにより、時刻（ｔ）におけるポテンシャル係数を求めることができる。さらに、このＰＵＣＫモデルを時間発展させることで、未来の市場価格がとる値を確率密度分布として得ることが可能である。
【００９１】
ｊ番目の粒子の時刻（ｔ＋１）における市場価格（ｔ＋１）は、式（１）を用いて、以下の式（３６）で表すことができる。
【数３８】
0006004284.tif 000039
【００９２】
図７は、実施の形態３にかかる粒子フィルタ法を適用したＰＵＣＫモデルの時間発展方法を示すフローチャートである。まず、実施の形態１にかかる時系列データの変化点検出方法で計算された時刻（ｔ）での粒子のパラメータ（ｂ（ｊ）ｉ（ｔ）、Ｍ（ｊ）（ｔ）、ｆ（ｊ）Ｐ（ｔ）を準備するステップＳ３１）。
【００９３】
次いで、実施の形態１にかかる時系列データの変化点検出方法で計算された時刻（ｔ）での真の市場価格Ｐ（ｔ）を準備する（ステップＳ３２）。
【００９４】
準備したパラメータ及び時刻（ｔ）での真の市場価格Ｐ（ｊ）（ｔ）を用いて、式（３６）を計算する。
【００９５】
そして、最新の処理の時刻（ｔ）が、ｔ＝ｔ＋Ｎに到達したかを判断する（ステップＳ３４）。最新の処理の時刻が最後のステップ（ｔ＝ｔ＋Ｎ）でなければ、時刻（ｔ）を１ステップ進めて（ステップＳ３５）、式（３６）での計算を再度実行する。最新の処理の時刻が最後のステップ（ｔ＝ｔ＋Ｎ）であれば、処理を終了する。
【００９６】
これにより、時刻（ｔ）でのパラメータ及び真の市場価格を用いて、式（３６）を時刻（ｔ＋Ｎ）まで時間発展させることができる。粒子フィルタ法を用いて推定した時刻（ｔ）でのｊ番目の粒子のＰＵＣＫモデルのパラメータセット（ポテンシャル係数ｂ（ｊ）ｉ（ｔ）、価格変動の中心を見積もるために必要なデータ数Ｍ（ｊ）（ｔ）、変動誤差ｆ（ｊ）Ｐ（ｔ））とすると、式（３６）を時刻（ｔ＋Ｎ）まで時間発展させた場合の時刻（ｔ＋Ｎ）での予想市場価格は、以下の式（３７）で表される。
【数３９】
0006004284.tif 000040
【００９７】
以上、各粒子について式（３７）を適用することにより、時刻（ｔ＋Ｎ）における市場価格Ｐ（ｊ）（ｔ＋Ｎ）の確率密度分布を得ることが可能となる。図８は、実施の形態３にかかる未来の時系列データがとる値の確率密度分布の予測方法により時刻（ｔ＋Ｎ）における市場価格の確率密度分布を求めた場合の確率密度分布の推移を示すグラフである。なお、図８では、Ｎ＝１００として、米ドル／日本円の為替レートを例として用いた。
図８では、式（３７）でＫ＝３の場合（すなわち、２次及び３次のポテンシャルが作用する）場合の時刻（ｔ）から時刻（ｔ＋Ｎ）までの確率密度分布について、第１四分位偏差ＱＤ１、中心値Ｖｃ及び第３四分位偏差ＱＤ３を表示した。なお、図８では、Ｎ＝１００として計算した。
【００９８】
図９は、実施の形態３にかかる未来の時系列データがとる値の確率密度分布の予測方法により時刻（ｔ＋Ｎ）における市場価格がとる値の確率密度分布を示すグラフである。なお、図９では、図８と同様に、Ｎ＝１００として、米ドル／日本円の為替レートを例として用いた。図９に示すように、時刻（ｔ＋Ｎ）での確率密度分布は、価格上昇側（分布曲線の頂部よりも右側の部分）での裾が広く、歪んだ分布を有することが確認できる。つまり、図９の例では、時刻（ｔ＋Ｎ）での市場価格は、上昇傾向の方がより強いと理解することができる。
【００９９】
以上のように、本実施の形態にかかる未来の時系列データがとる値の確率密度分布の予測方法によれば、３次以上のポテンシャル導入することで、価格変動の非線形挙動の影響を反映させることが可能となる。これにより、未来の市場価格の予測において、図９に示すような、左右非対称の歪んだ確率密度分布を生成することができる。これは、一般的なリスク評価などで用いる正規分布などの左右対称分布を用いる手法では実現できない、本実施の形態によりはじめて実現できる効果である。これにより、現実の金融市場等において生起する左右非対称の歪んだ確率密度分布を好適に予測することができ、より現実に即した時系列データがとる値の確率密度分布の予測方法を提供することができる。
【０１００】
以上、本実施の形態にかかる未来の時系列データがとる値の確率密度分布の予測方法によれば、粒子フィルタ法を用いて推定した時刻（ｔ）でのＰＵＣＫモデルのパラメータセットを基礎として、未来の市場価格がとる値の確率密度分布を得ることができる。また、式（３７）においては、Ｋ≧３の場合には３次以上のポテンシャルを導入できるため、市場価格変動の非線形な挙動を含めた予測が可能となる。これにより、より市場の実態に即した市場価格の予測が可能となる。
【０１０１】
本実施の形態にかかる未来の時系列データがとる値の確率密度分布の予測方法を為替や株式などの価格予測に適用すれば、為替や株式の価格変動のリスクを従来手法よりも高精度に見積もることができる。よって、従来の手法と比べて、より適正にリスクを見積もったオプションなどの金融商品の設計・開発に活用することができる。
【０１０２】
実施の形態４
次に、実施の形態４にかかる未来の時系列データがとる値の確率密度分布の予測プログラムについて説明する。図１０は、実施の形態４にかかる未来の時系列データがとる値の確率密度分布の予測プログラムが実行される未来の時系列データがとる値の確率密度分布の予測装置４００の構成を模式的に示すブロック図である。未来の時系列データがとる値の確率密度分布の予測装置４００は、実施の形態３にかかる未来の時系列データがとる値の確率密度分布の予測方法のステップＳ３１〜Ｓ３５を表現するアルゴリズムが記述されたプログラムが実行される装置である。図１０に示すように、未来の時系列データがとる値の確率密度分布の予測装置４００は、記憶部１、演算部２ａ、表示部３及びバス４を有する。記憶部１、表示部３及びバス４は、実施の形態２と同様であるので、説明を省略する。
【０１０３】
演算部２ａは、バス４を介して、記憶部１から、プログラム２０により算出されたパラメータ及び時系列データ及びプログラム４０を読み込む。演算部２は、読み込んだプログラム４０を実行し、未来の時系列データがとる値の確率密度分布の予測を行う。プログラム４０は、パラメータ読み込み部４１、データ読み込み部４２、算出部４３及びカウント部４４を有する。
【０１０４】
具体的には、パラメータ読み込み部４１は、プログラム２０で計算された時刻（ｔ）でのｊ番目の粒子のパラメータ（ｂ（ｊ）ｉ（ｔ）、Ｍ（ｊ）（ｔ）、ｆ（ｊ）Ｐ（ｔ））を読み込む。つまり、パラメータ読み込み部４１は、図７のステップＳ３１に対応する処理を行う。なお、プログラム２０で計算された時刻（ｔ）での粒子のパラメータは、例えば記憶部１に格納され、必要に応じてパラメータ読み込み部４１により読み出される。
【０１０５】
データ読み込み部４２は、プログラム２０で計算された時刻（ｔ）での真の市場価格Ｐ（ｊ）（ｔ）を読み込む。つまり、データ読み込み部４２は、図７のステップＳ３２に対応する処理を行う。なお、プログラム２０で計算された時刻（ｔ）での真の市場価格Ｐ（ｊ）（ｔ）は、例えば記憶部１に格納され、必要に応じてデータ読み込み部４２により読み出される。
【０１０６】
算出部４３は、図７のステップＳ３３に対応する処理を行う。すなわち、式（３６）を用いて、設定されている時刻よりも１ティック先の市場価格を算出する。
【０１０７】
カウント部４４は、図７のステップＳ３４及びＳ３５に対応する処理を行う。すなわち、最新の処理の時刻が最後のステップ（ｔ＝ｔ＋Ｎ）であるかを判定する（図７のステップＳ３４での処理）。そして、最新の処理の時刻が最後のステップ（ｔ＝ｔ＋Ｎ）でなければ、時刻（ｔ）を１ステップ進めて、算出部４３に帰還させる。これにより、式（３６）の計算が繰り返され、時刻（ｔ＋Ｎ）に到達した場合には、式（３７）に示す式と同等の計算件が得られる。最新の処理の時刻が最後のステップ（ｔ＝ｔ＋Ｎ）であれば、処理を終了し、算出結果に基づいて、時刻（ｔ＋Ｎ）における市場価格の確率密度分布を、表示部３に出力する。
【０１０８】
表示部３は、演算部２から出力される時刻（ｔ＋Ｎ）における市場価格の確率密度分布を、視認可能な状態で、例えば液晶表示画面上に表示する。この際、時刻（ｔ）から時刻（ｔ＋Ｎ）までの間の市場価格の確率密度関数を表示装置に逐次表示することも可能である。
【０１０９】
以上より、本実施の形態によれば、実施の形態４にかかる未来の時系列データがとる値の確率密度分布の予測プログラム及び予測装置を具体的に実現することができる。
【０１１０】
実施の形態５
次に、実施の形態５にかかる未来の時系列データがとる値の確率密度分布の他の予測方法について説明する。実施の形態３では、ＰＵＣＫモデルを時間発展させることで、未来の市場価格がとる値を確率密度分布として得ることが可能であることを説明した。つまり、時間発展の開始時点である時刻（ｔ）でのポテンシャル係数ｂ（ｔ）を用いて時間発展を行う。なお、本実施の形態では、ポテンシャルの次数ｉ、粒子の番号ｊについては、表示を省略している。
【０１１１】
実施の形態３では、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）の値が、将来の時刻（ｔ＋１）から（ｔ＋Ｎ）の間でも同様だと仮定して、将来を予測している。従って、実施の形態３では、短期的な価格変動の時間変化しか予測することができない。これに対し、本実施の形態のように、ポテンシャル係数ｂ（ｔ）の時間発展の方程式が得られれば、将来の時刻（ｔ＋１）から（ｔ＋Ｎ）の間のポテンシャル係数ｂ（ｔ）の変化を予測することができる。これにより、より正確な将来の時刻t+Nでの価格分布が得られる。従って、長期的な価格変動の時間変化を予測する場合、ポテンシャル係数ｂ（ｔ）の時間発展の方程式を推定することは重要である。例えば、時刻（ｔ＋１）でのポテンシャル係数ｂ（ｔ＋１）と時刻（ｔ）でのポテンシャル係数ｂ（ｔ）との関係を、以下の式（３８）で定義する。
【数４０】
0006004284.tif 000041
なお、右辺第１項（偏微分項）の関数Ｇ（ｂ）は、時刻（ｔ）までのポテンシャル係数の時間発展を記述する関数である。ｆｂ（ｔ）は、ノイズ項である。
【０１１２】
マイナス符号を除く右辺第１項（偏微分項）を、ポテンシャル係数の時間変化に作用するポテンシャルλとすると、式（３８）は式（３９）に変形できる。
【数４１】
0006004284.tif 000042
【０１１３】
式（３９）より、時刻（ｔ＋１）でのポテンシャル係数ｂ（ｔ＋１）は、以下の式（４０）で表される。
【数４２】
0006004284.tif 000043
式（４０）より、時刻（ｔ＋１）でのポテンシャル係数ｂ（ｔ＋１）は、係数（１−λ）が作用し、ノイズ項ｆｂ（ｔ）が加えられる。
【０１１４】
ここで、式（４０）で示すポテンシャル係数を用いて、例えば２次のポテンシャルが作用する場合の市場価格分布の予測について説明する。２次のポテンシャルが作用する場合の１ティックの間での真の市場価格の変動は、式（２）、（４）及び（５）より、以下の式（４１）で表すことができる。
【数４３】
0006004284.tif 000044
【０１１５】
図１１は、実施の形態５にかかる粒子フィルタ法を適用したＰＵＣＫモデルの時間発展方法を示すフローチャートである。図１１では、図７と比べて、ステップＳ３３をステップＳ５２に置換している。また、ステップＳ３２及びＳ３５とステップＳ５２との間に、ステップＳ５１が追加されている。
【０１１６】
ステップＳ５１は、式（４０）を用いて、ポテンシャル係数ｂ（ｔ＋１）を計算する。次いで、ステップＳ５２では、準備したパラメータ及び時刻（ｔ）での真の市場価格Ｐ（ｊ）（ｔ）を用いて、式（４１）を計算する。その他のステップについては、図７と同様であるので、説明を省略する。
【０１１７】
図１２は、現在の市場価格を数理モデルに代入して、粒子番号ｊ毎に時間発展させ、その経過を粒子毎に表示したグラフである。ここでは、例として、開始価格を１００とし、１４５０ティック後までの市場価格の分布を予測した。図１２に示すように、本実施の形態によれば、ポテンシャル係数を時間発展させることで、将来のポテンシャル係数の変動をも予測し、市場価格の長期トレンドをより正確に予測することが可能である。
【０１１８】
長期トレンドの予測は市場価格分布として与えられるので、ある将来の時点における価格変動についてのリスクの見積もりが可能である。これにより、本実施の形態によれば、将来の為替や株式の価格変動のリスクを実施の形態３よりも高精度に見積もることができる。よって、更に適正にリスクを見積もったオプションなどの金融商品の設計・開発に活用することができる。また、多数の金融機関が自己資産のリスクを計測するために用いているValue at Riskなどのリスク管理の指標を改善することができる。
【０１１９】
なお、本実施の形態で説明したポテンシャル係数の時間発展は例示に過ぎない。よって、ポテンシャル係数の時間発展は、式（４０）で示す例に限られない。
【０１２０】
実施の形態６
次に、実施の形態６にかかる未来の時系列データがとる値の確率密度分布の予測プログラムについて説明する。図１３は、実施の形態６にかかる未来の時系列データがとる値の確率密度分布の予測プログラムが実行される未来の時系列データがとる値の確率密度分布の予測装置６００の構成を模式的に示すブロック図である。未来の時系列データがとる値の確率密度分布の予測装置６００は、実施の形態５にかかる未来の時系列データがとる値の確率密度分布の予測方法のステップＳ３１、３２、Ｓ５１、Ｓ５２、Ｓ３４、Ｓ３５を表現するアルゴリズムが記述されたプログラムが実行される装置である。未来の時系列データがとる値の確率密度分布の予測装置６００は、未来の時系列データがとる値の確率密度分布の予測装置４００の変形例である。図１３に示すように、未来の時系列データがとる値の確率密度分布の予測装置６００は、図１０のプログラム４０をプログラム６０に置換した構成を有する。
【０１２１】
プログラム６０は、プログラム４０の算出部４３を算出部６２に置換し、更にポテンシャル算出部６１を追加した構成を有する。ポテンシャル算出部６１は、図１１のステップＳ５１に対応する処理を行う。すなわち、式（４０）を用いて、ポテンシャル係数を計算する。算出部６２は、図１０のステップＳ５２に対応する処理を行う。すなわち、式（４１）を用いて、設定されている時刻よりも１ティック先の市場価格を算出する。未来の時系列データがとる値の確率密度分布の予測装置６００のその他の構成は、未来の時系列データがとる値の確率密度分布の予測装置４００と同様であるので、説明を省略する。
【０１２２】
以上より、本実施の形態によれば、実施の形態５にかかる未来の時系列データがとる値の確率密度分布の予測プログラム及び予測装置を具体的に実現することができる。
【０１２３】
なお、本発明は上記実施の形態に限られたものではなく、趣旨を逸脱しない範囲で適宜変更することが可能である。例えば、半導体製造工程において装置から生じるゆらぎを伴う時系列データに対して本発明を適用することによって、工程の処理が安定的に進んでいるかどうかをリアルタイムで計測し、何らかの原因で不安定化した場合には速やかに異常を検出することが可能となる。この場合、例えば、市場価格を装置が出力する管理・計測データに置き換えることで、同様の解析を行うことが可能である。
【０１２４】
上述の実施の形態では、式（４）に示すように、時刻（ｔ）でのポテンシャル係数ｂ（ｔ）を２次以上のポテンシャルＷ（ｑ,ｔ）の係数として説明した。例えば、現実の金融市場などでは、２次ポテンシャルでは表現できない、更なる急激なトレンドの変動が有ることが知られている。そのため、時刻（ｔ）でのポテンシャルＵ（ｑ,ｔ）を２次に限定せず、高次ポテンシャル項を含む関数形を導入して高次ポテンシャル項の係数を推定することにより、急激なトレンド変動に追随することができる。特に、３次ポテンシャルを導入することにより、金融市場における暴騰又は暴落について、２次ポテンシャルと比べてより高速に追随することが可能である。また、暴騰（すなわち、トレンド上昇）又は暴落（すなわち、トレンド下降）の方向性についても、併せて解析することが可能となる。
【０１２５】
また、２次ポテンシャルモデルを用いて、時刻（ｔ）におけるｂ（ｔ）およびＭ（ｔ）の同時確率分布等を推定することにより、複数の時間スケールでのディーラーの異なる戦略が推定可能となる（マルチスケールＰＵＣＫモデル）。さらに、３次以上の高次ポテンシャル係数を推定することで、トレンドを追うディーラーの戦略、および、価格を一方的に変動させている集団的行動の可視化が可能となる（高次ポテンシャルを含むマルチスケールＰＵＣＫモデルの展望）。
【０１２６】
この出願は、２０１１年９月８日に出願された日本出願特願２０１１−１９６５１２及び２０１２年３月７日に出願された日本出願特願２０１２−５０８１９を基礎とする優先権を主張し、その開示の全てをここに取り込む。
【産業上の利用可能性】
【０１２７】
本発明は、為替市場や株式市場での市場価格の時間変動の分析又は予測、管理・計測データなどの他の時系列データの時間変動の分析又は予測に適用することができる。
【符号の説明】
【０１２８】
１  記憶部
２  演算部
３  表示部
４  バス
２０、４０、６０  プログラム
２１  初期設定部
２２  ＰＤＦ生成部
２３  適合度評価部
２４  再サンプリング部
２５  カウント部
４１  パラメータ読み込み部
４２  データ読み込み部
４３、６２  算出部
４４  カウント部
６１  ポテンシャル算出部
２００  時系列データの変化点検出装置
４００、６００  未来の時系列データがとる値の確率密度分布の予測装置
Ｓ１〜Ｓ６、Ｓ１１〜Ｓ１４、Ｓ２１〜Ｓ２４、Ｓ４１〜Ｓ４５、Ｓ５１、Ｓ５２  ステップ
【図１】
0006004284.tif 000047
【図２】
0006004284.tif 000048
【図３】
0006004284.tif 000049
【図４】
0006004284.tif 000050
【図５】
0006004284.tif 000051
【図６】
0006004284.tif 000052
【図７】
0006004284.tif 000053
【図８】
0006004284.tif 000054
【図９】
0006004284.tif 000055
【図１０】
0006004284.tif 000056
【図１１】
0006004284.tif 000057
【図１２】
0006004284.tif 000058
【図１３】
0006004284.tif 000059