JPH09233234A

JPH09233234A - 周期的なパターンの監視

Info

Publication number: JPH09233234A
Application number: JP8346854A
Authority: JP
Inventors: Ross Mitchell Cox; ミッチェルコックスロス
Original assignee: AT&T Corp
Current assignee: AT&T Corp
Priority date: 1995-12-29
Filing date: 1996-12-26
Publication date: 1997-09-05
Also published as: CA2191270A1; US5768284A; CA2191270C; EP0782313A3; EP0782313A2

Abstract

(57)【要約】【課題】本発明は、複数の顧客に対してサービスを提
供する電話ネットワークにおける記録の誤りを検出する
ための装置および方法に関する。【解決手段】その方法は、顧客に対するサービスの提
供に関連している電話ネットワークによって出力される
周期的なパターンを持っているような、情報の記録のス
テップと、所定の期間より長く持続する周期的パターン
における妨害を識別するステップと、識別された妨害が
しきい値を超えた時に警報を生成するステップとを含ん
でいる。また、本発明は以前の期間中に記録された情報
を利用して、まだ発生していない期間においてネットワ
ークによって出力されることになる情報の予測値を生成
するために、モデルを使うことをも特徴とし、そのモデ
ルは周期的なパターンにおける妨害を表している複数の
ステートを含んでいて、各ステートが持続時間の異なる
妨害を表している。

Description

【発明の詳細な説明】

【０００１】

【発明の属する技術分野】本発明は通信ネットワークに
おける呼出しの詳細記録などの周期的なパターンの監視
に関し、そして特に、ほぼリアルタイムで記録の誤りを
検出するためのシステムに関する。

【０００２】

【従来の技術、及び、発明が解決しようとする課題】電
話のシステムは提供されるサービスに対して顧客に適切
な料金請求をするために、電話のアクティビティについ
ての詳細を記録しなければならない。例えば、その記録
される情報としては、その呼出しの接続時刻、呼び出さ
れた番号および切り離し時刻などがある。それによっ
て、呼出しの持続時間および呼び出された番号に基づい
た料金を反映する請求書を発行することができる。

【０００３】呼出しの詳細記録監視（または「ＣＤＲ
Ｍ」）は、呼出しの詳細記録が正しく作られていること
を保証するために使われる。そのような記録を作る際の
誤りの結果、電話のサービスのプロバイダーにとって大
きな歳入損失が発生する可能性がある。これらの誤りは
呼出しの提供時の誤り、電話システム間での料金請求情
報の伝送における誤り、料金請求番号スクリーニング・
データベースにおける誤り、およびハードウェアおよび
ソフトウェアの誤りによって発生する可能性がある。そ
の誤りはそれらが発生した後、長い時間が経ってから下
流の料金請求システムにおいて発見されるのが普通であ
る。誤りを早期に検出することによって、問題を是正す
ることができ、歳入損失を最小にすることができる。従
来の技術の欠点は本発明によって克服される。本発明は
その好適な実施例において、呼出しの詳細記録のカウン
トの合計値に異常な変動がないかどうかを監視するこ
と、およびそのような変動が観察された時に警報を発す
ることによって、ほぼリアルタイムで電話のネットワー
ク内の記録の誤りを検出するためのシステムおよび方法
を提供する。

【０００４】

【課題を解決するための手段】１つの態様においては、
本発明は複数の顧客に対してサービスを提供する電話ネ
ットワークの中の記録の誤りを検出するための装置およ
び方法に関し、その方法は、顧客に対するサービスの提
供に関連している電話ネットワークからの出力情報を記
録し、所定の期間より長い期間持続する周期的なパター
ンにおける乱れを識別し、そしてその識別された乱れが
しきい値を超えた時に警報を生成することを含んでい
る。

【０００５】もう１つの態様においては、本発明は複数
の顧客に対してサービスを提供する電話ネットワークに
おける記録の誤りを検出するための装置および方法を含
み、その方法は、顧客に対するサービスの提供に関連し
ている第１の期間に電話のネットワークにおいて出力さ
れる情報を記録し、その情報は周期的なパターンを持っ
ているような、記録のステップと、第１の期間中に記録
された情報を利用して、まだ発生していない第２の期間
におけるそのネットワークによって出力される情報の予
測値を生成するためにモデルを使い、そのモデルは周期
的なパターンにおける乱れを表している複数のステート
を含んでいて、各ステートは異なる持続期間を持ってい
る乱れを表しているような、モデルの使用のステップ
と、第２の期間中の顧客に対するサービスの提供に関連
している電話ネットワークによって出力される情報を獲
得し、その偏差を生成するために第２の期間の間に得ら
れた情報と予測値とを比較するステップと、その偏差を
使って各ステートを更新するステップと、そのステート
の選択された１つがしきい値を超えた時に警報を生成す
るステップとを含んでいる。

【０００６】好適な実施例においては、そのステートの
うちの選択された１つは、他のステートより長く持続す
るステートである。その記録される情報は所定の期間中
にネットワーク内で切り離された呼出しに対するネット
ワーク内で接続された呼出しの比率であってよい。ま
た、その記録される情報は所定の期間にわたるネットワ
ーク内において行なわれた呼出しのボリュームであって
もよい。予測値を生成するためにモデルを使用するステ
ップは、第１の期間において記録された情報を利用する
ステップの他に、その第１の期間の前に記録された情報
を利用するステップを含んでいることが好ましい。

【０００７】その好適な方法は、第２の期間中に記録さ
れた情報およびその第２の期間より前に記録された情報
を利用し、まだ発生していない第３の期間に対するネッ
トワークによって出力されることになる情報の期待され
る内容を予測して利用するステップと、その第３の期間
において顧客に対するサービスの提供に関連している電
話ネットワークによって出力される情報を記録するステ
ップと、そしてその第３の期間に関係付けられている予
測のステップの結果と第３の期間の間に記録された情報
とを比較して第２の偏差を生成するステップとを含む。
その方法は各ステートを更新するために第２の偏差値を
使うステップと、そのステートの選択された１つがしき
い値を超えた時に警報を生成するステップとをさらに含
んでいてもよい。その方法のステップはネットワークに
よって提供される複数のサービスに対して行なわれるよ
うにすることができる。

【０００８】この方法の各ステップは、ネットワークに
よって最近出力されたデータと最近の出力データに対す
る予測値との比較に対応している偏差値を絶えず生成す
るために繰り返されることが好ましく、そのステートは
選択された更新済みのステートがしきい値を超えた場合
に、そのステートが繰り返し更新され、そして警報が発
せられるようにする。そのモデルはカルマン・フィルタ
を使って生成することができる。

【０００９】さらにもう１つの面においては、本発明は
複数の顧客に対してサービスを提供する電話ネットワー
クにおける記録の誤りを検出するための方法に関し、そ
の方法は、顧客に対するサービスの提供に関連している
電話ネットワークからの出力情報を記録し、その情報が
周期的なパターンを持っているような、記録のステップ
と、記録された情報に基づいて１つのモデルを生成し、
そのモデルは第１の期間の周期的パターンにおける乱れ
を表している第１のステートおよび第２の期間の周期的
パターンにおける乱れを表している第２のステートを持
ち、その第２の期間は第１の期間より長いようなモデル
生成のステップと、そして第２のステートにおける乱れ
がしきい値を超えた時に警報を生成するステップとを含
んでいる。

【００１０】上記の方法は、まだ発生していない期間に
対して電話ネットワークによって出力されることになる
情報の予測値を生成するためにそのモデルを使うステッ
プと、その期間中に電話ネットワークによって出力され
る情報を記録するステップと、そしてその予測のステッ
プの結果をその期間中に記録された情報と比較して偏差
値を生成するステップとをさらに含んでいてもよい。そ
の偏差値は各ステートを更新するために使うことができ
る。さらにもう１つの態様においては、本発明は複数の
顧客に対してサービスを提供する電話ネットワークにお
ける記録の誤りを検出するための方法を特徴とし、その
方法は、顧客に対するサービスの提供に関連している電
話ネットワークから出力される情報を記録し、その情報
が周期的なパターンを持っているような記録のステップ
と、第１の期間の周期的パターンにおける乱れを識別す
るステップと、第２の期間の周期的パターンにおける乱
れを識別し、その第２の期間が第１の期間より長いよう
な識別のステップと、第２の期間の乱れがしきい値を超
えた時に警報を生成するステップとを含んでいる。

【００１１】さらにもう１つの態様においては、本発明
は複数の顧客に対してサービスを提供する電話ネットワ
ークにおける記録の誤りを検出するための装置および方
法を特徴とし、その方法は、（ａ）第１の期間において
顧客に対するサービスの提供に関連した電話ネットワー
クによって出力される情報を記録するステップ、（ｂ）
その第１の期間において記録された情報に基づいてモデ
ルを生成し、そのモデルは第１の持続期間の情報におけ
る乱れを表している第１のステートおよび第２の持続期
間の情報における乱れを表している第２のステートのデ
ータを含んでいて、その第２の持続期間は第１の持続期
間より長いような、モデル生成のステップ、（ｃ）その
モデルを使って将来の期間に対するネットワークによっ
て出力されるべき情報の予測値を生成するステップ、
（ｄ）将来の期間の発生時に、顧客に対するサービスの
提供に関連している電話ネットワークによって出力され
る情報を記録するステップ、（ｅ）ステップ（ｃ）で生
成された予測情報とステップ（ｄ）で記録された情報と
を比較して偏差情報を生成するステップ、（ｆ）前のス
テップで生成された偏差に基づいて各ステートを更新す
るステップ、（ｇ）ステートのうちの選択された１つの
ステートがしきい値を超えた時に警報を生成するステッ
プ、および（ｈ）ステップ（ｃ）〜（ｄ）を繰返し実行
してステートのうちの選択された１つがしきい値を超え
た時は常に警報を絶えず生成するステップを含んでい
る。好適な実施例においては、ステートのうちの選択さ
れた１つ（ステップ（ｇ））は第２のステートである。

【００１２】

【発明の実施の形態】本発明の１つの実施例がＡＴ＆Ｔ
の４ＥＳＳ^TM交換機の内部で、これらの交換機に対する
トラフィック・パターンにおける偏差を解析することに
よって、呼出しの詳細記録の問題点を検出し、そして警
報を生成する環境の中で記述される。（４ＥＳＳ^TMの設
計の電子的交換システムについての詳細はＢｅｌｌＳ
ｕｓｔｅｍＴｅｃｈｎｉｃａｌＪｏｕｒｎａｌ（Ｂ
ＳＴＪ）の１９７７年９月第５６巻７号の１０１７ペー
ジに記載されている一連の記事の中で説明されており、
それは参照によってここに組み込まれている。

【００１３】４ＥＳＳ^TMは呼出しの詳細記録における誤
りを検出するために解析することができる各種のデータ
を生成する。例えば、ＭＬＳＳ（マシン・ロードおよび
サービスの集計）として知られている信号はその４ＥＳ
Ｓ^TMの中で接続された呼出しの数を表している。「トレ
ーサ」記録として知られている別の量は、完了した呼出
しの数の測度である。「接続された」呼出しは他の側に
おいて受信された呼出しである。「完了」呼出しは切り
離されたか、あるいは終了した呼出しである。接続され
た呼出しは切り離された呼出しと一対一に対応している
必要がある。というのは、接続された呼出しはすべて実
質的に切り離されなければならないからである。同様
に、以前に接続されていない呼出しが切り離されること
はあり得ない。従って記録の誤りがない場合、ＭＬＳＳ
およびトレーサの記録はそれぞれのタイミングが違って
いるだけとなる。従って、与えられた任意の時間に対す
るこれらの２種類の記録の比は１に近い値になることが
期待される。というのは、ほとんどの接続された呼出し
はそれらが起動された１時間の間に切り離され、そして
それらが起動されたあと１時間以内に完了した接続され
た呼出しはその次の時間において完了する新しい接続さ
れた呼出しの数だけ部分的にオフセットされることにな
るからである。

【００１４】通常の動作（記録が失われていない）の
間、比率の値は毎日および毎週の周期を持つパターンに
従う。このパターンは交換機によって、そしてサービス
によって変化する可能性がある（そして普通はそのよう
になる）。図１は４週間の期間にわたる比率の時系列の
一例である。そのパターンの規則性は明らかである。図
１の中の比率の値は、例えば、１日のうちの業務開始の
時間帯のように、接続されている呼出しの数が切り離さ
れている呼出しの数を超える時は１より大きくなる。こ
れと逆の状態、例えば、業務指向のサービスのための午
後の遅い時間帯に発生するように、切り離される呼出し
の数が接続される呼出しの数を超える時には、この比率
は１以下に落ちる。

【００１５】本発明の最初の実施例は比率時系列におい
て異常な変動を生じるデータの記録の中に潜んでいる誤
りを検出するために、図１に示されているような比率時
系列を監視するために使われる。その好適な実施例は現
在の時系列の知識（すなわち、その現在の値）およびそ
の時系列の歴史の知識に基づいて、時系列の予測値を生
成することによって監視を実行する。その時、この予測
値が実際のデータの値から差し引かれて誤差の項を形成
し、それがカルマン・フィルタの理論の骨組みを経由し
て、異常な乱れの存在が数学的に解析される。

【００１６】図２は比率の監視を実行するために使われ
る第１の実施例を示している。それは減算器１０を含ん
でおり、この減算器は監視される比率時系列ν_t の現在
値（すなわち、現在の時刻ｔにおける比率時系列の値）
およびその時系列の以前に計算された予測値ｆ_t を入力
として受け取る。時系列ν_t と予測値ｆ_t との間の差が
誤差信号ｅ_t として出力される。誤差信号ｅ_t は測定の
更新モジュール１２に入力され、そのモジュールはその
比率時系列の「状態空間モデル」のいくつかのステート
の現在値の推定値であるベクトル量ｓ＾_t を更新するよ
うに働く。その更新されたベクトル量ｓ＾_t はカルマン
・フィルタの理論を使って生成される。カルマン・フィ
ルタは１９６０年にＲ．Ｅ．カルマン（Ｋａｌｍａｎ）
によって発見された推定技法であり、機械的なシステム
の数学的モデリングの技法に習熟している人によってよ
く知られている。例えば、カルマン・フィルタはトラッ
キングおよびナビゲーションの分野（例えば、人工衛星
と関連した）において広く使われている。

