JPH0237402A - パラメータ推定方式 - Google Patents

Abstract

(57)【要約】本公報は電子出願前の出願データであるた
め要約のデータは記録されません。

Description

【発明の詳細な説明】 〔産業上の利用分野〕 本発明は、時系列モデルのパラメータを推定する方式に
関するものである。

〔従来の技術〕

従来のパラメータ推定法では、時々刻々のシステムの離
散人力周期に対応した時系列モデルを考え、そのモデル
の表現に必要なパラメータ全部を推定する考えをとって
いた。この従来のパラメータ推定法が記載された文献と
しては、日本自動詞′４Ｂ協会（現システム制御情報学
会）発行の「システムと制御、２５巻８号〜１２号、２
６巻２号」に「線形離散時間システムの同定手法」とし
て連載された「中溝、高好著」の文献がある。

〔発明が解決しようとする課題〕

しかるに、そのような従来の方法によると、システムの
特性、特にむだ時間と離散入力周期の兼ね合いによって
は、推定するパラメータの数が多すぎて計算に必要な記
憶容量と計算時間が増大し、推定が困難になる恐れがあ
った。

上記欠点は、次の原因から生ずる。

（ａ）　　推定しなければならないパラメータ数にシス
テムのむだ時間と比例する要素があり、このため、むだ
時間が長いシステムはど推定が困難になる。

山）　パラメータ推定法では、むだ時間をあらかじめ既
知と仮定しておく場合もあるが、この仮定は一般の応用
にあたっては適当でない。

本発明は、従来のパラメータ推定方式では応用が難しい
場合においても、■システムの離散入力周期とは異なる
周期を採用し、■推定するパラメータの低次元化を実現
することによって、従来のパラメータ推定法の成果をな
んら無駄にすることな（、より広い適用範囲を実現する
ことを目的とする。

〔課題を解決するための手段〕

このような目的を達成するために本発明は、時系列モデ
ルのパラメータを推定する方式において、推定する対象
に入力される離散入力の周期とは異なる周期を推定周期
として指定する推定周期設定手段と、推定周期間中に推
定する対象に入力された離散入力から移動平均を演算す
るための移動平均演算手段と、推定周期に同期してシス
テムの出力を取り込む同期手段と、推定周期間中に推定
する対象に入力された離散入力から移動平均を演算する
こと及び推定周期に同期してシステムの出力を取り込む
ことを推定周期ごとに繰り返して行なって得られたデー
タから従来のパラメータ推定法を用いてパラメータを推
定するパラメータ推定演算手段とを備え、指定する対象
を低次数でパラメータ推定するようにしたものである。

〔作用〕

本発明によるパラメータ推定方式においては、推定する
パラメータの数を低次元化できるので、計算に必要な記
憶量と計算量もそれに伴って減少し、従来の方法が適用
可能だった範囲では、より高速の演算が可能となり、従
来の方法では適用が困難の事例に対してもパラメータの
推定が可能になる。

〔実施例〕

本発明の主要な機能は成る期間中のデータを移動平均す
る機能と、成る周期ごとにデータを取り込む機能および
パラメータ推定機能からなる。これらの機能のうち、最
初の２つは〔発明が解決しようとする課題〕の欄でも述
べた従来の方法の適用範囲を広げるために用いられるも
のであり、パラメータ推定機能は従来から考えられてい
る方法のうち、不偏推定あるいは一致推定となるものな
らなんでもよい。

以下、本発明の一実施例について図面を用いて詳細に説
明する。第１図は、本発明によるマルチレートパラメー
タ推定方式を構成する要素とその結合を表わしたブロッ
ク系統図である。

第１図において用いられている記号と表現上の規約は次
の通りとする。

円弧に矢印のついたスイッチ入力の記号は連続信号を一
定周期でサンプリングするサンプル回路を意味し、その
周期は記号の側に英字で書かれた値（第１図ではＴ、Ｔ
２）とする。箱の中に「ＺＯＨＪと書かれている記号は
サンプル回路で、サンプリングされた値を保持するホー
ルド回路を意味する。ブロックとブロックとの間の矢印
の付いた線は矢印の向きに情報が伝達されることを意味
し、特にその情報が何であるかを線の側に書く。

