粒子滤波ブロックを使用した仿真软件でのパラメーターお金宝appよび状態の推定

この例では次を使用します:

この例では,系统辨识工具箱™の粒子滤波ブロックを使用する方法を示します。離散時間伝達関数のパラメ，タ，推定問題を再定式化し，状態の推定問題として再帰的に解きます。

はじめに

系统辨识工具箱には,非線形の状態推定を行うための仿真软件ブロックが3つあります。金宝app

粒子滤波:離散時間の粒子フィルタ，アルゴリズムを実装します。
扩展卡尔曼滤波器:1次離散時間の拡張カルマンフィルターアルゴリズムを実装します。
无味卡尔曼滤波:離散時間のアンセンテッドカルマンフィルターアルゴリズムを実装します。

これらのブロックは,さまざまなサンプルレートで動作する複数のセンサーを使用した状態推定をサポートします。これらのブロックの使用に関する一般的なワ，クフロ，は次のとおりです。

MATLABまたはSi金宝appmulink関数を使用してプラントとセンサ，の動作をモデル化します。
ブロックのパラメ，タ，を構成します。
フィルタ，をシミュレ，トして結果を解析し，フィルタ，の性能の信頼度を得ます。
ハ，ドウェアにフィルタ，を展開します。金宝appSimulink Coder™ソフトウェアを使用して，これらのフィルタ，のコ，ドを生成できます。

この例では粒子滤波ブロックを使用し,このワークフローの最初の2つの手順を示します。最後の2の手順にいては”次のステップの節で簡単に説明します。この例の目標は，離散時間伝達関数(出力誤差モデル)のパラメ，タ，を再帰的に推定することです。ここでは新しい情報が着信すると各タムステップでモデルのパラメタが更新されます。

拡張カルマンフィルターに関心がある場合,“複数のマルチレートセンサーをもつ非線形システムの状態の推定“の例を参照してください。アンセンテッドカルマンフィルタ，の使用手順は拡張カルマンフィルタ，の場合と似ています。

サンプルファ转换器ルフォルダ转换器をmatlabパスに追加します。

目录(fullfile (matlabroot,“例子”，“识别”，“主要”));

プラントのモデル化

大半の状態推定アルゴリズムは,タイムステップ間におけるプラントの状態の変化を記述するプラントモデル(状態遷移関数)に依存します。この関数は通常， $ x [k + 1] = f (x [k], w [k], u [k])美元で表されます。ここでx は状態、w はプロセスノイズ、u はたとえばシステム入力やパラメーターといったオプションの追加入力です。粒子フィルターブロックには、この関数を少し異なる構文で $ X [k + 1] = f {pf} u (X [k], [k])美元のように指定しなければなりません。違いは以下のとおりです。

粒子フィルターは多くの状態仮説(粒子)の軌跡を追うことによって機能し,ブロックはすべての状態仮説を一度に関数に渡します。具体的には，状態ベクトルに個の要素があり粒子の使用を選択した場合，は次元 $美元[N_s \;N_p]美元$ をも，その各列が状態仮説になります。
関数内の状態仮説に対するプロセスノ@ @ズの影響を計算します。ブロックでは、プロセスノイズの確率分布にいての仮定が行われず，入力としてを必要としません。

関数 f {pf}(…)美元美元はMATLAB编码器™の制限に準拠するMATLAB関数か,仿真软件功能ブ金宝appロックにすることができます。を作成した後，粒子过滤器ブロック内に関数名を指定します。

この例では，離散時間伝達関数のパラメ，タ，推定問題を，状態の推定問題として定式化し直します。この伝達関数は,離散時間プロセスのダイナミクスを表現するか,ゼロ次ホールドなど何らかの連続時間ダイナミクスに信号リコンストラクターを連動させたものを表現することができます。次の1次離散時間の伝達関数のパラメタ推定に関心があると仮定します。

$$ y [k] = \压裂{20问^ {1}}{1 - 0.7 - q ^ {1}} u e [k] + [k] $$

ここで美元$ y [k] はプラント出力， u [k]美元はプラント入力，美元$ e [k] は測定ノ电子邮箱ズ， $美元问^ {1}$$ は $美元问^ {1}u [k] = [t - 1]美元$ が成り立むだ時間演算子です。伝達関数を $\frac{n\;q^{-1}}{1+d\;q^{-1}}$ としてパラメ，タ，化します。ここで n美元と $ d $ は推定するパラメ，タ，です。伝達関数とパラメーターは,状態ベクトルの選択によって,必要な状態空間形式で複数の方法により表現することができます。1の選択肢として $$ x [k] = [\;y [k];d [k];n [k] \;]美元$ があり，ここで2番目と3番目の状態がパラメ，タ，の推定値を表します。したがって，伝達関数を同じく次のように記述することができます。