【００１７】本発明によって使われる状態空間モデルが
図３に示されている。図３は実際の（すなわち、物理的
な）システムではなく、図３の中のブロックは物理的な
構造に対応するのではなく、値ν_t を生成するために仮
定されるシステムのモデルに対応する。図３の中の唯一
の「実際の」要素は出力信号ν_t だけであり、それは比
率時系列の観測値を表す。図３に示されている他の「信
号」およびブロックはカルマン・フィルタの公式化のた
めに必要であり、その推定値を生成するモデルを含んで
いる。カルマン・フィルタの理論は推定のパラメトリッ
クな理論である。それはその推定値を算出するために、
その信号発生プロセスのパラメトリックな記述（モデ
ル）が必要である。ここで比率の監視のために使うモデ
ルは図３に示されているモデルである。このモデルは比
率の値ν_t が３つの他の値の合計であると仮定する。そ
の１つは比率の測定の週ごとの通常の変動を表している
項（図３の中のｓ_t ^norm 、そして２つの妨害プロセスｓ
_t ²およびｓ_t ³通常の変数ｓ_t ^norm は別の妨害変数ｓ_t ¹に
よってドライブされ、ｓ_t ¹はさらにノイズ変数ｎ_t ¹によ
ってドライブされる。妨害変数ｓ_t ²およびｓ_t ³はノイズ
変数ｎ_t ²およびｎ_t ³によってそれぞれドライブされる。

【００１８】入力ｎ_t ¹、ｎ_t ²およびｎ_t ³は、互いに統計
的に独立であり、従って、妨害変数ｓ_t ¹、ｓ_t ²、および
ｓ_t ³も統計的に独立である。長期の妨害変数ｓ_t ¹は通常
の比率トラフィック変数ｓ_t ^norm をドライブすることに
注意する必要がある。つまり、このことは長期の妨害変
数ｓ_t ¹は通常の変数ｓ_t ^norm のパターンに永久的に影響
する可能性があることを示している。少なくとも１つの
妨害プロセスと、通常のプロセスとの間のこのリンクが
ない場合、このモデルは通常のプロセスが時間によって
は変化しない固定のパターンを持つことを示す。長期の
妨害プロセスが通常のプロセスに影響することを許すこ
とによって、通常のプロセスがその時系列の中の永久的
な変動に合わせることができる。

【００１９】変数ｓ_t ^norm 、ｓ_t ¹、ｓ_t ²、およびｓ_t ³は
この状態空間モデルの中のステートである。以下に記述
されるように他のステートもあるが、それらはステート
ｓ_t ^norm 、ｓ_t ¹、ｓ_t ²、およびｓ_t ³を更新する目的のた
めだけに維持されなければならない補助的なステートで
ある。ｓ_t ¹のステート推定値ｓ＾_t ¹、すなわち長期の妨
害変数は以下にさらに説明されるように、警報を作り出
すために使われる量である。中期および短期の変数のス
テートの推定値は、この実施例において警報を作るため
には使われない、ただし、それらはｓ＾_t ¹を含む他のス
テートの推定値を更新するために必要である。図３の中
でｓ_t によって示されている完全な比率トラフィック変
動変数は３つのプロセス、ｓ_t ^norm 、ｓ_t ²、およびｓ_t ³
の和である。上記のように、変数ν_t は観測された比率
トラフィック変数、すなわち、ＣＤＲＭによって実際に
記録されるデータである。変数Ｎ_t ^obsは１つのノイズ変
数であり、それは他のノイズ入力ｎ_t ¹、ｎ_t ²、およびｎ
_t ³とは相関がなく、また、ゲインρ_t によって乗算され
る時に完全な比率トラフィック変動変数ｓ_t に加算され
て観測変数ν_t を形成する。Ｎ_t ^obsはｓ_t の観測値を乱
している観測ノイズとみなされる。

【００２０】図２に戻って、ステートの推定値ｓ＾_t
は、将来の期間ｔ＋１に対する予測モジュール１４の中
の時系列の各ステートの予測値を生成するために使われ
る。これは現在の値およびその履歴および他のステート
に基づいて各ステートの将来の値を外挿することによっ
て行なわれる。モジュール１４の出力は時刻ｔ＋１にお
ける各ステートに対する予測値を表している量ｓ￣_t+1
である。集計モジュール１６はｓ￣_t+1 における各種の
ステートに対する予測値の和をとって、時刻ｔ＋１にお
ける時系列に対する単独の予測値ｆ_t+1 を生成する。次
に、遅延モジュールが、時刻ｔ＋１が現在の時刻になる
までシステムを待たせ、以前に生成された予測値ｆ_t+1
を現在時刻ｔに対する予測値であるｆ_t に等しく設定す
る。上記のように、減算器１０はその後、時系列の予測
値ｆ_t と現在値ν_t との差を求める。

【００２１】時系列の現在値が利用できない場合（すな
わち、欠落している場合）、欠落データ検出器１８がｓ
￣_t+1 を時間更新モジュール１４の入力へ戻して次の期
間に対する予測値を形成する。言い換えれば、データν
_t が利用できない場合は比較のステップを実行すること
ができないので、システムは現在の予測値を使って次の
期間に対する予測値を生成する（すなわち、その外挿値
を次の期間に対して延長する）。

【００２２】上記のように、図１の中で示されているよ
うな推定値を生成するため、および記録の誤りが発生し
たと考えられる時を示すために使われる警報を生成する
ためにカルマン・フィルタの技法が使われる。次に、カ
ルマン・フィルタがこれらの結果を得るために本発明に
おいてどのように使われるかを説明する。推定のための
カルマン・フィルタの方法は、再帰的な推定の方法であ
る。すなわち、それは各時間において、保持されている
有限の数の量の関数として、比率トラフィック変数の更
新された推定値を計算する。維持されている量は比率の
推定値と同時に更新される。上記のようにカルマン・フ
ィルタの理論においては、予測値を生成するために「モ
デル」が使われる。本発明の好適な実施例の中で使われ
るモデルは図３に示されている。

【００２３】カルマン・フィルタのモデルの内部で更新
される量は「ステート」と呼ばれ、それらの推定値（カ
ルマン・フィルタによって更新され、維持される量）は
「ステート推定値」と呼ばれ、そしてカルマン・フィル
タのモデルは「状態空間モデル」と呼ばれる。観測変数
に対する統計的モデルの与えられた選択に対して、ステ
ートに対する可能な多くの選択がある。すなわち、状態
空間モデルに対する多くの選択がある。これらはすべて
観測されたデータのモデリングに適している。しかし、
本発明の目的は観測データそのものを推定または予測す
ることではなく、そのデータの中での何らかの特徴また
は偏差を推定または検出することであるので、１つのス
テートとしてこの機能を含んでいる状態空間モデルがカ
ルマン・フィルタ理論の観点から好ましい。これはこの
状態空間モデルに基づいたカルマン・フィルタが、その
観測変数の推定値と一緒にこの特徴の推定値を維持し、
更新することになるからである。

【００２４】カルマン・フィルタの方法によって使われ
る状態空間モデルは、複数のプロセスからの結果として
得られる変数をモデル化することができる。比率データ
に対して使われるモデルは、同時に発生する３種類の異
なるプロセスが重畳されたものである。１つは通常のト
ラフィック変動プロセスであり、他の２つは妨害プロセ
スである。この重畳プロセスは完全な比率トラフィック
・プロセスと呼ばれ、結果の変数は完全な比率トラフィ
ック変数と呼ばれる。さらに、以下に詳しく説明され
る、通常のトラフィック変動プロセスをドライブする１
つの妨害プロセスがある。

【００２５】この３種類の妨害プロセスは持続時間のレ
ベルが異なる妨害をモデル化する。これらの各妨害はス
テートによってモデル化される。カルマン・フィルタに
よって作られ、これらのステートに対応するステートの
推定値におけるしきい値設定によって、３種類の異なる
各妨害に対して３種類の異なる警報が作られる。本発明
の好適な実施例においては、１つのタイプの妨害だけが
警報のために使われ、他の２つの警報は単純に無視され
る（ただし３種類のすべてを使うこともできる）。

【００２６】実際の時系列とその予測時系列との間の大
きな差の理由が基本的に２種類ある。（１）呼出しの詳
細記録誤りを示す可能性のある時系列データにおける妨
害があり得る。そして（２）データにおける大きなラン
ダムな変動があり得る。本発明はこのうちの最初のもの
だけを警報することを求める。

【００２７】状態空間モデルはｓ_t ^norm 、ｓ_t ²、および
ｓ_t ³によって示される３つのステートの和として完全な
比率トラフィック変数を記述する。これらのステートは
通常の比率トラフィック変動プロセスおよび中期および
短期の妨害変数をそれぞれ表す。さらに第３の妨害変数
ｓ_t ¹および１９１の追加の補助ステートがあり、その補
助ステートはこれらの各ステートの個々の更新を記述し
ている１９５の式を書くために必要である。そのステー
トは単独の１９５次元のベクトルｘ_t のコンポーネント
であると考えるのが便利である。もっと一般的には、任
意の状態空間モデルにおけるステートの更新は行列の表
記で最も簡単に表される。この理由で、行列の代数学
は、状態空間モデル、カルマン・フィルタを記述するた
めに使われる自然な数学的言語である。従って、完全な
比率トラフィック変動プロセスに対して使われるモデル
についての以下の記述は行列の表記を採用する。

【００２８】上で説明したように、比率警報アルゴリズ
ムに対する状態空間モデルの完全な記述を与えるために
行列の表記が必要である。しかし、行列の表記はモデル
の入力−出力の関係、すなわち、或るノイズ入力ｎ_t ¹、
ｎ_t ²、およびｎ_t ³と、その出力（観察された）変数ν_t
との間の関係だけを記述するには不要である。入力−出
力の関係を示すことは比率トラフィック変数ν_t の統計
的な記述を与えることと等価であり、そしていわゆるボ
ックス・ジェンキンズ（Ｂｏｘ−Ｊｅｎｋｉｎｓ）の予
測理論の中で使われているもっとよく知られたＡＲＩＭ
Ａモデルによって完全に与えることができる（１９９４
年１１月３０日に申請されて、本発明出願の被譲渡人に
対して譲渡された米国特許出願第０８／３４６，８０８
号（ここで参照によって組み込まれている）参照）。予
測のためのＡＲＩＭＡのモデルおよびボックス・ジェン
キンズの方法はカルマン・フィルタによって使われてい
る状態空間モデルよりももっと広く知られているので、
ここでは比率のモデルの入力−出力の関係のＡＲＩＭＡ
の記述を示す。そこで、状態空間モデルの行列の記述は
次のようになる。

【００２９】長期の妨害変数ｓ_t ¹、中期の妨害変数
ｓ_t ²、および短期の妨害変数ｓ_t ³は入力ｎ_t ¹、ｎ_t ²、お
よびｎ_t ³についての３つのＡＲＩＭＡモデルの出力であ
る。これらのモデルは明示的に次の式によって与えられ
る。

【数１】 α_w 、α_d 、α_q 、およびβ_w 、β_d 、β_q は以下で与
えられる値の定数である。式１の中のｓ_t+1 ¹、ｓ_t+1 ²、
およびｓ_t+1 ³における時間のインデックス添字ｔ＋１は
ノイズ・プロセスｎ_t ¹、ｎ_t ²、ｎ_t ³における添字から１
時間進んだものであり、下記の同じモデルの状態空間記
述におけるこれらの変数における添字と一致している。
また、式１の中のＢは逆方向の時間シフト演算子であ
る。普通の比率トラフィック変数ｓ_t ^norm は、さらに次
の差分方程式によって長期の妨害変数ｓ_t ¹に関係付けら
れる。

【数２】

【００３０】これは単にｓ_t+1 ¹が通常のトラフィック変
数ｓ_t+1 ^norm と、その１週間前の値ｓ_t-167 ^norm との間
の差に等しいことを示しているだけである。式（２）は
通常のプロセスｓ_t+1 ^norm のＡＲＩＭＡの記述を限定す
るものではない。なぜなら、入力Ｓ_t+1 ^lは、時刻相関で
あり、ホワイトではないからである。しかし、式（２）
と式（１）を組み合わせることによって、ｓ_t ^norm はｎ
_t ¹の項で次のように表される

【数３】これはｎ_t ¹がホワイト・プロセスであるので、ｓ_t ^norm
のＡＲＩＭＡの記述である。式（１）および式（３）
は、変数ｓ_t+1 ^norm 、ｓ_t+1 ²、ｓ_t+1 ³は、３つのボック
ス・ジェンキンズ型のＡＲＩＭＡモデルの出力であるこ
とを示している。これらの式に対して、比率トラフィッ
ク・モデルの入力−出力の記述を完成させるために、次
のもう２つの式を追加する。

【数４】ここで＜Ｎ_t ^o＞は＜ｎ_t ¹＞、＜ｎ_t ²＞、および＜ｎ_t ³＞
とは独立のホワイト・ノイズ・プロセスであり、その標
準偏差がσ_t ^o2 （おそらく時間によって変わる）であ
る。式（１）、（４）、および（５）は比率トラフィッ
ク変動プロセスν_t を記述するために使われる状態空間
モデルの入力−出力の関係を要約している。このモデル
は次の３つの点で一般のボックス・ジェンキンズのＡＲ
ＩＭＡモデルとは異なっている。

【００３１】１．比率トラフィック変数ν_t の観測値
（推定アルゴリズムに対する入力）は、「信号」ｓ_t お
よび独立のノイズ・プロセスＮ_t ^oの和としてモデル化さ
れている。＜ｓ_t ＞は、ここでノイズを含んでいる観測
値＜ν_t ＞を使って推定しようとしている「真の」トラ
フィック変動プロセスと考えることができる。Ｎ_t ^oは観
測ノイズと呼ばれる。２．信号ｓ_t は３つのボックス・ジェンキンズ・モデ
ルの出力ｓ_t ^norm 、ｓ_t ²、ｓ_t ³、の和としてモデル化さ
れ、それぞれが異なるノイズ入力＜ｎ_t ¹＞、＜ｎ_t ²＞、
＜ｎ_t ³＞によってドライブされる。

【００３２】これで完全な比率トラフィック・プロセス
に対するＡＲＩＭＡモデルの項での記述が完了する。こ
のプロセスに対する状態空間モデルが次に記述される。
次の説明の中で、ｘ_t は状態空間モデルのステート・ベ
クトル（すなわち、そのコンポーネントがモデルのステ
ートである１９５次元のベクトル）を示す。ｘ_tのコン
ポーネントのうちの４つが上記の４つのステートｓ_t
^norm 、ｓ_t ¹、ｓ_t ²およびｓ_t ³である。残りの１９１個
のコンポーネントは、これらの４つのステートを更新す
る目的で厳密に維持されている補助ステートである。こ
れらの補助ステートのほとんどは、ｓ_t ^norm およびｓ_t ¹
の単なる遅延されたバージョンである。これらのステー
トの定義は状態空間モデルそのものによって暗黙に与え
られている。

【００３３】式（１）、（２）および（４）によって記
述されるプロセスに対する状態空間モデルは次の式によ
って与えられている。

【数５】ここで行列Ｃ、Ａ、Ｆは以下の表１で示されている次元
および非ゼロ要素が散在している行列である。

【表１】

【００３４】ここでＮ_t ³＝［ｎ_t ¹、ｎ_t ²、ｎ_t ³］^T は、
コンポーネントがｎ_t ¹、ｎ_t ²、ｎ_t ³である列ベクトルで
ある。式（６）の中の第２の式はｘ_t のコンポーネント
ｘ_t ⁱ；ｉ＝１、．．．、１９５の項で次のように書く
ことができる。