次に、第１図に番号をつけたブロック（＝要素）につい
て説明する。要素１は、そのパラメータを推定しようと
する対象（＝システム）である。この要素は１人力１出
力からなり、ここでは入力を記号Ｕ、出力を記号Ｙで表
わす。なお、システムには離散化されたＵ、が入力され
るものとする。

要素２は、本発明におけるマルチレートパラメータ推定
での推定周期を与える推定周期設定手段であり、推定を
行なう以前にあらかじめ適当な値を設定しうる機能を備
えている。この手段で設定された値Ｔ２は、要素３　（
移動平均演算手段）と要素４（同期手段）に伝達され、
それぞれの要素での時間的動作を決定する。

要素３は、要素２から与えられた値Ｔ２の周期ごとに、
システムに入力された時々刻々の入力ｕＪから演算、例
えば移動平均演算を行ない、その演算値Ｘｋを出力する
手段である。

要素４は、要素２から与えられた値Ｔ２の周期ごとに、
要素３からの出力Ｘｋが入力されるのと同期してシステ
ムの出力Ｙｋを取り込む手段である。この要素によって
同期した信号として取り込まれた（Ｘｋ、Ｙｍ）が要素
５．６に出力される。

要素５は、要素６（パラメータ推定演算手段）でパラメ
ータ推定演算を行なう際のモデルの次数を決定する時系
列モデル次数設定手段である。第１図の例では、時系列
の次数として過去の出力の影響を考慮する数としてｎ′
、過去の入力の影響を考慮する数としてｍ“を要素６に
設定している。

この次数は手動で決定することもできるが、手動によら
ない場合はこの次数に基づき要素６が推定したパラメー
タが再び要素５にフィードバックされ、その結果から情
報量基準ＡＩＣあるいはＦ検定法等の方法を用いて試行
的に決定することもできる。このＡＩＣ，Ｆ検定等の方
法については、社団法人計測自動制御学会から発行され
ている「自動制御ハンドブック基礎編」のｐ、５８４〜
ｐ、５８５に説明がされている。

要素６は、要素５から伝達された時系列モデルの次数に
基づいて要素４から伝達される（Ｘｋ、Ｙｋ）から時系
列モデルのパラメータを推定する手段である。パラメー
タ推定の方法としては従来から用いられている方法でよ
い。推定されたパラメータは要素５にフィードバックさ
れ、時系列モデルの次数の決定に用いられる。

最後に、要素５で最終的に決定した時系列モデルの次数
に基づき要素６が推定したパラメータの結果が格納され
る。なお、Ａは要素２，３．４を含む本発明に特に関わ
る手段である。

第２図は、第１図においてシステムの入出力信号の様子
を表わしたものである。第２図のｒｔＪは時間を意味し
、ｊ’、　　ｊ−１，ｊ−２，・・・等の記号は入力信
号Ｕがサンプル回路によってサンプルされた時刻を表わ
す。このサンプル間の時間はすべて一定でその値はｒＴ
Ｊである。また、Ｕｊ−１＋　ｕｊ−２＋・・・等の記
号のついた長方形の図形はサンプル回路によってサンプ
ルされたＵのホールド回路を経た出力を表わし、各々の
図形の高さが出力値を意味する。また、Ｙｊ、Ｙｊ−１
＋　’／ｊ−２゜・・等はＵがサンプル回路に入力され
た時点でシステムの出力値をサンプリングした値を表わ
し、その下に曲線で描かれた軌跡は、システムの出力値
がサンプリングされる時刻以外も含めた時刻でのシステ
ムの出力値の変化の様子を模式的に表現したものである
。このよう−に第２図は、時々刻々変化する離散入力に
よってシステムの出力値が連続的に変化する様子を表わ
したものである。