$$ x[k+1] = \left[ c \开始{数组}{}& # xA;-x_2[k] x_1[k] + x_3[k] u[k] \\ x_2 [k] \ \ & # xA;x_3 [k] \ \ & # xA;结束\{数组}& # xA;\] $$

測定ノ@ @ズ項美元$ e [k] はセンサ，モデリングで処理されます。この例では上記の式を,計算の効率性のためにベクトル化された形式を用いてMATLAB関数で実装します。

类型pfBlockStateTransitionFcnExample

函数xNext = pfBlockStateTransitionFcnExample (x, u) % pfBlockStateTransitionFcnExample粒子状态转换函数%过滤器,估算参数%的输出,一阶,离散时间传递% % %函数模型输入:% x -粒子,NumberOfStates-by-NumberOfParticles矩阵% u -系统输入,一个标量输出% %:% xNext——预测粒子,与相同的维数作为输入x %实现状态转换函数% xNext = [x (1) * (2) + x (3) * u;% x (2);% x (3)];%的向量化形式(所有粒子)。xNext = x;xNext (1) = bsxfun (@times x (1:) - x (2:)) + x (3:) * u;%添加一个小的过程噪声(相对于每个状态的预期大小)，以%增加粒子多样性xNext = xNext + bsxfun(@times，[1;1飞行;1 e 1], randn(大小(xNext)));结束

センサ，モデリング

粒子滤波ブロックには,各状態仮説の尤度(確率)を計算する測定尤度関数を指定しなければなりません。この関数は $L[k] = h_{pf}(X[k]，y[k]，u[k])$ の形式をとります。 L [k]美元は N_p美元要素のベクトルで，は選択する粒子の数です。内の $m ^ {th}$ 要素は，美元$ X [k] 内の $m ^ {th}$ 粒子(列)の尤度です。美元$ y [k] はセンサ，測定値です。 u [k]美元はオプションの入力引数で，これは状態遷移関数の入力と異なるものを指定できます。

センサ，はこの例の最初の状態を測定します。この例では,実際の測定値と予測された測定値の誤差がガウス分布に従うと仮定していますが,尤度の計算には任意の確率分布またはその他の方法を使用できます。 $美元$ h_ {pf}(…)$ を作成し，粒子过滤器ブロック内に関数名を指定します。

类型pfBlockMeasurementLikelihoodFcnExample；

function likelihood = pfblockmeasurementlikehoodfcnexample(粒子，测量)% pfblockmeasurementlikehoodfcnexample测量似然函数用于粒子过滤器% %测量是第一个状态% %输入:%粒子-一个numberofstats -by- numberofparticles矩阵% measurement -系统输出，一个标量% %输出:% likelihood -一个含有NumberOfParticles元素的向量，其中第n个%元素是第n个粒子的可能性%#编码原%预测测量yHat =粒子(1，:);根据实际测量值%与预测测量值%之间的误差计算每个粒子的似然性，假设误差分布在多元正态分布中，%平均值为零，方差为1。求对应的概率%密度函数e = bsxfun(@minus, yHat, measurement(:)');%误差数为1;Mu = 0;%平均西格玛=眼睛(数量的测量);%方差测量enterprod = dot((e-mu)， Sigma \ (e-mu)， 1);c = 1/√((2*pi)^numberOfMeasurements * det(Sigma));似然= c * exp(-0.5 * measurementErrorProd); end

フィルタ，構築

粒子过滤器ブロックを推定用に構成します。状態遷移関数と測定尤度関数の名前、粒子の数およびその初期分布を指定します。

ブロックダ@ @アログの[システムモデル]タブで，次のパラメ，タ，を指定します。

状態遷移

[関数]に状態遷移関数pfBlockStateTransitionFcnExampleを指定します。関数名を入力して[適用]をクリックすると，ブロックは関数に追加の入力があることを検出して入力端子“StateTransitionFcnInputs”を作成します。この端子にシステム入力を接続します。

初期化

[粒子数]に10000を指定します。通常は粒子の数が多いほど優れた推定が得られますが，計算コストは増加します。
[分布]に高斯を指定して，多変量ガウス分布から粒子の初期セットを取得します。その後，[平均]に[0;0;0]を指定します。これは、3 つの状態があり、これが最適な推定であるためです。[共分散]にDiag ([10 5 100])を指定して3番目の状態の推定の分散(不確かさ)が大きく,最初の2つは分散が小さいことを指定します。この粒子の初期セットは,可能な真の状態をカバーするのに十分な幅に広がっている(分散が大きい)ことが重要です。

測定 1

[関数]に，測定尤度関数の名前pfBlockMeasurementLikelihoodFcnExampleを指定します。

サンプル時間

ブロックダ@ @アログの下部にある[サンプル時間]に1を入力します。状態遷移関数と測定尤度関数の間にサンプル時間の相違がみられる場合や、サンプル時間の異なる複数のセンサーを使用する場合、これらを[ブロック出力，マルチレ，ト]タブで構成できます。

粒子フィルタ，は，尤度の低い粒子を削除し，尤度の高いものを使って新しい粒子をシ，ドします。これは[リサンプリング]グル，プの下のオプションによって制御されます。この例では既定の設定を使用します。

既定では，ブロックは状態仮説の平均値をその尤度で加重したものだけを出力します。すべての粒子や重みを表示したり，状態推定の別の抽出方法を選択するには，[ブロック出力，マルチレ，ト]タブのオプションを使用します。

シミュレ，ションと結果

簡単なテストとして，ホワ@ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @プラントからの入力およびノイズを含む測定値が粒子滤波ブロックに渡されます。次のS金宝appimulinkモデルはこのセットアップを表します。

open_system (“pfBlockExample”）

システムをシミュレ，トして，推定されたパラメ，タ，と実際のパラメ，タ，を比較します。

sim卡(“pfBlockExample”) open_system (pfBlockExample /参数估计的）

プロットには分子と分母の実際のパラメ，タ，と，その粒子フィルタ，推定値が表示されています。推定値は10タescムステップ後に真値にほぼ収束します。初期状態推定が真値から遠い場合でも収束が得られます。

トラブルシュ，ティング

ここではアプリケーションで粒子フィルターが予想どおりに機能しない場合に備えて,潜在的な実装の問題とトラブルシューティングのアイデアをいくつか示します。

通常，粒子フィルタ，のトラブルシュ，ティングを行う場合，粒子とその重みのセットを調査します。これらを取得するには，ブロックダ[ブロック出力，マルチレ，ト]タブにある[すべての粒子を出力]と[输出所有权重]を選択します。

最初の正常性チェックとして,状態遷移関数と測定尤度関数がシステムの動作を適切に捉えていることを確認します。システムのシミュレーションモデルがある(よってシミュレーション内の実際の状態にアクセスできる)場合,実際の状態でフィルターを初期化してみることができます。その後,状態遷移関数によって実際の状態の時間伝播が正確に計算され,また測定尤度関数によってこれらの粒子で高い尤度が計算されるかどうかを検証できます。

粒子の初期セットは重要です。シミュレションの開始時に少なくとも一部の粒子が高い尤度をもことを確認してください。実際の状態が状態仮説の初期の広がりの外にある場合,推定値が不正確になったり発散したりする可能性があります。

状態推定の精度が最初は良好であるが,時間とともに劣化する場合,粒子が縮退したり粒子の精度が悪化したりする問題が起きている可能性があります[2]。粒子の縮退は粒子の分布が広すぎる場合に発生し,粒子の精度悪化はリサンプリングの後に粒子が塊を形成するために発生します。粒子の縮退は，直接リサンプリングの結果，粒子の精度悪化にながります。この例で使用されている状態遷移関数内での人工的なプロセスノイズの追加は,実用的なアプローチの1つです。このような問題の解決方法については多くの文献が公開されており,アプリケーションによってはさらに体系的な方法を利用できる可能性もあります。[1]，[2]は有用な参考資料です。

次のステップ

状態推定の検証:シミュレーションでフィルターから予想どおりの性能が得られるようになったら,通常は幅広いモンテカルロ法シミュレーションを使用してこの性能をさらに検証します。詳細にいては，金宝appSimulinkでのオンラeconlン状態推定の検証を参照してください。粒子过滤器ブロックダ安大安大市アログの[ランダム性]グル，プのオプションを使用してこれらのシミュレ，ションを実行できます。
コードの生成:粒子滤波ブロックは,仿真软件编码器™ソフトウェ金宝appアを使用したCおよびc++コードの生成をサポートします。このブロックに提供する関数は,MATLAB编码器™ソフトウェア(MATLAB関数を使用してシステムをモデル化している場合)および仿真软件编码器ソ金宝appフトウェア(仿真软件函数ブロックを使用してシステムをモデル化している場合)の制限に準拠しなければなりません。

まとめ

この例では,系统辨识工具箱の粒子滤波ブロックを使用する方法を示しました。ここでは離散時間の伝達関数のパラメーターを再帰的に推定し,新しい情報の着信に合わせて各タイムステップでパラメーターを更新しました。

close_system (“pfBlockExample”, 0)

Matlabパスからサンプルファ转换器ルフォルダ转换器を削除します。

rmpath (fullfile (matlabroot,“例子”，“识别”，“主要”));

参考文献

西蒙，丹。最优状态估计:卡尔曼，H∞和非线性方法。约翰·威利父子出版社，2006年出版。

[2]杜塞，阿诺，亚当·m·约翰森。《粒子滤波和平滑教程:15年后》非线性滤波手册12.656-704(2009):3。

参考

扩展卡尔曼滤波器|无味卡尔曼滤波器|粒子滤波