【数６】

【００３５】式（７）と（８）とを連結することによっ
て次の式が得られる。

【数７】従って、ｘ_t+1 ²⁵ は式（１）の変数ｓ_t+1 ¹と同じであ
る。さらに、式（１０）と（９）とを連結することによ
って次の式が得られる。

【数８】従って、ｘ_t+1 ¹⁹³は式（３）の変数ｓ_t+1 ^norm と同じで
ある。ｘ_t+1 の残りの２つのコンポーネントｘ_t+1 ¹⁹⁴お
よびｘ_t+1 ¹⁹⁵は式（１）の変数ｓ_t+1 ²とｓ_t+1 ³と同じで
あることが容易に分かる。従って、システム（６）の信
号ｓ_t+1 に対する最終の式は次のようになる。

【数９】

【００３６】そして、システム（６）の信号ｓ_t+1 は式
（４）の信号ｓ_t+1 と同じである。すなわち、この２つ
のシステムの入力−出力の関係は同じである。式（６）
は信号ｓ_t の状態空間記述を与える。観測されたトラフ
ィック変数ν_t の式のシステムは次の観測ノイズＮ_t ^oが
加わることだけが式（６）と異なっている。

【数１０】これは式（１）、（３）、（４）および（５）で与えら
れているシステムの状態空間記述である。式（１）の
中、および表１の中に現われている定数α_w 、α_d 、α
_q 、β_w 、β_d 、およびβ_q の実際の値は次のようにな
っている。

【数１１】

【００３７】α_w 、α_d 、α_q はすべて０．２^1/x+1 の
形になっている。ここでｘはそれぞれ１日の中の時間
数、１週間の中の日数および４分の１日の中の時間数で
ある。定数α_d の値は式（７）の中の演算子［１-α_d
Ｂ］^-1のインパルス応答を２４サンプル＝１日の後、そ
の初期値の２０％まで減衰するように選ばれている。α
_qに対しても同様である。α_w の値は式（８）の中の
［１−α_w Ｂ²⁴］^-1のインパルス応答が７×２４＝１６
８サンプル＝１週間の後、その初期値の２０％にまで減
衰するように選ばれている。β_d の値はβ_d ［１−α_d
Ｂ］^-1のインパルス応答のエネルギーＥが次の式で与え
られているように選ばれている。

【数１２】これは１に等しい。これによってｘ_t+1 ¹の分散が１に等
しくなり、さらにそれによってｘ_t+1 ²⁵ の分散が式
（８）から分かるようにほぼ１に等しくなる。同様なコ
メントがβ_w およびβ_q に対しても適用される。妨害プ
ロセスｘ_t ²⁵ ＝ｓ_t ¹、ｘ_t ¹⁹⁴＝ｓ_t ²、およびｘ_t ¹⁹⁵＝ｓ
_t ³の分散が互いに等しくなり、それぞれの観測されたト
ラフィック変数ν_t に及ぼす影響が等しくなるようにす
ることが望ましい。これによって変数ν_t における妨害
の観測値が、その大きさによるのではなく、それがどの
程度長く持続するかのベースで厳格に長期、中期、また
は短期として分類されることになる。

【００３８】状態空間モデル（１１）に対するカルマン
・フィルタは、それが完全なトラフィック変動変数ｓ_t
についてのその推定値ｓ＾_t ＝ｘ＾_t ¹⁹³＋ｘ＾_t ¹⁹⁴＋ｘ
＾_t ¹ ⁹⁵を更新するのと同時に、３種類の妨害ステートに
対する推定値ｘ＾_t ²⁵ 、ｘ＾_t ¹⁹⁴およびｘ＾_t ¹⁹⁵を更新
することになる。この監視のアルゴリズムは持続時間の
最も長い妨害のうちのステートｘ_t ²⁵ の推定値であるｘ
＾_t ²⁵ の正規化されたバージョンが、しきい値を超えた
時は常に警報を生成する。これが行なわれる正確な方法
について以下に説明する。

【００３９】比率の監視に対するカルマン・フィルタの
方法と、ボックス・ジェンキンズの方法との間の問題形
成における３つの主な違いをまとめると次のようにな
る。１．カルマン・フィルタの方法は観測された比率の値ν
_t+1 が信号ｓ_t+1 そのものではなく、観測ノイズ＜Ｎ
_t+1 ^o＞によって乱された信号、すなわち、ν_t+1＝Ｓ_t+1
＋Ｎ_t+1 ^oであると仮定している。２．カルマン・フィルタのモデルにおける「信号」ｓ
_t+1 は３つのボックス・ジェンキンズ型のモデルの出力
の和ｓ_t+1 ＝ｓ_t+1 ^norm ＋ｓ_t+1 ²＋ｓ_t+1 ³としてモデル
化される。３．カルマン・フィルタの方法において警報のために使
われる変数は予測の誤りｅ_t+1 ではなく、そのモデルの
内部ステートｘ＾_t+1 ²⁵ の推定値ｘ＾_t+1 ²⁵ である。

【００４０】本発明によって使われる警報を作るための
方法は、上記のように、完全な比率トラフィック変動プ
ロセスのためのモデルを選定し、このモデルを使ってそ
の観測された比率の値に関係付けられた１つまたはそれ
以上の量を推定する方法である。警報はこれらの推定値
のどれかまたはいくつかの組合せが或るしきい値を超え
た時に生成される。上記のように、警報のための推定値
は、ｓ＾_t ¹＝ｘ＾_t ²⁵のスケールされたバージョン、長
期の妨害変数のステートの推定値となる。ｓ＾_t ¹をスケ
ーリングするために使われる定数は、標準偏差の推定
値、あるいはこの推定値の分散の平方根である。本発明
のカルマン・フィルタ理論の使用目的は、推定値ｓ＾_t ¹
およびその分散の両方を計算することである。

【００４１】次に、状態空間モデル（１１）からステー
ト・ベクトルｘ_t の推定値、そして特に、このベクトル
のコンポーネントｘ_t ²⁵ を得る方法について説明する。
ｘ_tは次の式によって更新される。

【数１３】ｘ＾_t はこれに似た式によって更新することができそう
である。しかし、信号ノイズＮ_t ^oの値は観測不可能であ
り、式（１２）を使ってｘ＾_tを更新するために観測さ
れた量によって置き換えられなければならない。Ｎ_t ^oを
式（１２）の中で置き換えるために選ばれる候補は、予
測ｅ_t+1 ＝ν_t+1 −ＣＡｘ＾_t である。というのは、そ
れは観測可能であり、推定値ｘ＾_t がｘ_t にどれだけ近
いかという測度を提供するからである。しかし、予測の
誤差ｅ_t+1 はスカラであり、そして次元の首尾一貫性の
ために、それは式（１２）の中に置き換えるには、ベク
トルＫ_t+1 によって乗算されなければならない。これを
行なうと、式（１２）は次のようになる。

【数１４】

【００４２】これは式（１２）に類似している便利なク
ラスの推定式を提供する。ただし、式（１３）の中のｘ
＾_t+1 はＦＮ_t ^oではなく、Ｋ_t+1 （ν_t+1 −ＣＡｘ＾
_t ）の項によってドライブされている。しかし、推定値
ｘ＾_t+1 は列ベクトルＫ_t+1 が指定されるまではまだユ
ニークには決定されない。カルマン・フィルタの理論は
Ｋ_t+1の選択において使う規則を提供する。この理論の
中心の結果は、観測値ν_t、ν_t-1 、．．．、が与えら
れた時の線形システム（１１）のステート・ベクトルｘ
_t の最小二乗法による推定値ｘ＾_t が、カルマン・フィ
ルタの更新式（付録Ａ参照）によって与えられるベクト
ル＜Ｋ_t ＞のシーケンスによる更新式（１３）を満足す
ることである。

【００４３】ステートの推定値ｘ＾_t を更新するための
式の他に、カルマン・フィルタの理論は理論的な誤差の
共分散行列Ｐ_t の現在値を計算するための再帰型の式も
提供する。（これらの式は付録Ａの中で（Ａ．３）、
（Ａ．４）、および（Ａ．６）として現われる。）こ
の行列は状態空間モデル（１１）が正確であるという仮
定の下で、状態誤差ベクトルｘ＾−ｘ_t の共分散行列に
等しい。量Ｐ_tは観測されたデータν_t 、ν
_t-1 、．．．とは無関係であり、そしてモデルによって
のみ決定される。これはステート誤差ベクトルｘ_t −ｘ
＾_t の実際の共分散行列の推定値ではあり得ない理由を
説明している。また、付録Ａの中の式（Ａ．３）、
（Ａ．４）、および（Ａ．６）から導かれるようなカル
マン・フィルタのアルゴリズムの実際の実行の前に、Ｐ
_t をオフラインで計算できる理由をも説明している。

【００４４】誤差共分散行列Ｐ_t ^25,25の２５番目の対角
要素は、この場合も状態空間モデル（１１）が正確であ
るという仮定の下で、警報変数ｘ＾_t ²⁵ の分散である。
警報のためのしきい値を作るため、あるいは等価的に
は、ｘ＾_t ²⁵ の分散の値が１になるように正規化するた
めに、ｘ＾_t ²⁵ に実際の分散が必要である。ｘ＾_t ²⁵ の
実際の分散、あるいは実際の誤差共分散行列を求めるに
は、Ｐ_tはデータに依存している量によってスケーリン
グされなければならない。これを行なう方法を説明する
前に、Ｐ_t およびＫ_t が、モデルのパラメータ、特に、
ノイズの分散σ^o2およびσⁱ²およびＰ_t のＰ_toの初期値
そのものにどのように依存するかをもう少し詳細に説明
する。

【００４５】付録Ａの中の式（Ａ．３）、（Ａ．４）、
および（Ａ．６）から検証されるように、σ^o2、σⁱ²、
およびＰ_toのすべてを同じ定数ｃで乗算することはＰ_t ^o
をｃで乗算することの効果はあるが、Ｋ_t の値を変化さ
せない。一方、Ｋ_t の値は２つのノイズ分散の比ρ² ＝
σⁱ²／σ^o2に依存する。この比はステートの推定値ｘ＾
_t の過渡的な挙動を決定する。ρ² ＜＜１であった場
合、Ｋ_t ＝０およびカルマン・フィルタは変化している
比率のトラフィック・パターンに対してゆっくりと適応
するが、それ自身ゆっくり変化しているパターンに対す
る非常に正確な定常状態の推定値を与える。

【００４６】一方、ρ² ＞＞１であった場合、このフィ
ルタはトラフィック・パターンの変化に対して迅速に適
応するが、ゆっくりと変化しているパターンに対する推
定値よりも不正確な値を与える。好適な実施例において
採用される方法は、比ρ² を「固定する」ことである。
比ρ² を固定化するための手順は次の通りである。先ず
最初に、データのパターンにおける永久的な変化に対し
てこのアルゴリズムがどの程度迅速に適応するようにし
たいかを前の経験に基づいて決定する。「適応する」と
いうことの意味は、そのアルゴリズムは警報を生成して
停止し、新しいパターンを確立された通常のトラフィッ
ク・パターンとしてみなすことを意味する。次に、その
データのパターンの中に永久的な変化（例えば、ステッ
プ）を導入された１つの範囲のデータについて、そのア
ルゴリズムを実行することによって、ρ² のどの値がこ
の適応のレートを生じるかを実験的に求める。好適な実
施例は現在約４週間でステップの変化に対して適応す
る。ρ² の固定化によって、このモデルの中でσ^o2＝
１、そしてσⁱ²＝ρ² と取り、Ｐ_t を計算するためにこ
のモデルに対してカルマン・フィルタの更新式を使うこ
とによって、そのアルゴリズムの所望の適応のレートが
与えられる。これによって、次の式が得られるように、
乗算の定数γ_t ²までのＰ_t の正しい値が得られる。

【数１５】

【００４７】式（１４）の右辺は定義によって、このス
テート・ベクトル推定値の実際の共分散行列である。式
（１４）からγ_t ²を求めるために、誤差ｅ_t ＝ｘ＾_t −
ｘ_tの少なくとも１つのコンポーネントを、観測された
比率のトラフィックのデータから直接に推定できるよう
にしたい。というのは、これによって真の誤差共分散行
列Ｅ［ｅ_t ｅ_t ^T］の対応するコンポーネントを推定する
ことができるからである。しかし、真のステート・ベク
トルｘ_t のコンポーネントを観測しないので、これは不
可能である。しかし、ε＝ν_t −Ｃｘ＾_t によって与え
られる「適合誤差」εを観測する。そしてこの量の分散
が次の式で求められるということを使って、εの観測値
とＰ_t についての知識とを組み合わせて、γ_t ²を導くこ
とができる。

【数１６】

【００４８】次に、γ_t ²および式（１４）を使ってＥ
［（ｘ＾_t −ｘ_t ）（ｘ＾_t −ｘ_t ^T］のコンポーネント
を完全に求めることができる。この行列の（２５、２
５）番目のコンポーネントはσ² _et ²⁵ ＝Ｖａｒ［ｘ＾_t
²⁵ ］であり、これは警報の統計的なｘ＾_t ²⁵ の分散で
あり、従って、次の式によって与えられる正規化された
分散β_tはその分散が１となり、警報を作るための適切
な量となる−付録Ｄ参照。

【数１７】

【００４９】ここで、完全な比率警報システムにするた
めに核のカルマン・フィルタ理論に対して行なった変更
および追加について説明する。特に、以下に説明される
のは、完全な比率警報アルゴリズムによって解決されな
ければならない次の条件が、カルマン・フィルタの理論
のコンテキストの中でどのように満足されるかというこ
とである。・飛び離れた値のリジェクト・ステップの取扱い・誤差分散の推定・警報の形成・初期化・欠落データの扱い

【００５０】先ず最初に飛び離れた値のリジェクション
について考える。最も普通の推定技法は飛び離れた値に
対して明示的な備えを行なっていないが、飛び離れた値
のリジェクションは既に調査中の母集団において実行さ
れていると仮定し、そして結果のサンプルのすべての測
定値が含められ、そして推定値を作り出すために等しく
重み付けられていることを仮定する。従って飛び離れた
値のリジェクションはこの技法を適用する前に実行され
る時「アドオン」として扱われる。これは時間または空
間のいずれかの定期的な間隔において分布しているサン
プルの測定値を必要とするボックス・ジェンキンズの予
測法のような方法に対して１つの問題を提供する可能性
がある。これらの方法によって、除外される極端な値は
何らかの他の値によって置き換えられなければならな
い。しかし、カルマン・フィルタによって使われる状態
空間モデル（１１）は極端な値を明示的に考慮するのに
十分に一般的であるので、それらを仮定の値で置き換え
ずにｘ＾_tの推定値を作り出すために使うことができ
る。