従来のパラメータ推定法においては、Ｕがサンプル回路
によってサンプリングされる時刻でのシステムの入出力
値を使ってパラメータ推定をしようとするシステムをｆ
ｌ１式のごとく記述し、未知のａＩ＋　ａ！＋”’＋　
　ａｎｒｂｏ＋　ｂ＋＋”’＋　　ｂｍを推定する手法
を考えていた。

ｙ（ｊ）＋ａ＋ｙ（ｊ−１）＋・＝＋ａ＋５ｙ（ｊ−ｎ
）＝ｂｏｕ（ｊ）＋ｂ＋ｕｃｊ４）＋−＋　ｂ、ｕ（ｊ
−ｍ）＋ｅ　（Ｄ　　・−−（１）ここで、（ｕ（Ｄ）
　はシステムの入力データ列、（ｙ（ｊ））はシステム
の出力データ列、（ｅ（ｊ））は測定不可能量である。

（１）式は過去のシステムの入力と出力の線形結合が現
在のシステムの出力値を決定していることを意味するが
、むだ時間が存在するようなシステム等ではかなり以前
の入力まで考慮しなければならず、このようなシステム
のパラメータ推定が困難になっていた。また、むだ時間
がそれほど長くないシステムの場合であっても、Ｕのサ
ンプリング（離散化）周期がパラメータ推定に都合のい
い値になっているとは限らず、この場合離散化周期をパ
ラメータ推定に都合のいい周期に設定しておかなければ
ならなかった。

ここで、パラメータ推定の工業分野への応用を考えてみ
るとき、推定される個々のパラメータが必要であるとさ
れることはあまりない。それよりも、推定されたパラメ
ータに基づき工業的に意味のあるもとのシステムの情報
を復元できることが重要である。

ところで、上のような従来の方法の欠点は個々のシステ
ムの離散入力に異なる係数を考えてパラメータ推定を行
なっている点に原因している。システムの出力値に個々
の入力が影響を及ぼしていることは間違いないが、近い
時刻でのシステムの離散入力値同士はシステムの出力に
同じくらいの影響を持っているとも考えられる。そこで
、時系列Ｕの整数倍の周期Ｔ２を考え、成る時刻とそれ
よりもＴ２時間だけ前の時刻の間に有る時系列Ｕの個々
の値が同じくらいの影響を持って成る時刻でのシステム
の出力に関与しているものと考える。

このことを数式を使って説明する。周期Ｔの時系列Ｕに
対してその整数倍の周期Ｔ２の時系列Ｘを表のようにＴ
２周期ごとに同期させ（最初に同期する時点は任意）、
時系列Ｘの個々の値Ｘｉ＋（ｋ＝１．２，３．　・・・
）を時系列Ｕの個々の値ｕｊ−Ｔｇ＋１ｕｊ−ア２＋２
．・・・、ｕｊからに、ｊの値にはよらないＦを用いて
、 Ｘ！ｌ＝Ｆ　　（ｕｊ−丁２＋ｌ＋　　ｕｊ−ＴｉＺ＋
　　　Ｈ＋　　＋、　　ｕｊ）・・・・（２） のように求める（具体的なＦについては後述する）。

表においては、下線を付けた時点で上下の時系列が同期
している。

このように時系列Ｕから時系列Ｘに変換される様子を示
したのが第３図である。この例では、Ｔ２＝３Ｔである
としている。ただし、関数Ｆについてはこのままでは計
算できない。そこで、このＦの一例を次にあげる。

■ｕｊ−Ｔ２．１＋　　”ｊ−１２＊Ｚ＋　　”　’ｒ
　　ｕｊの単純平均を計算し、それをＸｋとする。

■ｕｊ−Ｔ２＋ｌ＋　　ｕｊ−Ｔ！＋２＋　　”　’　
＋　　ｕｊの重み付は平均を計算し、それをＸｋとする
。

いずれの場合も時系列Ｕの移動平均演算として考えるこ
とができるが、■の単純平均の場合を例にとり関数Ｆの
計算式を（３）式に示す。