【００５１】本発明のこの利点は重要である。というの
は、飛び離れた値と同様に極端な値はデータのパターン
における永久的な変化を示す可能性があるからである。
これらの変化は推定のアルゴリズムによって実質的に注
意されなければならない。式（１１）を使って極端な値
をモデル化する手段は極端な値が観測される時間ｔに対
する観測ノイズＮ_t ^oの分散σ^o2を増加させることによ
る。カルマン・フィルタは、観測値ν_t 、ν
_t-1 、．．．および状態空間モデルおよび観測値の分散
およびステート・ノイズ・プロセスの知識に基づいてｘ
_t の推定値ｘ＾_t の最適値を各時間ｔにおいて生成する
ことを思い出されたい。さらに、そのモデルおよびノイ
ズの変動は時間と共に変わることが許されている。これ
は非常に大きな柔軟性を提供するが、上記の実施例にお
いては、時間と共に変化することを許しているモデルの
唯一のコンポーネントは観測ノイズの瞬時変動だけであ
る。任意の時間において、この値を任意の値に設定する
ことは自由である。特に、１つの極端な値が観測された
時刻ｔ＋１においてそれを大きな値に設定した場合、フ
ィルタは観測された比率の値ν_t+1 と観測ノイズに対す
る予測値Ｃｘ￣_t+1 との間にほとんどの食違いを割り当
て、そのステートの推定値の更新を計算する際にその観
測値を軽い値で重み付ける。これは付録Ａの中の更新式
から明確に分かる。式（Ａ．４）の中のσ_t+1 ^o2 の増加
によって、Ｋ_t+1 のコンポーネントのサイズが減少し、
それによってさらに新考案のプロセス（定義については
付録Ａ参照）の現在値ν_t+1 （＝ν_t+1 −Ｃｘ￣_t+1 ）
がｘ＾_t+1 に対する更新の式（Ａ．５）の中でより小さ
く重み付けられることになる。

【００５２】飛び離れた値のリジェクションの好適な方
法に対する正確な式は付録Ｂに与えられている。特に、
付録Ｂはこの目的に対して使われる「極端な値」の定義
以外に、極端な値が観測される時間における観測ノイズ
の変動を増加させるために使われる式を提供する。ステ
ップの取扱いは、どの比率監視アルゴリズムも観測され
たデータ・パターンにおける永久的な変化に対して適応
できるようにするための必要条件である。ステップの取
扱いの条件は飛び離れた値のリジェクションの条件と衝
突する。というのは、前者は推定のアルゴリズムが極端
なデータの値に注目することを要求し、一方、後者はそ
れが極端な値に注意を払わないことを要求するからであ
る。

【００５３】ステップの取扱いは本発明において別のス
テップとして実装されない。代わりに、飛び離れた値の
リジェクションが扱われる方法の中に組み込まれてい
る。前に説明されたように、極端な値が飛び離れた値で
あって、その推定のプロセスから除外されるかどうかに
関してイエスまたはノーの判定を行なわずに、値の新し
いパターンが持続する時間が長いほど自動的に増加する
量によって重み付けられる。これによって、フィルタは
ステップおよび他のパターンの変化に対して実質的に適
応することができる。また、カルマン・フィルタの方法
はパターンの変化の後のデータに対して新しいモデルを
適合させる必要はない。というのは、それは固定された
モデルについて監視するのに十分頑健であるからであ
る。

【００５４】適合誤差の変動の推定値は、γ_t ²の推定値
と等価である。というのは、これらの２つは式（１５）
によって関係付けられているからである。これらの量の
１つの推定値が、カルマン・フィルタのアルゴリズムに
よって算出される理論的な誤差共分散行列Ｐ_t の値をス
ケーリングするために必要である。さらにこれは警報変
数ｘ＾_t ²⁵ の実際の（スケーリングされない）分散の推
定値を算出するために必要であり、正規化された警報変
数β_t （分散が１である）が付録Ｄの中の式（Ｄ．１）
に従って形成されるために必要である。

【００５５】ノイズにおける一定のオフセットを推定す
るために単純なカルマン・フィルタを使ってγ_t ²を推定
する。その式が付録Ｃに示されている。極端な値が存在
しない場合、このフィルタは係数が一定である単純なフ
ィルタへ収れんする。定常状態に収れんするγ_t ²の推定
機能として単純なフィルタではなく、カルマン・フィル
タとしての推定機能を実装する理由は、この推定値が比
率データの中の極端な値の場合に爆発してしまうことを
避けるためである。これは極端な値を欠落データとして
扱うことにより、そしてそれらをγ_t ²の推定値から完全
に消去することによって行なわれる。これについての詳
細は付録Ｃの中で説明されている。警報の形成を考え
て、付録ＡおよびＢの中の式に従って推定値ｘ＾_t+1 ²⁵
が形成された後、そしてスケーリングの定数γ_t+1 が付
録Ｃの中の式に従って形成された後、前者は後者によっ
てスケーリングされ、警報変数β_t+1 を作り出す。β
_t+1 の大きさがしきい値を超えた場合、警報が発せられ
る。このための式が付録Ｄに示されている。

【００５６】ここでシステムの初期化を考える。このア
ルゴリズムの最初の繰返しの間に、付録Ａ、Ｂ、および
Ｃの中に現われている更新式（Ａ．１）〜（Ａ．６）、
（Ｂ．１）〜（Ｂ．４）、および（Ｃ．１）〜（Ｃ．
６）の中のｔ＝ｔ₀ および、これらの式の右辺に現われ
る量ｘ＾_t0、Ｐ_t0、γ_t0 ² およびｐ_t0の量が初期化され
なければならない。これらの量に対する初期値は次の通
りである。ｘ＾_t0＝０Ｐ_t0＝１０００・Ｉ γ_t0 ²＝０ｐ_t0＝１０００ｘ＾_t0およびγ_t0 ² の値は多かれ少なかれ任意である
が、Ｐ_t0およびｐ_toは大きい値（＞＞１）に選定され、
システムが最初の観測された比率の値ν_t0+1を大きく重
み付け、観測されたデータに対して速やかに（２、３時
間程度で）観測されたデータに適応するように選定され
る（付録Ａの中の式（Ａ．３）、（Ａ．４）、および
（Ａ．５）参照）。同様に、ｐ_t0はγ_t0がそのデータに
迅速に適応するように大きく選定される。また、付録Ｂ
の中の式（Ｂ．２）および付録Ｃの中の式（Ｃ．４）〜
（Ｃ．６）を取り囲んでいる“ｉｆ”ブロックの中の条
件から分かるように、飛び離れた値のリジェクションが
オフされ、γ_t ²がこの最初の繰返しに対してオンされ
る。初期化は実行の最初においてのみ行なわれる。シス
テムがスタートした後、再度初期化する必要はない。

【００５７】欠落データの扱いが次に説明される。比率
時系列の値は各時間において実数値である。しかし、Ｃ
ＤＲＭはトレーサの値に対するＭＬＳＳの比として解釈
するのが難しい（それは普通非常に大きな値であるの
で）値を時々渡すことがある。これらの値が大きくなる
原因となっている飛び離れた値のリジェクションまたは
欠落データの処理について、ＣＤＲＭがそれ自身のフォ
ームを実装できると考えられる。本発明の観点から、Ｃ
ＤＲＭがこれらの欠落値を充填しなかったが、代わり
に、これらの時間に対応している比率を負の数に等しく
設定することによって、その時間がデータを含んでいな
いことをその情報と一緒に渡していることが好ましい。
これはカルマン・フィルタが働いて、欠落データあっ
て、データ・ストリームの中に仮定の値を挿入すること
によって得られるものは何もないような時系列において
さえも最適だからである。

【００５８】本発明は欠落しているデータのある時系列
について使えるように変更することができる。比率の値
が利用できる任意の時間に対するこのアルゴリズムの操
作は、以前にあったものと同一となる。データが欠落し
ている残りの時間に対して、そのアルゴリズムは次のよ
うに変更されなければならない。カルマン・フィルタの
更新式（付録Ａ）の中で、ステップ（Ａ．２）および
（Ａ．４）は省略され、ステップ（Ａ．５）および
（Ａ．６）は次で与えられる式（Ａ．５’）および
（Ａ．６’）によって、それぞれ置き換えられる。

【数１８】

【００５９】ステップ（Ｂ．１）は省略され、式（Ｂ．
２）を取り囲んでいる“ｉｆ”ブロックの中の条件が偽
であると考えられるので、△σ_t+1 ^o2 が０に等しく設定
される。最後に、ステップ（Ｃ．１）および（Ｃ．２）
が省略されて、 γ_t+1 ²＝γ_t ² という式が式（Ｃ．２）、（Ｃ．３）および（Ｃ．４）
を囲んでいる“ｉｆ”ブロックの代わりに、更新中のγ
_t+1 ²に対して使われる。本発明は呼出しの詳細記録およ
び特に比率の時系列の監視のコンテキストにおいて説明
されてきた。しかし、前の説明は単に本発明の一実施例
に過ぎず、本発明は他の環境においても広く応用され
る。この分野の技術に熟達した人には、多くの変更およ
び適用があることは明らかである。例えば、ＣＤＲＭの
環境において、本発明は記録の誤りを検出するために他
のタイプのデータを分析するために使うことができる。
いわゆる「ボリューム」データは上記のカルマン・フィ
ルタの方法によって解析することができる。

【００６０】さらに、本発明はＣＤＲＭ以外の環境にも
適用される。一例として、本発明は電話回線を違法にト
ラップすることによって、あるいは盗んだテレフォンカ
ードを使うことによってサービスを盗むことから生じ
る、呼出しのボリュームの異常な増加について警告する
ことによって、電話のネットワークにおける不正の検出
を行うために適用される可能性がある。

【００６１】以下の付録Ａ−Ｇと称する部分は、本明細
書の一部を構成するものとして、本文中から参照される
記述である。

【００６２】付録Ａ：カルマン・フィルタの更新式本論文の中で使われる信号発生モデルは式（１１）で与
えられている。

【数１９】ここで次の関係がある。

【数２０】カルマン・フィルタはこのモデルに関連している次の４
つの量を更新することに関与する。ステートの推定値：ｘ＾_t＝Ｅ［ｘ_t｜ν，ν_t-1，．．．］ステートの予測値：ｘ￣_t＝Ｅ［ｘ_t｜ν_t-1，ν_t-2，．．．］ステートの推定誤差の共分散：Ｐ_t＝Ｅ［（ｘ＾_t−ｘ_t）（ｘ＾_t−ｘ_t）^T］ステートの予測誤差の共分散：Ｈ_t＝Ｅ［（ｘ￣_t−ｘ_t）（ｘ￣_t−ｘ_t）^T］厳密に言うと、条件付きの期待値としてｘ＾_t およびｘ
￣_t のための上記の式は信号および観測ノイズ・プロセ
ス＜Ｎｓ_t ＞および＜Ｎ_t ^o＞がガウス分布である場合に
のみ成立する。しかし、ｘ＾_t およびｘ￣_t は常に最小
二乗法での推定値となる。

【００６３】時間ｔ＋１において、新しい比率の測定値
ν_t+1 が利用できるようになり、上記の量は次のシーケ
ンスですべてのｔ≧ｔ₀ ：に対して更新される。

【数２１】式（Ａ．２）は、新考案のシーケンスν_t を定義してい
る。これはセクション２で述べたように、過去の観測値
ν_t 、ν_t-1 ・・・とは独立のホワイト・ノイズのプロ
セスであり、従って、ｘ￣_t+1 とは無関係である。ν
_t+1 は時間ｔ＋１において観測値ν_t+1 と共に利用でき
るようになっているステートｘ_t+1 に関する新しい情報
を含んでいる。従って、カルマンのゲイン行列Ｋ_t+1 の
値は、この新しい情報を重み付ける方法を決定すること
になる。ν_t+1 が飛び離れた値である場合、Ｋ_t+1 を小
さい値にしたい。これについてはこれを行なうための装
置について説明した。

【００６４】ステップ（Ａ．１）および（Ａ．３）は、
次の事実を使うことによって非常に効率的にプログラム
することができる。すなわち、Ａは、ステート推定ベク
トルｘ＾_t にＡを乗算することがｘ＾_t のコンポーネン
トのほとんどを左へ１単位シフトする効果を持つこと、
すなわち、１、２５、１９３、１９４、１９５以外のす
べてのｉに対して（Ａｘ＾_t）ⁱ＝ｘ＾_t ⁱ⁺¹である。そこ
で、ｘ＾_t とｘ￣_t を表している論理的な配列がメモリ
の中で同じ物理的な配列の場所を占める場合、それぞれ
ｘの別々の部分行列をポイントしている２つのポインタ
を設け、少数の追加の乗算および加算を加えることによ
って、上記のステップ（Ａ．１）を「イン・プレース
で」実行することができる。これを行なうためのＣライ
クなコードをいくつか以下に示す。

【００６５】

【数２２】

【００６６】この中のｆｏｒループは、本発明の監視ア
ルゴリズムのメイン・ループである。それは変数「ｈｏ
ｕｒ」（時間のインデックスを表している）をインクリ
メントし、そして付録Ａ〜Ｄの中のすべての式全体を通
してループする。式（Ａ．１）に対応しているステップ
が上記の行（１）〜（７）に与えられている。ここでｘ
は式（Ａ．１）の中のｘ＾_t およびｘ￣_t+1 の両方を表
している配列である。もっと正確には、ループの（ｔ＋
１）番目の繰返し（ｈｏｕｒ＝ｔ）において、そして式
（１）の実行の前に、ｘはｘ＾_tの値を含んでおり、一
方、式（７）の実行の後に、ｘはｘ￣_t の値を含んでい
る。式（１）の実行の前に、整数「ｐｔｒｌ」は長期の
妨害変数ｓ_t ^l＝ｘ_t ²⁵ のステートの推定値ｘ＾_t ²⁵ の値
を含んでいるｘの配列の内部の場所をポイントする。そ
して整数「ｐｔｒ２」は、通常のトラフィック変動変数
ｓ_t ^norm ＝ｘ_t ¹⁹³のステートの推定値ｘ＾_t ¹⁹³の値を含
んでいるｘの中の場所をポイントする。上記の式（１）
から分かるように、ｐｔｒｌはループを通過するたびに
１だけインクリメントされる。ただし、上限（＝２４）
に等しいループに続くパスにおいては、その下限（＝
１）へ回り込む。同様にｐｔｒ２は上限（＝１９２）に
等しいパスに続くパスまでインクリメントされ、その
後、下限（＝２５）へ回り込む。

【００６７】上記のステップ（Ａ．３）の「イン・プレ
ース」実装は少し複雑であるが、行列の乗算を使って直
接実装するのに比べてより効率的である。ステップ
（Ａ．３）を実行するためのＣライクなコードを以下に
示す。

【数２３】

【００６８】この中のｆｏｒループは上記と同じであ
り、変数ｐｔｒ１およびｐｔｒ２も同じである。上記の
偽（８）〜（２１）は行（１）〜（７）の後に挿入され
るべきである。これらの行はステップ（Ａ．３）を実行
する。ここでＰは（Ａ．３）における論理行列Ｐ_t およ
びＨ_t+1 の両方を表している１つの行列（ポインタの配
列に対するポインタ−行列のＣバージョン）である。ル
ープの（ｔ＋１）番目の繰返し（ｈｏｕｒ＝１）におい
て、そして上記の（８）の実行の前に、ＰはＰ_tの値を
含んでおり、一方（２１）の実行の後、ＰはＨ_t+1 の値
を含んでいる。行列Ｐ_t およびＨ_t+1 は、対称行列であ
り、従ってＰ［ｉ］［ｊ］はｉ≧ｊに対して宣言され定
義されるだけで済む。従って、メモリの使用量を減らす
ために、Ｐの部分対角要素に対してのみメモリを割り当
てる。これは上記のループの前にある文の中で見ること
ができる。^* 行（８）〜（２１）の中で、インデックス
ｉが常にメモリの割り当てられた部分をポイントしてい
ることに注意しなければならない。従って、例えば、行
（１１）および（１２）は、ループのための次の単独行
にマージすることはできない。

【数２４】というのは、ｉ＞ｐｔｒ１に対して、Ｐ［ｐｔｒ１］
［ｉ］はメモリの割り当てられていない部分を参照する
からである。また、行（１４）、（１５）および（１
６）を単独のループとして書くことができず、（１８）
および（１９）も同様である。

【００６９】ステップ（Ａ．１）および（Ａ．３）の
「イン・プレース」実装のために、物理的配列ｘの内部
の与えられた値をポイントしているインデックスは、一
般に論理配列ｘ＾_t またはｘ￣_t の内部の同じ値をポイ
ントしているインデックスとは同じではない。すなわ
ち、一般に１〜１９２の範囲のｉの値に対してループの
トップにおいて、ｘ［ｉ］≠ｘ＾_t ⁱ⁺¹である。同様に、
この同じ範囲内のｉおよびｊ（ｊ≦ｉ）の値に対するル
ープのトップにおいてＰ［ｉ］［ｊ］≠Ｐｔ^i+1、j+1
である。代わりに、物理的配列に対するインデックス
と、論理的配列に対するインデックスとの間の対応は、
次のように定義されているｉ＝０，．．．、１９４での
インデックスされたシーケンスＰｉの項で次のように与
えられる。

【数２５】次に、ループのトップにおいて、０〜１９２の間のｉお
よびｊ（ｊ≦ｉ）に対して、ｘ［Ｐ_i ］＝ｘ＾_t ⁱ⁺¹そし
てＰ［Ｐ_i ］［Ｐ_j ］＝Ｐ_t ^i+1、j+1である。このことを
知っている必要はない。プログラムはインデックス間の
正しい関係を自動的に維持する。また、本発明のプログ
ラムの中の物理的配列の任意のインデックスと、それに
対応している、それが表す論理配列との間にも同じ対応
関係が成立する。特に、上記の式（Ａ．４）〜（Ａ．
６）の中に現われる論理的配列Ｋ_t+1 およびＰ_t+1 を表
している物理的配列に対して成立する。

【００７０】付録Ｂ：飛び離れた値のリジェクションに
対する補助方程式次の式は、飛び離れた値のリジェクションを実装するた
めに更新の式（Ａ．１）〜（Ａ．６）に対して追加され
る式である。（Ａ．３）と（Ａ．４）との間に挿入する

【数２６】

【００７１】付録Ａの中で以前に定義されていなかった
上記の変数は本発明のモデル（１１）の項において次の
ように解釈される。新考案のシーケンスの分散： σ_tν²＝Ｅ［ν_t ²］極端な値のための観測ノイズ分散における増加： △σ
_t ^o2 観測ノイズの分散： σ_t ^o2 警報変数に対するスケーリング定数γ_t ² 式（Ｂ．１）は新考案のシーケンスの分散σ_t+1ν²を更
新し、そして式（Ｂ．３）、（Ｂ．４）と組み合わせら
れたｉｆブロック（Ｂ．２）は、観測された比率の値ν
_t+1 が極端な値であった場合に、観測ノイズＮ_t+1 ^oと新
考案のシーケンスとの分散を増加させる。これは上記の
ようにステートの推定値ｘ＾_t+1に及ぼす現在の観測値
ν_t+1 の影響を減少させることによって、一種の飛び離
れた値のリジェクションを実行する。スケーリング定数
γ_t ²に対する更新の式が付録Ｃの中で与えられている。

【００７２】付録Ｃ：γ_t の推定のための補助式次の式はγ_t ²を推定するために、付録Ａの中の式（Ａ．
１）〜（Ａ．６）および付録Ｂの中の式（Ｂ．１）〜
（Ｂ．４）に対して追加された式である。（Ａ．６）の
後に挿入する

【数２７】式（Ｃ．１）は、次のｉｆブロックと一緒に、次の状態
空間モデルによる単純なカルマン・フィルタを使ってγ
_t+1 ²を更新する。

【数２８】ここでαγの値は０．２０／（１６８＋１）である。式
（Ｃ．５）の中に現われる量を、式（１１）の中に現わ
れている量のどれとも混同してはならない。

【００７３】式（Ｃ．１）は付録Ａの中の式（Ａ．３）
に対応し、式（Ｃ．２）は式（Ａ．４）に対応し、式
（Ｃ．３）は式（Ａ．５）に対応し、そして式（Ｃ．
４）は式（Ａ．６）に対応する。ステップ（Ｃ．２）〜
（Ｃ．４）はｉｆ文の内部にあり、この文の条件が満足
された場合にのみ実行される。この条件はｆｉｒ誤差ｅ
_t+1 ^fが大き過ぎないこと、あるいは等価的に、観測され
た比率の値ν_t+1が大き過ぎないことをチェックする。
常にこのケースであってステップ（Ｃ．２）〜（Ｃ．
４）が常に実行される場合、式（Ｃ．２）の中のゲイン
Ｋ_t+1 はｔ→∞に対して１−αへ収れんし、定常状態の
式（Ｃ．３）において次のようになる。

【数２９】これは入力が（ν_t+1／σ_t+1ν）²であって、出力がγ
_t+1 ²である単極のフィルタである。

【００７４】付録Ｄ：警報の形成のための補助式付録Ａ．Ｂ、およびＣの中の更新ステップがすべて時間
ｔ＋１に対して実行された後、スケールされた警報変数
β_tが次のように計算される。

【数３０】そして警報信号は次の条件が成立する場合に与えられて
いる。

【数３１】ここで、ＴＨＲＥＳＨＯＬＤは警報のためのしきい値を
表している固定定数である。ＴＨＲＥＳＨＯＬＤは、交
換機／サービスの組合せによって異なると考えられる。

【００７５】付録Ｅ：一般的なカルマン・フィルタ監
視アルゴリズムこの申請書の本文および前の付録の中で、比率監視アル
ゴリズムについて説明し、その実装の詳細のいくつかを
示した。このセクションでは、一般のカルマン・フィル
タ監視（ＧＫＦＭ）アルゴリズムについて完全に記述す
る。このアルゴリズムは実行時に構成ファイルから１つ
またはそれ以上の時系列に対するＡＲＩ（自動縮退統合
化）モデル（正確な定義は後で示される）を読み、この
モデルを状態空間モデルへ変換し、対応しているカルマ
ン・フィルタによってこの状態空間モデルのステートを
推定し、これらのステートの１つまたはそれ以上をあら
かじめ設定されたしきい値と比較し、これらのしきい値
が超過されていた場合、警報を生成する。この申請書の
本文の中で詳細に説明されている比率のモデルによって
構成された場合、このアルゴリズムはそこで説明された
アルゴリズムと機能的に等価である。

【００７６】Ｅ．１モデルのクラスΨ このセクションでは、本発明のアルゴリズムが受け付け
るモデルのタイプを指定する。これをクラスΨと呼ぶ。
クラスΨは次の規則を使って単独チャネルの相互に相関
のないノイズ・プロセスｎⁱ 、ｉ＝１，．．．；Ｅ［ｎ
_t ⁱｎ_s ^j］＝δ_tsδ_ijのシーケンスから生成することがで
きる（シンボルδ_nmは、ｎが≠ｍの場合０、そしてｎ＝
ｍの場合１であるように定義されている）。１．ｎⁱ∈Ψ 各ｉに対してプロセスｎⁱ はΨに属している。２．ベクトル・ジョインの下での閉止Ｎチャネルのプロセスｘ＝（ｘ¹ ，．．．，ｘ^N ）およ
びＭチャネルのプロセスｗ＝（ｗ¹ ，．．．．，ｗ^M ）
が両方共Ψに属しており、それらが共有しているノイズ
・コンポーネントがない場合、すなわち、ｘⁱ ＝ｗ^j ＝
ｎ^k となるようなｉ、ｊ、ｋがない場合、次のＮ＋Ｍチ
ャネルのプロセスはΨに属する。 η＝（ｘ¹ ，．．．，ｘ^N ，ｗ¹ ，．．．，ｗ^M ）３．並列および直列の実現の下での閉止Ｎチャネルのプロセスｘ＝（ｘ¹ ，．．．，ｘ^N ）がΨ
に属している場合、そしてこのφ（ｚ）＝１−Σ_i=1ρ
α_iｚⁱ が、任意のＡＲ多項式（係数が一定で１である
多項式）である場合、コンポーネントがηⁱ＝ｘⁱ、ｉ＝
１，．．．，ＮであるＮ＋１チャネルのプロセスη＝
（η¹ ，．．．，η^N+1 ）および次の式で表されるη
^N+1 はΨに属する。

【数３２】

【００７７】ここでβ_iは任意の実数である。β_iをゲイ
ン係数と呼ぶ。例えば、次の式のシステムを考える。

【数３３】このシステムは規則１〜３を逐次適用することによって
作ることができる、Ψに属している複数チャネルηを定
義する。規則

【００７８】

【数３４】このプロセスも、Ψに属しているプロセスの場合と同様
に、方向付けられた非巡回型のグラフ、すなわち、ＤＡ
Ｇで表現されている。（Ｅ．２）のＤＡＧは、図４に示
されている。

【００７９】図４の中で現われている楕円および円は、
そのＤＡＧのノードであり、矢印は方向付けられたエッ
ジである。ノードおよび方向付けられたエッジは或る種
の追加の情報も含んでいる。例えば、ｓ² からｓ¹ へ向
かっている方向付けられたエッジはゲイン係数β₂ の値
を含んでおり、ｓ² のノードは自動縮退演算子（１−Ｂ
²⁴）^-1を含んでいる。このＤＡＧが式（Ｅ．２）からど
のように作られるかを説明する。式（Ｅ．２）を形成す
る際に規則１を適用するたびに、リーフ・ノードが生
じ、規則３を適用するたびに、非リーフ・ノードが発生
する。方向付けられたグラフのリーフはそれから出てい
る方向付けられたエッジが出ていないノードである。同
様に、ルートは方向付けられたエッジが入って来ないノ
ードである。

【００８０】リーフ・ノードに対して付加されている名
前は、それに対応しているノイズ・プロセス（規則１に
おけるｎⁱ ）の名前であり、非リーフ・ノードに対して
付加されている名前は式（Ｅ．１）におけるｚ_N+1 の役
割りを演じている単独チャネル・プロセスの名前であ
る。さらに、どの非リーフ・ノードにも式（Ｅ．１）に
おけるφの役割りを演じているＡＲ演算子の値が付加さ
れている。ただし、φがトリビアルな場合、すなわち、
φ≡１である場合を除く。その時は、この演算子を例え
ば、上記のο１およびο２のように明示的には書かな
い。非リーフ・ノードには、式（Ｅ．１）の右辺の和の
中に現われるすべての非ゼロ係数β_i に対してそこから
出ている１つの方向付けられたエッジがある。その係数
β_iに対応している方向付けられたエッジは式（Ｅ．
１）の中のｘⁱ の役割りを演じている単独チャネル・プ
ロセスに対するノードに対して出ている。

【００８１】この方向付けられたグラフのルートは変数
ο¹ およびο² であり、リーフは変数ｎ¹ 、．．．、ｎ
⁶ である。Ψの中の任意のプロセスに対応しているＤＡ
Ｇのルートを「観測変数」と呼ぶ。というのは、それら
を物理的なシステムにおける測定された、あるいは観測
されている変数と考えるからであり、そしてそのリーフ
を「ノイズ変数」と呼ぶ。というのは、それらは物理的
なシステムに対する妨害であると考えるからである。リ
ーフまたはループのいずれでもないノードを「ステート
変数」と呼ぶ。監視のアルゴリズムにおける第１のステ
ップは、観測変数の過去および現在の値の知識からステ
ート変数の値および標準偏差を推定することである。こ
の問題に対してカルマン・フィルタの理論を適用するた
めに、Ψの中のプロセスに対する状態空間表現を見つけ
ることができなければならない。これについては次のセ
クションで説明する。

【００８２】Ｅ．２ Ψの中のプロセスに対する状態空
間モデル式（付録Ａも参照）の中で、比率監視アルゴリズムに対
して使われる状態空間方程式の形式を示した。ＧＫＦＭ
（ｇｅｎｅｒａｌＫａｌｍａｎｆｉｌｔｅｒｍｏｎ
ｉｔｏｒｉｎｇ（汎用カルマン・フィルタ監視））アル
ゴリズムに対して使われる状態空間方程式は類似してい
る。

【数３５】

【００８３】式（１１）と式（Ｅ．３．ａ．ｂ）との間
の違いに注意されたい。式（Ｅ．３．ａ、ｂ）の中のシ
ステム行列Ａ_t 、Ｆ_t 、およびＣ_t は時間のインデック
スｔに依存しており、式（Ｅ．３．ａ）の中の観測ノイ
ズＮ_t ^oには、ここではゲイン行列σが乗算されており、
その観測値ｏ_t は、式（１１）の中のν_t と違ってベク
トルである可能性がある。ο_t のコンポーネントは、例
えば、図４の中のο１およびο２のような個々の監視さ
れているデータ・ストリームとなる。式（Ｅ．３．ａ、
ｂ）のタイプの状態空間モデルをクラスΨの中の或る種
のプロセスに対して構築する。それについてＧＫＦＭの
アルゴリズムを適用したい。この状態空間モデルが或る
種の望ましい特性を持つようにするために、ＤＡＧにつ
いていくつかの条件を課す必要がある。その課される条
件は次の通りである。

【００８４】１．そのＤＡＧは少なくとも１つのステ
ート変数、すなわち、ルートまたはリーフのいずれでも
ないノードを含んでいなければならない。２．そのＤＡＧは連結されていなければならない（方
向付けられていないグラフとして見られた時）。３．各ノードからノイズ変数（リーフ）へ向かう方向
付けられたエッジは、多くとも１つでなければならな
い。４．各観測ベクトル（ルート）は、ノイズ変数（リー
フ）へ向かう正確に１つの方向付けられたエッジを持っ
ていなければならない。さらに、各観測変数に対して付
加されているＡＲ演算子φは、トリビアル、すなわち、
φ≡１でなければならない。５．ノイズ変数に向かっているエッジが２つ以上あっ
てはならない。図４のＤＡＧがこれらの特性を満足することは明らかで
ある。そこで、条件１〜４を満足している１つのＤＡＧ
があると仮定し、そのノードに対する状態空間モデルを
構築することを試みる。ルートでなく、非リーフ・ノー
ドのすべて、すなわち、ステート変数ｓ^kに対応してい
るすべてのノードは次の形の式となる。