Ｘｉ＋＝　　（ｕＪ−ｒｚ＋ｒ”　　ｕｊ−丁ｚ＋ｚ＋
　　”　　’　　”　　ｕｉ）／（Ｔ２／Ｔ）　　・・
・・（３） 以上が第１図で要素３を中心とした要素２，３４の説明
であったが、従来の方法で（１）式の記述されたモデル
が本発明のマルチレートパラメータ推定法でどのように
記述され直したかを（４）式に示す。

ｙ（ｋ）＋ａ＋’ｙ（ｋ−１）＋・＋ａｎ’　　ｙ（ｋ
−ｎ’）＝　ｂｏＸ（ｋ）＋ｂ＋’Ｘ（ｋ−１）＋”・
＋ｂｌＩ’　　Ｘ（ｋ−ｍ’）＋　ｅ（ｋ）・・・・（
４） ここで、（Ｘ（ｋ））はシステムの移動平均入力データ
列、（ｙ（ｋ）ｌ　はシステムの出力データ列、（ｅ（
ｋ））は測定不可能量である。（１）式の表現と（４）
式の表現の違いはモデルの表現に必要なパラメータ数で
ある。例えば、（１１式でｍ＝８として推定しなければ
ならなかったシステムの事例では第３図のようにシステ
ムの入力の変換を行なった場合、ｍｏ　＝３ですむこと
がわかる。（１）式のｎと（４）式のｎ“に関しては明
確な比較ができないが、少なくともｎ゛がｎを越えるこ
とはない。

以上がマルチレートパラメータ推定法の基本的考え方で
ある。考え方を明瞭にするために第１図の要素と関連づ
けて言えば、要素２では推定のための周期Ｔ２を決定す
る。要素３は、第２図に示した新たな時系列Ｘ、を（２
）式に基づいて計算する。

要素４は、Ｘ、が生成されるのと同期してＹｋを取り込
み、同期信号列（ｘ＋＝、ｙｋ）を生成する。要素５は
、同期信号列（Ｘｋ、Ｙｋ）を基本データとして（３）
式でａＩＺ　　ａ２Ｚ　・・’＋　　ａ、ｌ’　＋　　
ｂＯ’＋　　ｂ・・・、ｂ、”を推定する際の次数を決
定する。

要素６については要素４，５より与えられた情゛報を基
にしてパラメータを推定する従来の方法で不偏推定ある
いは一致推定となるものならばなんでもよく、最小二乗
推定、重みつき最小二乗推定、補助変数法、最尤推定法
等を用いることができる。

以上説明したマルチレートパラメータ推定法を、工業分
野におけるプロセス解析の手法として有効であることを
確認するため行なった具体的な一実施例について説明す
る。

ここでは、第１図における工業分野との対応を次のよう
に考える。要素１のシステムは、プロセス解析を行なう
プロセスそのものに相当する。プロセスへの入力信号Ｕ
としてはプロセス操業中のコントローラからの操作出力
信号とし、プロセスからの出力信号Ｙは制御信号を用い
る。

このような対応のもとで、実際のプロセスのかわりにあ
らかじめその特性（伝達関数）のわかっているモデルを
用いてマルチレートパラメータ推定を行なった。（５）
式に示す伝達関数がモデルとして採用したプロセスの伝
達関数である。このプロセスに実際の現場で用いられて
いるデジタルコントローラの操作信号更新周期と同じ０
．５秒の間隔でプロセスに離散信号を入力した。なお、
この間隔は今までの説明で用いた記号によってＴ＝０．
５と表わせる。

この他にも、第１図で説明した各種の設定パラメータに
ついては本実施例の場合、手動入力で次のように設定し
た。

Ｔ２＝１０．ｎ　　＝２．ｍ　　＝４−・・・（６１さ
らに、（２）式の関数Ｆとしては単純平均をとる（３）
式を用いた。