【００８５】

【数３６】ここで、β_k,i はＤＡＧの中でｓ^k から他のノードｓⁱ
へ向かっているエッジに付加されているゲイン係数であ
り、γ_k はｓ^k からｎ^k へ向かう方向付けられたエッジ
に対して付加されているゲイン係数である（あるいは、
そのような方向付けられたエッジがない場合はゼロであ
る）。

【００８６】この式は次のように書くことができる。

【数３７】ここで、Ｍはグラフの中のステート変数の個数であり、
Ｄ_k およびα_k,i 、ｉ＝０，．．．，Ｄ_k ^-1 は、それぞ
れ次の多項式の次数および非定数係数である。

【数３８】

【００８７】式（Ｅ．４）を使ってｋ番目のステート変
数の値ｓ_t ^kを更新しようとした場合、問題が生じる可能
性がある。というのは、この式の右辺の第２項の和は、
他のステート変数の既に更新されている値ｓ_t ⁱに関係
し、元の式が正しい順序になっていない限り、これらの
値は利用できないからである。このため、プロセス
（Ｅ．２）のＤＡＧ表現は式のシステムとしての表現よ
り優れている。というのは、前者によって、変数ｓ
¹ 、．．．、ｓ⁶ 、ο¹ 、ο² の依存関係を直ちに読み
取ることができるからである。式（１．２）のオリジナ
ル・リストからこの依存性の情報を抽出することは、よ
り難しく、前回の操作時に見つかった式の右辺として、
左辺が現われる式を探して繰り返しリストを操作する必
要がある。上記の式（Ｅ．２）のように、式が既に深さ
優先の順序（すなわち、ＤＡＧの中の順序と首尾一貫し
ている順序）になっている場合、この操作は一度行なわ
れるだけで済み、あるいは等価的に、ｋ＝１，．．．，
Ｍの順序で式（Ｅ．４）全体を通して実行することによ
ってステート変数の値ｓ_t ^kをすべて更新することができ
る。従って、Ψの中の任意のプロセスに対して深さ優先
の順序を見つけられることが重要である。これは次の２
つのステップで行なうことができる。

【００８８】１．先ず最初に、式のリスト全体を通し
て単独の走査を実行し、最初に遭遇した時、各変数に対
して新しいノードを割り当て、各式の右辺に現われてい
るノードを左辺に現われているすべてのノードに対して
連結することにより、進行と同時に対応しているＤＡＧ
を動的に構築する。２．次に、このＤＡＧの深さ優先のサーチを行い、そ
のノードから戻る時に各変数に対する定義式を書き出
す。

【００８９】その結果の式は、次にそれぞれの左辺の順
次更新に適した深さ優先の順序で書き出される。１およ
び２の実行のためのアルゴリズムは、方向付きのグラフ
の構築およびサーチのための標準のアルゴリズムに基づ
くことができる。そしてＧＫＦＭのアルゴリズムはこれ
らのアルゴリズムを使ってΨの任意のメンバーに対する
式を深さ優先の順序に配置する。ここで、式（Ｅ．４）
によってステート変数ｓ_t ^kの値の更新へ戻り、これらの
式が深さ優先の順序になっていることを仮定する。とい
うのは、それらが深さ優先の順序になっていない場合
に、それらをそのように配置する方法が分かったからで
ある。ｉ≧ｋでの式（Ｅ．４）の中の係数β_k,i はゼロ
になり、この式は次のように書き換えることができる。

【数３９】

【００９０】式（Ｅ．５．ａ）をｋ＝１，．．．，Ｍの
順序で通して実行することによって、すべての変数ｓ
¹ ，．．．，ｓ^M が更新される。次に、その観測された
変数を、次の式に従って更新することができる。

【数４０】ここで、変数ｓⁱ 、ｉ＝１，．．．，Ｍおよびｎ^M+l は
ゲイン係数ｃ_l,i 、ｉ＝１，．．．，Ｍおよびσ_l、ｌ＝
１，．．．，Ｌをそれぞれ持つ方向付けられたエッジ
によってο^l が連結されるノードである。式（Ｅ．３．
ｂ）は時間の添字ｔにおける一次式である。このことは
式（Ｅ．３．ｂ）がｘ_t の項でｘ_t を定義し、一方、式
（Ｅ．５．ａ）は、ｍａｘ₁ ≦_k≦_M ｛Ｄ_k ｝の順序で
あることを意味する。式（Ｅ．５．ａ）を状態空間の形
式に変換するために、元のステートｓ_t ¹、．．．、ｓ_t ^M
に対して次のように追加のステートを付加することによ
って１次式のシステムに変換しなければならない。Ｓ_k
＝Σ_j=1 ^kＤ_j がＤ_k の部分和であると定義し、次の変数
を定義する。

【数４１】

【００９１】これは１からＳ_M ＝Σ_i=1 ^MＤ_i までのすべ
てのｉに対してｘ_t ⁱを定義する。この式（Ｅ．５．ａ）
は、ｘ_t ⁱの項で次のように書き直すことができる。

【数４２】そのとき、観測ベクトルに対する式は、次のようにｘ_t ⁱ
の項で書くことができる。

【数４３】

【００９２】式（Ｅ．６．ａ）はｔの一次式であり、そ
してそれぞれの右辺は、ｋ＝１，．．．，Ｓ_M の順序で
評価されてステート・ベクトルｘ_t ＝（ｘ_t ¹，．．．，
ｘ_t ^SM ）の更新された値を作り出すことができる。それ
らは完全には式（Ｅ．３．ｂ）の形式ではない。という
のは、それらはステート・ベクトルｘ_t のコンポーネン
トにおいて再帰的であるからである。すなわち、更新さ
れたコンポーネントｘ_t ⁱのいくつかに対する式はｊ＜ｉ
での更新されたコンポーネントｘ_t ^jを含むからである。
しかし、ｘ_t ^jはこれらの式の右辺から簡単に消去され
て、式（Ｅ．３．）の形式のシステムを作ることができ
る。これが行なわれると、このシステムに対する行列Ａ
_t 、Ｆ_t 、およびＣ_t は、時間インデックスｔには依存
しなくなる。

【００９３】式（Ｅ．６．ａ）は、ｘ_t-1 →ｘ_t の更新
演算が結局は任意のｉ≠Ｓ_j 、ｊ＝１，．．．，Ｍに対
するベクトルｘ_t-1 のコンポーネントｘ_t-1 ⁱの単純な置
き換え、あるいは等価的には、これらのコンポーネント
のインデックスの並べ換えになることを示している。解
体された各インターバル［Ｓ_k-1 ＋１，Ｓ_k ］に及ぼす
この並べ換えの効果は、このインターバルの中でインデ
ックスを左へ１単位だけ巡回的に左へシフトすることで
ある。これを並べ換えｌ（左シフトに対して）と呼ぶ。
ｌに対する明示的な式は次のようにになる。

【数４４】

【００９４】任意の整数ｍに対して、ｌ^m （ｌのｍ番目
の繰返し）は、各インターバル［Ｓ_k-1 ＋１，Ｓ_k ］に
おけるｍ単位だけの左巡回シフトである。これは負数ｍ
だけ左に巡回シフトすることは｜ｍ｜単位だけ右に巡回
シフトすることであると解釈する場合、負の整数に対し
ても成立する。従って、ｌ^-1をｒ（右シフトに対して）
と定義する。式（Ｅ．６．ａ）の中のステート・ベクト
ルｘのコンポーネントのシフトは、時間の掛かる演算で
あり、それは行列Ａ_t 、Ｆ_t 、およびＣ_t が時間のイン
デックスｔに依存することが許される異なる状態空間モ
デルを選定することによって消去することができる。こ
の目的で、次の関係によって新しいステート・ベクトル
ξ＝（ξ¹ ，．．．，ξ^SM）のコンポーネントを定義す
る。

【数４５】

【００９５】次に式（Ｅ．６．ａ）はξ_t の項で書かれ
て次のようになる。すべてのｉ≠ｒ^t（Ｓ_k）、ｋ＝
１，．．．，Ｍに対して、

【数４６】これは新しいステート・ベクトルξのコンポーネントξ
ｉ、ｉ≠ｒ^t（Ｓ_k）がエポックｔ−１からエポックｔま
で不変であることを示している。これは論理配列ξ_t-1
およびξ_t を保持するために使われる同じ物理的配列に
おけるコンピュータの実装において、式（Ｅ．７）を全
く消去できることを意味している。残りの式（Ｅ．６．
ａ）は、ここでξ_t の項で書かれて次のようになる。

【数４７】

【００９６】ｋ＝１，．．．，Ｍに対する式（Ｅ．８）
は、ξ_t-1 からξ_t までのステート・ベクトルを更新す
るための式の完全なセットを含んでいる。その式のセッ
ト（Ｅ．８）を式（Ｅ．３．ｂ）の形の行列式として書
き直すことができる。固定された各ｎに対して、そして
ｋ＝１，．．．，Ｍに対して、１つのＳ_M ×Ｓ_M 行列_k
Ａ＝［_k Ａ_ij］および１つのＳ_M ×Ｍ行列_k Ｆ＝［_k Ｆ
_ij］を次の式によって定義する。

【数４８】そして、ｋ＝０，．．．，Ｍに対して、ベクトル_k ξ≡
（_k ξ¹，．．．，_kξ^SM）を次の式によって定義する。

【数４９】

【００９７】そのとき、式（Ｅ．８）がｋ＝
１，．．．，Ｍに対して次のように書けることは明らか
である。

【数５０】ここで、Ｎ_t ＝（ｎ_t ¹，．．．，ｎ_t ^M）は、コンポーネ
ントが式（Ｅ．８）の中に現われているノイズ変数に等
しいベクトルである。Ｎ_t は、プロセスまたはステート
・ノイズと呼ばれる。 ₀ξ＝ξ_t-1 、そして _Mξ＝ξ_t
なので、次の式が得られる。

【数５１】

【００９８】これはシステム行列が次の式によって与え
られるステート変数に対する（Ｅ．３．ｂ）の形式の１
つの式を与える。

【数５２】観測された変数ο^l に対する式をξおよびｎⁱ の項で次
のように書くことができる。

【数５３】

【００９９】次の式によってＬ×Ｍの行列Ｃ_t ＝
［Ｃ_ij］を定義した場合、

【数５４】そしてベクトルο_t を次の式によって定義し、

【数５５】そしてσが対角行列であるとすると、

【数５６】観測ベクトルο_t に対する式は、次のように行列の形式
で書くことができる。

【数５７】ここでＮ_t ^o＝（ｎ_t ^M+1，．．．，ｎ_t ^M+l）は、そのコン
ポーネントが式（Ｅ．１３）に現われているノイズ変数
に等しいベクトルである。この式は式（Ｅ．３．ａ）の
形の式である。要約すると、このプロセス（Ψに属する
と仮定されている）に対する状態空間表現は次の式で表
されることが分かる。

【数５８】ここで、システム行列Ａ_t 、Ｆ_t 、Ｃ_t 、およびσは式
（Ｅ．１１）、（Ｅ．１２）、（Ｅ．１４）、（Ｅ．１
５）によって与えられ、そこではステート・ベクトルに
対するξの代わりに元の表記に戻している。というの
は、システム（Ｅ．１６）はこれ以降で説明する唯一の
ベクトルであり、式（Ｅ．１６．ｂ）の中のｘと式
（Ｅ．６．ａ）の中のｘとの間に混同があってはならな
いからである。

【０１００】Ｅ．３カルマン・フィルタの更新式付録Ａにおいて、単独の観測のケース（Ｌ＝１）に対す
るカルマン・フィルタの更新式を示した。ここで、複数
の観測のケース（Ｌ＞１）に対するカルマン・フィルタ
の更新式を示し、それらを２つの部分、すなわち、時間
更新部分および測定値更新部分に分割する。これはＧＫ
ＦＭアルゴリズムの中で実装されている方法に対応す
る。最適のステート推定値ｘ＾_t およびその誤差共分散
行列Ｅ［（ｘ＾_t −ｘ_t ）（ｘ＾_t −ｘ_t ）^T ］に対す
る複数の観測値のケースに対する従来のカルマン・フィ
ルタの式は次の式によって与えられる。時間更新式

【数５９】測定値更新式

【数６０】

【０１０１】ｘ￣_t+1 はο_t 、ο_t-1 、．．．が与えら
れた場合のｘ_t+1 の線形最小二乗法による推定値であ
り、そして、Ｈ_t+1 ＝Ｅ［（ｘ￣_t −ｘ_t ）（ｘ￣_t −
ｘ_t ） ^T ］は、その誤差共分散行列である。Ｋ_tはカル
マン・ゲイン行列である。これは観測された１ステップ
先行の予測誤差ο_t −Ｃ_t ｘ￣_t に掛かり、１ステップ
先行の予測値に対して加えられなければならない補正項
を与え、現在の推定値ｘ＾_t が得られる。

【０１０２】ＧＫＦＭのアルゴリズムにおいて、Ｐ_t そ
のものではなく、Ｐ_t の平方根（すなわち、Ｐ_t ＝ＬＬ
^T であるような行列Ｌ）を更新し、ゲイン行列Ｋ_t をこ
の平方根の項で表す。（対称行列Ｐは、それが平方根を
持っている場合、そしてその場合だけ非負有界であ
る。）Ｐ_t の平方根を更新するような、ここで考えて
いるアルゴリズムは、平方根アルゴリズムと呼ばれ、そ
の詳細はビールマンＧ．Ｊ．の「ディスクリート順次推
定のための因数分解法（Ｆａｃｔｏｒｉｚａｔｉｏｎ
ＭｅｔｈｏｄｓｆｏｒＤｉｓｃｒｅｔｅｓｅｑｕ
ｅｎｔｉａｌＥｓｔｉｍａｔｉｏｎ）」Ａｃａｄｅｍ
ｉｃＰｒｅｓｓ，Ｉｎｃ．ＮｅｗＹｏｒｋ，１９７
７年の中で説明されている。この平方根に対して利用で
きる多くの選択肢のうち、所要メモリが最小であるもの
はＰ_t の方形下側三角形（またはチョレスキー）ファク
タＬである。というのは、メインの対角要素の上にある
Ｌの要素は０であり、記録される必要がないからである
（要素がＬ_ijである行列Ｌは、すべてのｉ＜ｊに対して
Ｌ_ij＝０である場合、下側三角形である。この定義は正
方行列以外に、非正方行列を含む）Ψの中のプロセスに
対する超対角要素の数は１／２Ｓ_M （Ｓ_M −１）であ
り、そして周期的なＡＲＩプロセスに対してはＳ_M＝Σ
_i=1 ^MＤ_i は大きくなり得るので、これは重要な考慮事項
である。この理由で、ＧＫＦＭのアルゴリズムはＰ_t お
よびＨ_t のチョレスキーの平方根を更新する。次に、そ
の方法を示す。

【０１０３】Ｅ．４ＧＫＦＭ時間更新アルゴリズムＥ．４．１説明ＧＫＦＭのアルゴリズムは、Ｈ_t+1の下側三角形チョレ
スキー・ファクタＬ〜を形成する。

【数６１】これは、最初にＨ_t+1 を形成する必要なしに、Ｐ_t の下
側三角形のチョレスキー・ファクタＬを直接形成する。
Ｈ_t+1 ＝Ａ_t Ｐ_t Ａ_t ^T＋Ｆ_t Ｆ_t ^Tであり、そして、Ｐ_t
＝ＬＬ^T であることから、次の式を書くことができ
る。

【数６２】これは、行列［Ａ_t ＬＦ_t ］が下側三角形でもなく正
方形でもない（その次元はＳ_M ×（Ｓ_M ＋Ｌ）である）
ことを除いて式（Ｅ．１９）の形をしている。式（Ｅ．
２０）をＬ〜が正方形の下側三角形になるように式
（Ｅ．１９）の形式に変換する必要がある。次の式で表
されるような（Ｓ_M ＋Ｌ）×（Ｓ_M ＋Ｌ）の直交行列Ｏ
（ＯＯ^T ＋Ｉ）を見つけることができると仮定する。

【数６３】ここで、Ｌ〜は正方形の下側三角形であり、０はＳ_M ×
Ｌのゼロ行列に等しい。そうすると、次の式が得られ
る。

【数６４】そして、Ｌ〜はＨ_t+1 のチョレスキー・ファクタとな
る。そこで、Ｈ_t+1のチョレスキー・ファクタを見つけ
る問題を、式（Ｅ．２１）が成立するような直交行列Ｏ
を見つける問題に帰せしめたことになる。これを行なう
ためのアルゴリズムの概要を示す前に次の２つの定義に
ついて説明する必要がある。

【０１０４】ギブンズのローテーション−ギブンズのロ
ーテーションＧ（ｊ，ｋ，θ）、ｊ＜ｋは、ｊ番目およ
びｋ番目の座標ベクトルによって張られた部分空間上で
θラジアンだけ平面的に回転するように働き、そして残
りの座標ベクトルによって張られた部分空間上でのアイ
デンティティとして働く行列であると定義される。明ら
かにＧ（ｊ，ｋ，θ）は、行列式が１である正方形の直
交行列である。ヘッセンバーグの行列−（下側の）ヘッ
センバーグ行列Ｍ＝［Ｍ_ij］は、そのメインの対角ｉ＝
ｊの上側にある非ゼロ要素はｊ＝ｉ＋１の要素だけであ
るような、すなわち、ｊ＞ｉ＋１のすべてのｊに対して
Ｍ_ij＝０であるような行列である。

【０１０５】Ｐ_t のチョレスキー・ファクタＬが与えら
れて、Ｈ_t+1 のチョレスキー・ファクタＬ〜を見つける
ためのアルゴリズムは次のようになる。₀Ｌ＝Ｌに設定
する。ｋ＝１，．．．，Ｍに対して、１． _kＡ・_k-1Ｌを作る。ここで_kＡはコンポーネント
が式（Ｅ．９）によって与えられている行列である。２．ギブンズの行列Ｓ_k −ｒ^t （Ｓ_k ）によって後か
ら乗算することによって、_kＡ・_k-1Ｌを下側のヘッセン
バーグの形式_kＬ_Hに変換する。３．ギブンズの行列Ｓ_k −ｒ^t （Ｓ_k ）−１を後から
乗算することによって、_kＬ_Hを下側の三角形状_kＬ_Tに変
換する。４．ギブンズの行列Ｓ_M −ｒ_t （Ｓ_k ）＋１によって
後から乗算することによって、［_kＬ_T _kＦ］を下側の
三角形状［_kＬ０］に変換する。ここで _kＦはコンポ
ーネントが式（Ｅ．１０）によって与えられている行列
である。Ｌ〜＝ _MＬに設定する。行列Ｌ〜はこれでＨ
_t+1 の下側三角形（チョレスキー）ファクタになる。

【０１０６】このアルゴリズムは次の事実によって有効
に使える。Ｉ．直交行列の積は直交行列である。特に、ステップ
２および３におけるギブンズの行列の積は直交行列であ
り、そしてステップ４におけるギブンズの行列の積もそ
うである。 II．Ｏ１がステップ２および３におけるギブンズの行
列の積であり、次の式が成立し、

【数６５】ここで_kＬ^Tは正方形の下側三角形行列であり、そしてＯ
₂がステップ４におけるギブンズの行列の積であって、
次の式が成立する。

【数６６】ここで_kＬが正方形の下側の三角形行列である場合、次
の式は式（Ｅ．２１）の形式の式である。

【数６７】というのは、［Ｏ₁］・Ｏ₂が直交行列だからである。Ｉ

【０１０７】III. 次の式によってｋ＝０，．．．，Ｍ
に対する行列_kＰを誘導的に定義する。

【数６８】そうすると、ａ．ステップ４の中に現われる_kＬは_kＰのチョレスキ
ー・ファクタである。すなわち、_kＰ＝_kＬ・_kＬ^Tであ
る。ｂ． _MＰ＝Ｈ_t+1 III.ａは式（Ｅ．２２）が式（Ｅ．１７．ｂ）および
（Ｅ．２０）から得られるのと同じ方法で、式（Ｅ．２
４）および（Ｅ．２３）から得られる。III.ｂはｊ＜ｋ
に対して_kＦ（_jＦ）^T ＝０であるという事実（セクショ
ン３の条件５から得られる）を使って、ｋについての誘
導によって式（Ｅ．２４）から得られる。この時間更新
アルゴリズムの完全な疑似コードのリスティングが次の
サブセクションに与えられている。

【０１０８】Ｅ．４．２ＧＫＦＭの時間更新部分に対
する疑似コードのリスティング。図５はＧＫＦＭの時間更新部分に対する疑似コードのリ
スティングである。図５において、非ゼロの自動縮退係
数α_k,i がα_k,0 だけであるｋの値に対してステップ２
および３がスキップされていることに注意されたい。こ
れらのｋの値に対して、行列_kＡ・_k-1Ｌは既に下側の三
角形であり（行列 _kＡがそうであるので）、従ってステ
ップ２および３は不要である。これは図４に示されてい
るプロセスに対して、ｋ＝１，．．．，６のすべてのｋ
に対して成立する。また、ステップ４はγ_k ＝０となる
ｋの値に対してスキップされている。ｋのこれらの値に
対して、行列Ｆ_t ^kは０である。これは図４の中のｋ＝３
および４に対して成立する。

【０１０９】Ｅ．４．３ＧＫＦＭ測定更新アルゴリズ
ム前のセクションにおいて、カルマン・フィルタの時間更
新方程式（Ｅ．１７．ａ）および（Ｅ．１７．ｂ）を実
装するための平方根タイプのアルゴリズムを提示した。
このアルゴリズムは、ステートの推定値ｘ＾_t を更新す
るための完全なアルゴリズムを作り出すために、測定更
新式（Ｅ．１８）に対する平方根アルゴリズムと組み合
わせられなければならない。この測定更新アルゴリズム
は、右辺に与えられているＨ_t のチョレスキー・ファク
タが与えられて、式（Ｅ．１８）の左辺に現われている
Ｐ_t のチョレスキー・ファクタを見つけ、このＨ_t の項
でＫ_t を計算しなければならない。このケースＬ＝１
（単独観測ケース）において、これを行なうための標準
のアルゴリズム、すなわち、カールソンのアルゴリズム
がある。それはビールマン、Ｇ．Ｊ．の「ディスクリー
ト順次推定に対する因数分解の方法（Ｆａｃｔｏｒｉｚ
ａｔｉｏｎＭｅｔｈｏｄｓｆｏｒＤｉｓｃｒｅｔ
ｅＳｅｑｕｅｎｔｉａｌＥｓｔｉｍａｔｉｏ
ｎ）」、（ＡｃａｄｅｍｉｃＰｒｅｓｓ，Ｉｎｃ．，
ＮｅｗＹｏｒｋ，１９７７）（参照によってここに組
み込まれている）のページ８１に記述されており、ＦＯ
ＲＴＲＡＮのリスティングが同じ本のページ１０２に示
されている。

【０１１０】Ｌ＞１のケースをＬ＝１に縮小するための
標準のアルゴリズムも存在し、それは順次アルゴリズム
と呼ばれる。このアルゴリズムは、チュイ、Ｃ．Ｋ．お
よびチェン、Ｇ．の「リアルタイム・アプリケーション
でのカルマン・フィルタ（ＫａｌｍａｎＦｉｌｔｅｒ
ｉｎｇｗｉｔｈＲｅａｌＴｍｅＡｐｐｌｉｃａ
ｔｉｏｎｓ）」（ＳｐｒｉｎｇｅｒＶｅｒｌａｇ，Ｎ
ｅｗＹｏｒｋ，１９９１）（参照によってここに組み
込まれている）のセクション７．１、ページ９７〜１０
２に示されている。ＧＫＦＭアルゴリズムの測定更新部
分は、これらの２つのアルゴリズムを組み合わせること
によって得られる。この組み合わせられたアルゴリズム
は、チュイ他の本の中のセクション７．３、ページ１０
５〜１０６に示されている。ただし、チュイ他によって
与えられているステップ（ｉｉ）（ｂ）に対するカール
ソンのアルゴリズムがＧＫＦＭアルゴリズムによって置
き換えられている。

【０１１１】時間更新ステップにおいて、Ｍ個の式
（Ｅ．２４）に対して共分散更新式（Ｅ．１７．ｂ）を
変形し、その各式を効率的に更新できたことを思い出さ
れたい。同様に、この順次アルゴリズムは複数観測（Ｌ
＞１）の測定更新式（Ｅ．１８）をＬ個の単独観測測定
更新式に変形する。さらに正確には、Ｃ_t ＝［ｃ¹ ・・
・ｃ^L ］^T とし、次の方式に従って、ｌ＝１，．．．，
Ｌに対してベクトルのシーケンス_lｋ，_lｘ＾および行列
_l Ｈを定義する。

【数６９】この時、順次アルゴリズムは次の事実に基づいている。

【数７０】

【０１１２】これはチュイにおいて定理７．１として証
明されている。量σ_l、ｌ＝１，．．．，Ｌは、行列σの
対角要素である。式（Ｅ．２５．ｌ）は式（Ｅ．１８）
の単独観測（Ｌ＝１）のケースに対する形式と同じであ
り、従って、カールソンのアルゴリズムは _lＨのチョレ
スキー・ファクタ以外に、式（Ｅ．２５．ｌ）の左辺に
ある量_lｋおよび_lｘ＾を計算するために適用できる。こ
れをｌ＝１，．．．，ＬのすべてのＬに対して行なう
と、式（Ｅ．２６）および（Ｅ．２７）による結果はス
テートの推定値ｘ＾_t の更新された値以外のＨ_t のチョ
レスキー・ファクタが得られることになる。順次アルゴ
リズムが効果的にＨ_t およびｘ＾_t を計算できることの
他に、このアルゴリズムの重要事項は、このアルゴリズ
ムによって飛び離れた値のリジェクションおよび欠落デ
ータの処理が各チャネルに別々に適用できることであ
る。このことは次のサブセクションにおいて明らかとな
る。そこではＧＫＦＭアルゴリズムの疑似コード・リス
ティングの残りの部分を示す。

【０１１３】Ｅ．４．４完全なＧＫＦＭアルゴリズムここで、上で開発したアルゴリズムを完全な汎用の監視
アルゴリズムにする。このアルゴリズムのステップは次
のようになる。エポックｔにおいて：１． γ_temp＝０Ｌγ＝０と設定する２．ＧＫＦＭアルゴリズムの時間更新部分を使ってＰ
_tのチョレスキー・ファクタＬをＨ_t+1のチョレスキー・
ファクタＬ〜に変換する。ｌ＝１，．．．，Ｌに対し
て、３．ｌ欠落データの処理−ｏ_t ^lが欠落している場合、
ループのインデックスｌをｌ＋１にインクリメントし、
このループを継続する（ステップ３．ｌ＋１から始め
て）４．ｌチャネルｌについて飛び離れた値のリジェクシ
ョンを実行する−下記参照。５．ｌカールソンのアルゴリズムを使って、ステート
の推定値_lｘ＾および式（Ｅ．２５．ｌ）における _lＨ
のチョレスキー・ファクタの更新値を計算する。６． γ_tをγ_t+1に更新する。７．警報変数を作る。

【０１１４】上記で、ステップ２および５．ｌについて
説明した。ステップ３．ｌは自明である。次に（４．
ｌ、６、および７）について詳しく説明する。ステップ
４．ｌは付録Ｂに示されている単独観測の場合と同様で
ある。以下の式（Ｂ’．１〜４）は、付録Ｂの中の式
（Ｂ’．１〜４）に対応している。ステップ４．ｌ式（Ｂ’１〜Ｂ’４）

【数７１】

【０１１５】同様に、ステップ６は、付録Ｃに示されて
いる単独の観測の場合に類似している。以下の式
（Ｃ’．１−４）は、付録Ｃの中の（Ｃ．１−４）に対
応している。ステップ６式（Ｃ’１）〜（Ｃ’４）

【数７２】

【０１１６】ここで、ステップ７について説明する。付
録Ｄにおいて、比率監視アルゴリズムに対する警報変数
β_t+1 は、ステート推定値ベクトルｘ＾_t+1 の１つのコ
ンポーネントｘ＾_t+1 ²⁵ の正規化されたバージョンであ
ることを知った。同様に、ＧＫＦＭのアルゴリズムの中
で、ステート推定値ベクトルｘ＾_t+1 のコンポーネント
の線形の組合せの正規化されたバージョンとして、１つ
またはそれ以上の警報変数を定義することができる。そ
のような警報変数β_t の一般的な定義は次の通りであ
る。

【数７３】ここで、ａ_iは任意の実数であってよい。警報（β_t に
対応しているタイプの）は、次の時に発行される。｜β_t｜＞ＴＨＲＥＳＨＯＬＤこれは付録Ｄの中の式と同じである。ここでＴＨＲＥＳ
ＨＯＬＤは警報のためのしきい値を表している固定の定
数である。勿論、異なる値のＴＨＲＥＳＨＯＬＤを選定
して警報の異なるレベルの厳格さを表すことができる。