この事例の場合、従来の方法（（１）式）でパラメータ
を推定しようとすると、ｍとしては少なくとも１２個が
必要になるが、本発明によるマルチレートパラメータ推
定方式では必要なｍ゛は６個となっており、さらにｎ、
ｎ’　についても明らかにｎ゛　≦ｎとすればよいから
、推定に必要なパラメータの数は減っていることがわか
る。

次に、このようにパラメータの数を減らしても工業分野
での応用に十分段にたつことを示す。（６）の条件によ
ってマルチレートパラメータ推定を行なった結果は次の
（８）式のようになった。ただし、ここで第１図の要素
６のパラメータ推定演算手段としては周知の最小二乗法
を使った。

Ｇ　（Ｚ−’）　＝Ｍ／Ｎ　　・・・・（８）Ｍ　＝　
０．００５３６２Ｚ−’　＋０．０３５４６９Ｚ−２＋
０．００６０９９Ｚ−’＋０．００００６３Ｚ−’ Ｎ　＝　１−１．２５６６４７２−　’＋０．３２３７
７５２−”（５）式と（８）式では、プロセスを表現し
ている領域が前者は連続時間、後者は離散時間であり異
なるので、比較が簡単にはできない。そこで、（５１，
（８）式のボード線図（工業分野では広く用いられてい
る）を描いたのが第４図、第５図である。第４図が（５
）式のボード線図、第５図が（８）式のボード線図であ
る。これらのボード線図では、ゲイン曲線が実線で位相
曲線が点線で示されているが、両図を比較すると、ゲイ
ン曲線、位相曲線ともほとんど一致しており、工業分野
での応用には充分利用可能な精度で求まっていることが
わかる。

〔発明の効果〕

以上説明したように本発明によるパラメータ推定方式は
、推定する対象に入力される離散入力の周期とは異なる
周期を推定周期として指定し、推定周期間中に推定する
対象に入力された離散入力から移動平均を演算し、推定
周期に同期してシステムの出力を取り込み、上記演算と
取り込みを推定周期ごとに繰り返して行なって得られた
データから従来のパラメータ推定法を用いてパラメータ
を推定し、指定する対象を低次数でパラメータ推定する
ことにより、推定するパラメータの数を低次元化できる
ので、計算に必要な記憶量と計算量もそれに伴って減少
し、従来の方法が適用可能だった範囲では、より高速の
演算が可能となり、従来の方法では適用が困難の事例に
対してもパラメータの推定が可能になる効果がある。

また、推定の計算方法は従来から考えだされている方法
を用いることができるので、従来のパラメータ推定で得
られた資産を継承することができ、また本発明と従来方
式との比較を容易に行なえる効果がある。

さらに、移動平均を演算しているので、高周波のノイズ
を低減化することができ、高信頼性のパラメータ推定が
可能になる効果がある。

【図面の簡単な説明】

第１図は本発明によるパラメータ推定方式の一実施例を
示すブロック系統図、第２図はシステムの入出力信号の
様子を表わすタイムチャート、第３図はシステムの入力
の変換を行なった場合を示すタイムチャート、第４図お
よび第５図はボード線図である。 １・・・推定する対象、２・・・推定周期設定手段、３
・・・移動平均演算手段、４・・・同期手段、５・・・
時系列モデル次数設定手段、 ６・・・パラメータ推定演算手段。

Claims (1)

【特許請求の範囲】
1. 時系列モデルのパラメータを推定する方式において、推
   定する対象に入力される離散入力の周期とは異なる周期
   を推定周期として指定する推定周期設定手段と、推定周
   期間中に推定する対象に入力された離散入力から移動平
   均を演算するための移動平均演算手段と、推定周期に同
   期してシステムの出力を取り込む同期手段と、推定周期
   間中に推定する対象に入力された離散入力から移動平均
   を演算すること及び推定周期に同期してシステムの出力
   を取り込むことを推定周期ごとに繰り返して行なって得
   られたデータから従来のパラメータ推定法を用いてパラ
   メータを推定するパラメータ推定演算手段とを備え、指
   定する対象を低次数でパラメータ推定することを特徴と
   するパラメータ推定方式。