【０１１７】付録Ｆ：Ｒ７．０ＣＤＲＭＫＦ監視
アルゴリズムに対するディレクトリ構造図１は、ＣＤＲＭカルマン・フィルタ監視アルゴリズム
に関連しているすべてのファイル（脚注１）を含んでい
るディレクトリ・ツリーを示している。ディレクト
リ．．／ｂｉｎは、実行可能なプログラムＫａｌｍａｎ
を含んでいる。ディレクトリ．．／ｓｒｃは、ソース・
コードのファイルを含んでいる。これらのファイルのう
ち２つを除くすべてのファイルはＣで書かれており、拡
張子＊．ｃが付いている。残りの２つのファイル、すな
わち、ｒｅａｄｌ．ｌおよびｒｅａｄｙ．ｙはｌｅｘお
よびｙａｃｃでそれぞれ書かれている。また、ディレク
トリ．．／ｓｒｃは、Ｋａｌｍａｎをコンパイルするた
めのメイクファイルであるｍａｋｅｆｉｌｅをも含んで
いる。ディレクトリ．．／ｈｄｒは、ヘッダー・ファイ
ルを含んでおり、それらは拡張子が＊．ｈであり、コン
パイル時にソース・コード・ファイルの中に含められ
る。最後に、ディレクトリ．．／ｐａｒｍは１つまたは
それ以上のサブディレクトリを含んでいる。各サブディ
レクトリは、この監視アルゴリズムに対する完全な一組
のコンフィギュレーション・ファイルを含んでいる。こ
のアルゴリズムは、Ｋａｌｍａｎが呼び出される時に、
−ｐのコマンド行オプションに対する引き数として、こ
れらのサブディレクトリの１つのパス名を指定すること
によって構成される。．．／ｐａｒｍの異なるサブディ
レクトリは、異なるタイプのデータを監視するために必
要なコンフィギュレーション・ファイルを含んでいる。
例えば、．．／ｐａｒｍは、サブディレクトリ．．／ｐ
ａｒｍ／ｒａｔｉｏｓおよび．．／ｐａｒｍ／ｖｏｌｕ
ｍｅｓを含むことができ、それぞれが比率およびボリュ
ームの時系列を監視するために必要なコンフィギュレー
ション・ファイルを含んでいる。３つのコンフィギュレ
ーション・ファイルが、この監視アルゴリズムを実行す
るために必要であり、それらは−ｐコマンド行オプショ
ンに対する引き数として現われる．．／ｐａｒｍのどれ
かのサブディレクトリの中に存在していなければならな
い。

【０１１８】これらのファイルおよびその目的は次の通
りである。ｋａｌｍａｎ．ｄ−入力および出力のデータ・ファイル
を指定する。ｋａｌｍａｎ．ｍ−カルマン・フィルタに対するモデル
を指定する。ｋａｌｍａｎ．ｐ−或る種の他のパラメータの値を指定
する。これらの３つのファイルは、ユーザ編集によるＡＳＣＩ
Ｉファイルである。さらに、第４のファイル、ｋａｌｍ
ａｎ．ｉがあり得る。このファイルが．．／ｐａｒｍの
１つのサブディレクトリの中に存在する時、それは或る
１つの時点における監視アルゴリズムの完全なステート
のダンプを含んでいるバイナリ・ファイルである。この
ファイルｋａｌｍａｎ．ｉは、ｋａｌｍａｎがコマンド
行オプション−ｘで実行されている時に生成され、それ
以降の実行時に、そのアルゴリズムを初期化するために
使うことができる。

【０１１９】また、図１は、監視コードの中で定義され
ている外部シンボルもリストしている。．．／ｈｄｒの
中にある各ヘッダー・ファイルの下に、そのファイルの
中で定義されている外部変数がリストされており、．．
／ｓｒｃの中の各ソース・ファイルの下には、そのファ
イルの中で定義されている関数およびサブルーチンがリ
ストされている。（脚注２）脚注脚注１入力および出力のデータ・ファイルは例外。こ
れらのファイルに対するポインタはコンフィギュレーシ
ョン・ファイル．．／ｐａｒｍ／＊／ｋａｌｍａｎ．ｄ
の中で与えられている。脚注２ｍａｋｅでｋａｌｍａｎをコンパイルすると、
ｒｅａｄｌ．ｌおよびｒｅａｄｙ．ｙにおいてｌｅｘお
よびｙａｃｃが実行されて、２つのソース・ファイル、
すなわち、ｒｅａｄｌ．ｃおよびｒｅａｄｙ．ｃが作ら
れる。これらのファイルは図１の中にはリストされてい
ないｙｙｌｅｘ（）およびｙｙｐａｒｓｅ（）などの追
加のファイルを含んでいる。ｍａｋｅはコンパイルした
後、ｒｅａｄｌ．ｃおよびｒｅａｄｙ．ｃを取り除く。

【０１２０】

【数７４】

【図面の簡単な説明】

【図１】４週間の期間にわたって記録されたＡＴ＆Ｔの
４ＥＳＳ^TM交換機に関連している比率時系列の一例であ
る。

【図２】本発明の１実施例に従って比率時系列を監視す
るために使われるシステムのブロック図である。

【図３】図２の実施例に従って各種のステートを生成す
るために使われる数学的モデルのブロック図である。

【図４】好適な実施例に関連して使われる方向付きの非
周期的グラフである。

【図５】好適な実施例の一般的なカルマン・フィルタ監
視アルゴリズムの部分を時間的に更新する部分に対する
疑似コードのリスティングである。

【図６】本発明の好適な実施例の監視方法に関連したデ
ィレクトリのツリーである。

Claims

【特許請求の範囲】

【請求項１】複数の顧客に対してサービスを提供する
電話ネットワークにおける記録の誤りを検出するための
方法であって、前記顧客に対するサービスの提供に関連している前記電
話ネットワークによって出力される情報を記録し、前記
情報は周期的なパターンを持っているような、記録のス
テップと、前記周期的パターンの中であらかじめ決められた期間よ
り長い間持続する妨害を識別するステップと、前記識別された妨害がしきい値を超えると、警報を生成
するステップとを含む方法。
【請求項２】複数の顧客に対してサービスを提供する
電話ネットワークにおける記録の誤りを検出するための
方法であって、前記顧客に対するサービスの提供に関連している第１の
期間において、前記電話ネットワークによって出力され
る情報を記録し、前記情報は周期的なパターンを持って
いるような、記録のステップと、モデルを使って、前記第１の期間の間に記録された前記
の情報を利用して、まだ発生していない第２の期間にお
いて前記ネットワークによって出力されることになる情
報の予測値を生成し、前記モデルは前記周期的パターン
における妨害を表している複数のステートを含み、前記
各ステートは持続時間が異なる妨害を表しているよう
な、モデルを使うステップと、前記第２の期間の間に前記顧客に対する前記サービスの
提供に関連している前記電話ネットワークによって出力
される情報を得るステップと、前記第２の期間の間に得られた情報と、前記予測値とを
比較して偏差値を生成する比較のステップと、前記偏差を使って前記の各ステートを更新するステップ
と、前記ステートのうちの選択された１つがしきい値を超え
た時に警報を生成するステップとを含む方法。
【請求項３】前記ステートのうちの前記選択された１
つが前記他のステートのどれよりも持続時間が長いこと
を特徴とする、請求項２に記載の方法。
【請求項４】前記の記録された情報が所定の期間に渡
って前記ネットワークの中で切り離された呼出しに対す
る前記ネットワークの中の接続された呼出しの比率であ
ることを特徴とする、請求項２に記載の方法。
【請求項５】前記記録された情報が、所定の期間にお
いて前記ネットワークの中で行なわれた呼出しのボリュ
ームであることを特徴とする、請求項２に記載の方法。
【請求項６】予測値を生成させるためのモデルを使う
前記のステップが、前記第１の期間において記録された
情報を利用するステップに加えて、前記第１の期間より
前に記録されていた情報を利用するステップを含んでい
ることを特徴とする、請求項２に記載の方法。
【請求項７】前記第２の期間の間に記録された前記情
報および前記第２の期間より前に記録された情報を利用
して、まだ発生していない第３の期間に対する前記ネッ
トワークによって出力されることになる情報の期待され
る内容を予測するステップと、前記第３の期間において前記顧客に対する前記サービス
の提供に関連している、前記電話ネットワークによって
出力される情報を記録するステップと、前記第３の期間に関連している予測の前記ステップの結
果と、前記第３の期間の間に記録された情報とを比較し
て、第２の偏差を生成するステップとをさらに含む、請
求項２に記載の方法。
【請求項８】前記各ステートを更新するために前記第
２の偏差を使うステップと、前記ステートの選択された１つがしきい値を超えた時に
警報を生成するステップとをさらに含む、請求項７に記
載の方法。
【請求項９】前記ステップが前記ネットワークによっ
て提供されている複数のサービスに対して実行されるこ
とを特徴とする、請求項２に記載の方法。
【請求項１０】前記の各ステップが絶えず繰り返され
て、前記ネットワークによって最近出力されたデータと
前記の最近出力されたデータに対する予測データとの比
較に対応している偏差値を生成し、そして前記ステート
が繰り返し更新され、選択された更新済みのステートが
しきい値を超えた場合に、警報が生成されることを特徴
とする、請求項２に記載の方法。
【請求項１１】前記モデルを生成するステップにおい
て、前記モデルがカルマン・フィルタを使って生成され
るステップをさらに含む、請求項２に記載の方法。
【請求項１２】複数の顧客に対してサービスを提供す
る電話ネットワークにおける記録の誤りを検出するため
の方法であって、前記顧客に対するサービスに関連している前記電話ネッ
トワークによって出力される情報を記録し、前記情報が
周期的なパターンを持っているような、記録のステップ
と、前記記録された情報に基づいて１つのモデルを生成し、
前記モデルは第１の期間の前記周期的パターンにおける
妨害を表している第１のステートと、第２の期間の前記
周期的パターンにおける妨害を表している第２のステー
トを持ち、前記第２の期間が前記第１の期間より長いよ
うな、モデル生成のステップと、前記第２の期間における妨害がしきい値を超えた時に警
報を生成するステップとを含む方法。
【請求項１３】前記モデルを使って、まだ発生してい
ない期間に対する前記電話ネットワークによって出力さ
れることになる情報の予測値を生成するステップと、前記期間の間に前記電話ネットワークによって出力され
る情報を記録するステップと、前記予測のステップの結果と、前記期間の間に記録され
た情報とを比較して偏差値を生成するステップとをさら
に含む、請求項１２に記載の方法。
【請求項１４】前記各ステートを更新するために、前
記偏差値を使うステップをさらに含む、請求項１３に記
載の方法。
【請求項１５】複数の顧客に対してサービスを提供す
る電話ネットワークにおける記録の誤りを検出するため
の方法であって、前記顧客に対するサービスの提供に関連している前記電
話ネットワークによって出力される情報を記録し、前記
情報が周期的パターンを持っているような記録のステッ
プと、第１の期間の前記周期的パターンにおける妨害を識別す
るステップと、第２の期間の前記周期的パターンにおける妨害を識別
し、前記第２の期間が前記第１の期間より長いような、
識別のステップと、前記第２の期間の前記妨害がしきい値を超えた時、警報
を生成するステップとを含む方法。
【請求項１６】複数の顧客に対してサービスを提供す
る電話ネットワークにおける記録の誤りを検出するため
の方法であって、（ａ）第１の期間において前記顧客に対するサービスの
提供に関連している前記電話ネットワークによって出力
される情報を記録するステップと、（ｂ）前記第１の期間において記録された前記情報に基
づいてモデルを生成し、前記モデルが第１の持続時間の
前記情報における妨害を表している第１のステートおよ
び、第２の持続期間の前記情報における妨害を表してい
る第２のステート・データを含んでいて、前記第２の持
続期間が前記第１の持続期間より長いような、モデル生
成のステップと、（ｃ）将来の期間に対して前記ネットワークによって出
力されることになる情報の予測値を生成するために前記
モデルを使用するステップと、（ｄ）前記顧客に対する前記サービスの提供に関連して
いる前記電話ネットワークによって出力される情報を、
前記将来の期間の発生時に記録するステップと、（ｅ）ステップ（ｃ）において生成された予測値と、ス
テップ（ｄ）において記録された情報とを比較するステ
ップと、以前のステップにおいて生成された前記偏差値に基づい
て前記各ステートを更新するステップと、（ｇ）前記ステートのうちの選択された１つがしきい値
を超えた時、警報を生成するステップと、（ｈ）ステップ（ｃ）〜（ｇ）を繰り返し実行して、前
記ステートのうちの前記選択された１つが前記しきい値
を超えた時は常に絶えず警報を生成するステップとを含
む方法。
【請求項１７】前記ステートのうちの前記選択された
１つが前記第２のステートであることを特徴とする、請
求項１６に記載の方法。
【請求項１８】複数の顧客に対してサービスを提供す
る電話ネットワークにおける記録の誤りを検出するシス
テムであって、前記顧客に対するサービスの提供に関連している、前記
電話ネットワークによって出力される情報を記録し、前
記情報が周期的パターンを持っているような、情報の記
録装置と、所定の期間より長い間持続する前記周期的パターンにお
ける妨害を識別するための装置と、前記識別された妨害がしきい値を超えた時に警報を生成
するための警報発生器とを含むシステム。
【請求項１９】複数の顧客に対してサービスを提供す
る電話ネットワークにおける記録の誤りを検出するため
のシステムであって、前記顧客に対するサービスの提供に関連している第１の
期間において前記電話ネットワークによって出力された
情報を記録するための情報記録装置と、前記第１の期間の間に記録された情報を利用して、まだ
起こっていない第２の期間における前記ネットワークに
よって出力されることになる情報の予測値を生成し、前
記モデルは前記周期的パターンにおける妨害を表してい
る複数のステートを含み、前記各ステートが、持続期間
の異なる妨害を表しているモデルと、前記第２の期間の間に、前記顧客に対する前記サービス
の提供に関連している前記電話ネットワークによって出
力される情報を得るための装置と、前記予測値を前記第２の期間の間に得られた情報と比較
して、偏差値を生成し、前記偏差値が前記の各ステート
を更新するために使われるようにするための比較器と、前記ステートの選択された１つがしきい値を超えた時に
警報を生成するための警報発生器とを含んでいる装置。
【請求項２０】前記ステートのうちの前記選択された
１つが、前記他のステートのどれよりも長く持続するス
テートであることを特徴とする、請求項１９に記載の装
置。
【請求項２１】前記の記録された情報が、所定の期間
において前記ネットワーク内で切り離された呼出しに対
する前記ネットワーク内で接続された呼出しの比率であ
ることを特徴とする、請求項１９に記載の装置。
【請求項２２】前記の記録された情報が、あらかじめ
選択された期間において前記ネットワーク内で行われた
呼出しのボリュームであることを特徴とする、請求項１
９に記載の装置。
【請求項２３】前記予測値が、前記第１の期間におい
て記録された情報を利用することに加えて、前記第１の
期間の前に記録された情報を利用して生成されることを
特徴とする、請求項１９に記載の装置。