制御系解析基礎 - [yamakita-lab.ctrl.titech]Home

システムモデルと確率過程
東京工業大学機械制御システム専攻
山北昌毅
よく用いられるシステムのモデル
（計算機で信号をサンプルする際）
AR(Auto Regressive)モデル
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n ),
(  i  0, i  0,
, n  1)
MA(Moving Average)モデル
y ( k )   n 1e( k  1) 
  0 e( k  n ),
( i  0, i  0,
, n  1; u (i )  e(i ), i  k  1,
, k  n)
ARX(exogenous)モデル
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n )   n 1u ( k  1) 
  0 u ( k  n )  e( k )
ARMAモデル
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n )   n 1e( k  1) 
  0e(k  n )
ARMAXモデル
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n )   n 1u ( k  1) 
  0u ( k  n )   n 1e( k  1) 
  0e(k  n )
工学系で対象にされるシステムは連続系のシステム。本当にこれで表現できるの？
第１回講義の内容
1.
2.
3.
4.
5.
6.
7.
8.
9.
10.
11.
状態空間表現とＬＴＩシステム
ＬＴＩシステムと伝達関数・可観測正準系
ＬＴＩシステムの解とシステムの離散化
可観測正準系と入出力モデル
Ｚ変換・パルス伝達関数・シフトオペレータ・微分オペレータ
式誤差モデル・出力誤差モデル
伝達関数とマルコフパラメータ
確率過程
エルゴード性
確率変数の収束
推定量の性質
状態空間表現とLTIシステム
u(t )
y(t )
d
システム
 x  f ( x, u, t ) 状態方程式
 dt
（ベクトル値関数の一階の常微分方程式）

 y  h ( x, u , t )
観測方程式
x(t )  R n , u (t )  R m , y (t )  R p
 x1 (t ) 
 x (t ) 
LTI（Linear Time Invariant) システム
2




x (t ) :
 x  Ax  Bu   






xn (t ) 





 y  Cx  Du 

 n
 A1 j x j 
 j 1
 n
 A2 j x j 
 j 1

 n

 Anj x j 
 j 1

B1 j u j 
j 1

m

B2 j u j 
j 1



m

Bnj u j 
j 1

m
m
 n

  C1 j x j   D1 j u j 
j 1
 j 1

m
 n

  C2 j x j   D2 j u j 
j 1
 j 1



 n

m


A
x

B
u
  pj j  pj j 
j

1
j

1


ＬＴＩシステムのブロック線図表現
D
u
x
B


A
x
C

y
状態空間表現（１）

１．厳密な数学モデル
(ml 2  I )    mgl sin( )  
２．近似数学モデル
(ml 2  I )    mgl  
状態空間表現（２）
 x1   
x    :  
 x2   
１．厳密な数学モデル
x2
 x1    

x  

2
x
(


x

mgl
sin(
x
)


)
/
(
ml

I
)
2
1
 2    

y    x1
状態空間表現（３）
２．近似数学モデル
x2
 x1    

x  

2
x
(


x

mglx


)
/
(
ml

I
)

2
1
 2   

y    x1

0
1
0

  x1  



x  




2
2
2
 mgl / (ml  I )  / (ml  I )   x2  1 / (ml  I ) 


 y  1 0  x1   [0]
 x 

 2

ラプラス変換の性質
（６）合成則
（１）線形性
L  f * g  L  f  L  g
L c1 f1  c2 f 2   c1 L  f 1  c2 L  f 2 
t
ただし、 f * gを  f (t   ) g ( )dとする。
（２）相似則
0
L  f (at ) 
1 s
F   , F ( s ) : L  f 
a a
（畳み込み積分）
（７）反転公式
（３）推移則
L  f (t   )  e  s L  f 


L  f (t   )  e s L  f    f (t )e  st dt
0
（４）積分法則
L
 f ( )d   1s L  f 
t
0
（５）微分法則
 df 
L    L  f (1)   sL  f   f (0)
 dt 
n 1
L  f ( n )   s n L  f    f ( k ) (0) s n 1 k
k 0
1 c 
F ( s )e st ds

c

2
（８）初期値の定理
f (t ) 

f (0 )  lim sF ( s )
s 
（９）最終値の定理
f ()  lim sF ( s )
s 0
ただし、 sF ( s )は閉複素右半面で解析的
（１０）パーセバルの定理

2
1 
f
(
t
)
dt

F ( j ) d 
0
2 
ただし、 F ( s )はs  jで定義されるものとする
2
良く使うラプラス変換と逆変換
f
L f 
 (t )
1
U (t )(step func.)
eat
t n eat
1
s
1
sa
n!
( s  a) n 1
f
L f 
sin t

s2   2
s
cos t
s2   2
 cos   ( s  a) sin 
at
e sin(t   )
( s  a)2   2
n!
n
t
s n 1
ラプラス変換の利用法
時間領域での表現
時間領域での表現
畳み込み積分
ラプラス変換
逆ラプラス変換
周波数領域での表現
周波数領域での表現
掛け算
状態空間表現と伝達関数行列
 x  Ax  Bu

 y  Cx  Du

( A, B, C , D)が定数行列であるとき、
時不変線形( LTI : Li near Ti me I nvar i ant ) システムという
x(0)  0の下で両辺を Lapl ace変換すると
sX ( s )  AX ( s )  BU ( s )  X ( s )  ( sI  A) 1 B
Y ( s )  CX ( s )  DU ( s )  C ( sI  A) 1 BU ( s )  DU ( s )
: H ( s )U ( s),
H ( s ) : C ( sI  A) 1 B  D
可観測正準系（１）
ＬＴＩシステムの座標変換
x  Tx ( xは新しい状態変数ベクトル）を用いて
xを用いてシステムを表現することをシステムの
座標変換という
Tx  ATx  Bu  x  T 1 ATx  T 1 Bu
 x  Ax  Bu



 y  Cx  Du  y  CTx  Du
 y  CTx  Du
 x  Ax  Bu

, ( A : T 1 AT , B : T 1B, C : CT )
 y  Cx  Du
可観測正準系（２）
[ Fact ] 完全可観測な1入力１出力システムは以下のように
座標変換することができる

0 0

1 0


 x  0 1







 y   0 0
0
0
0
0
0
1
0
a0 
 b0 
 b 
a1 
 1 
x 
u



 an 1 
bn  2 
 bn 1 
an 
1 x  du
連続系の時間領域の入出力表現
y  xn
y  xn  xn 1  an 1 xn  bn 1u  xn 1  an 1 y  bn 1u
y  xn  2  an  2 y  an 1 y  bn  2u  bn 1u
y ( n 1)  x1  a1 y 
 an  2 y ( n 3)  an 1 y ( n  2)  b1u  bn  2u 
 bn 1u ( n  2)
 an  2 y ( n  2)  an 1 y ( n 1)  b0u  b1u 
 bn 1u ( n 1)
y ( n )  a0 y  a1 y 
y ( n )  an 1 y ( n 1)  an  2 y ( n  2) 
 a1 y  a0 y  bn 1u ( n 1)  bn  2u ( n  2) 
 b1u  b0u
（注意：微分方程式表現で入力の微分項があっても実現には微分器は不要）
両辺を y (0)  y (0) 
( s n  an 1s n 1 
 y ( n 1) (0)  0の下でLapl ace変換すると
 a1s  a0 )Y ( s )  (bn 1s n 1 
bn 1s n 1   b1s  b0
Y ( s)
: H ( s )  n
U (s)
s  an 1s n 1   a1s  a0
 b1s  b0 )U ( s )
状態方程式の一般解
物理システムは連続系のシステムとして表現されることが普通
 x(t )  Ax(t )  Bu (t ), x(0)  x0

 y(t )  Cx(t )  Du (t )
t
x(t )  e x0   e
At
0
A ( t  )
Bu ( )d （一般解）
2
3
(
At
)
(
At
)
At
e : I  At 

 ...
2
3!
d At
A3t 2
( At )2
2
( e )  0  A A t 
... A ( I  At 
...)  Ae At
dt
2!
2
（遷移行列）
一般解の証明
公式
t 
d t
h(t , )d  h(t , t )  
h(t , )d

0
0
dt
t

t
d
d At
x(t ) 
e x0   e A(t  ) Bu ( )d
0
dt
dt

t
 Ae x0  A e A( t  ) Bu ( )d  Bu (t )
At

0
t

 A e x0   e A(t  ) Bu ( ) d  Bu (t )
At
0
 Ax(t )  Bu (t )
入出力関係の畳み込み積分表現
t
x(t )  e At x0   e A(t  ) Bu ( )d
0
y (t )  Cx(t )  Du (t )
x0  0と仮定すると
t
t
0
0
y (t )   Ce A(t  ) Bu ( )d  Du (t )   H (t   )u ( )d
H (t ) : Ce At B  D (t )
A2t 2 A3t 3
 C ( I  At 

  ) B  D (t )
2!
3!
t2
t3
2
3
 CB  CABt  CA B  CA B   D (t )
2!
3!
t2
t3
: H1  H 2t  H 3  H 4   H 0 (t )
2!
3!
H iを連続系のマルコフパラメータという
マルコフパラメータと伝達関数の関係
H ( s ) : C ( sI  A) 1 B  D
 C ( I  A / s ) 1 / sB  D
 C ( I / s  A / s 2  A2 / s 3 
1
)B  D
2
3
1
1
2 1
 CB    CAB    CA B   
s
s
s
1
2
3
1
 D 
s
0
1
1
1
1
: H1    H 2    H 3     H 0  
s
s
s
s
H iは i個の積分器を通った後の重みを表す
n
 1  t
L  n 1  
 s  n!
1
2
t
L1  H ( s )  CB  CABt  CA2 B 
2!
 D (t )
0
システムの離散化（１）
u(t )
u (k )
0
12 3
k
0
ＺＯＨ
Ｄ／Ａ
y(t )
T2T 3T
0
t
プラント
y (k )
T2T 3T
t
0
12 3
Ａ／Ｄ
計算機
計算機でy(k)を観測してu(k)を決定することになる
k
システムの離散化（２）
kT
x(kT )  e A( kT ) x0   e A( kT  ) Bu ( )d
0
( k 1)T

x((k  1)T )  e A(( k 1)T ) x0 
e A(( k 1)T  ) Bu ( )d
0
kT
( k 1)T
0
kT
 e AT e A( kT ) x0  e AT  e A( kT  ) Bu ( ) d 
( k 1)T

 e AT x( kT ) 
e A(( k 1)T  ) d Bu (kT )
kT
 x( kT )  u ( kT )
T
 : e AT ,  :  e A d B
0

e A(( k 1)T  ) Bu ( ) d
システムの離散化（３）
 x  Ax  Bu

 y  Cx  Du
＋
T
x (( k  1)T )   x ( kT )  u ( kT )
y ( kT )  Cx ( kT )  Du (kT )
T
 : e AT ,  :  e A d B
0
離散時間系の時間領域の入出力モデル

 o 
0
 0 

1 0

  


1



 1 
 Ao   1 ...
  x(k  1)  A0 x(k )  b0u (k )
.  , b0 = 





.



n2 

 n  2   y (k )  c0 x(k )


  n 1 

1  n 1 

 co   0 0 ... 0 1 （可観測正準系）
y (k )  xn (k )
 xn 1 (k  1)   n 1 xn (k  1)   n 1u (k  1)
 xn 1 (k  1)   n 1 y (k  1)   n 1u (k  1)
 xn  2 (k  2)   n  2 xn (k  2)   n  2u (k  2)   n 1 y (k  1)   n 1u (k  1)
 xn  2 (k  2)   n  2 y (k  2)   n  2u (k  2)   n 1 y (k  1)   n 1u (k  1)
 x1 (k  n  1)  1 y (k  n  1)  1u (k  n  1) 
  n 1 y (k  1)   n 1u (k  1)
 0   0 y (k  n)   0u (k  n)  1u (k  n  1) 
  n 1 y (k  1)   n 1u (k  1)
  n 1 y (k  1)   n  2 y ( k  2) 
  0 y (k  n)   n 1u (k  1) 
  0u ( k  n)
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n )   n 1u ( k  1) 
  0u ( k  n )
AR(Auto Regressive)モデル
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n ),
(  i  0, i  0,
, n  1)
MA(Moving Average)モデル
y ( k )   n 1e( k  1) 
  0 e( k  n ),
( i  0, i  0,
, n  1; u (i )  e(i ), i  k  1,
, k  n)
ARX(exogenous)モデル
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n)   n 1u ( k  1) 
  0 u ( k  n)(  e( k ))
ARMAモデル
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n )   n 1e( k  1) 
  0e(k  n )
ARMAXモデル
y ( k )   n 1 y ( k  1)   n  2 y ( k  2) 
  0 y ( k  n )   n 1u ( k  1) 
  0u ( k  n )   n 1e ( k  1) 
  0e(k  n )
Ｚ変換・パルス伝達関数、シフト・微分オペレータ
数列{y(k),k=0,1, }のＺ変換

Y ( z ) :  y (k ) z  k
k 0
関数y(t)(t  0)のラプラス変換

Y(s)=  y (t )e  st dt
0
パルス伝達関数＝ゼロ状態での入力のＺ変換伝達関数＝ゼロ状態での入力のラプラス変換と
と出力のＺ変換の比
安定性：極が複素単位円内に存在する
出力のラプラス変換の比
安定性：極が複素左半平面に存在する
q  z , q 1  z 1 (初期関数値を無視）
q  s, q 1  s 1
qy (k )  y (k  1), q 1 y (k )  y (k  1)
sy (t )  y (t ), s 1 y (t )   y ( )d
t
0
式誤差モデル・出力誤差モデル（ＢＪモデル）（１）
式誤差モデル：ダイナミックスの前に外乱が入る構造
（式の誤差として外乱が入る構造）
C ( z ) w (t )
B ( z )u (t )
1
A( z )
y (t )
y ( k )   n 1 y ( k  1) 
 1 y ( k  n  1)   0 y (k  n)   n 1u (k  1) 
  0 u (k  n)  w(t )
y ( k )   n 1 y ( k  1) 
 1 y ( k  n  1)   0 y (k  n)   n 1u (k  1) 
  0 u (k  n)  w(t )
出力誤差モデル：ダイナミックスの後に外乱が入る構造
（出力に誤差が入る構造）
w (t )
B ( z )u (t )
1
A( z )
y (t )
B ( z )u (t )
1
A( z )
C ( z ) w (t )
1
D( z )
y (t )
BJ（Box Jenkins)モデル
式誤差モデルはＢＪモデルの特殊な場合！
D( z )  A( z )
パラメータ同定用モデル
ARＸモデル
y (k )   T ( k  1)  e(k )
 T ( k  1) : [ y ( k  1), y ( k  2),
 T : [ n 1 ,
,  0 ,  n 1 ,
, y ( k  n), u ( k  1),
, u ( k  n)]
, 0 ]
 (k  1)を回帰ベクトルと呼ぶ
NARXモデル
y ( k )  h ( y ( k  1), y ( k  2),
, y ( k  n ), u ( k  1),
, u ( k  n ))
NARMAXモデル
y ( k )  h ( y ( k  1), y ( k  2),
, y ( k  n ), u ( k  1),
, u ( k  n ), e ( k  1),
, e ( k  n ))
パルス伝達関数とマルコフパラメータ（１）
y (k )   n 1 y (k  1)   n  2 y (k  2) 
  0 y (k  n)   n 1u (k  1) 
  0u ( k  n)
初期条件ゼロの下で両辺をＺ変換する
ただし、数列y (k )のZ変換は

Z  y (k ) : Y ( z ) :  y (k ) z  k
k 0

Z  y (k  1)   y (k  1) z
k
k 0

  y (k ') z
 k ' 1
k '1

zy (0)  zy (0)   y (k ') z  k '1
k '1



k ' 
  zy (0)  z  y (0)   y (k ') z    zy (0)  z  y (k ') z  k '   zy (0)  zY ( z )
k '1
k ' 0



Z  y (k  1)   y (k  1) z
z
k 0
n
  n 1 z n 1 
k

 y (k  1) z
k
k 1

  y (k ') z  k '1 z 1Y ( z )
k ' 0
 1 z   0  Y ( z )   bn 1 z n 1 
bn 1 z n 1  b1 z  b0
Y ( z)

: H ( z )
U ( z ) z n   n 1 z n 1   1 z   0
b1 z  b0 U ( z )
パルス伝達関数とマルコフパラメータ（２）
x(k  1)  x(k )  u (k )
y (k )  Cx(k )  Du (k )
u (0)  1, u (k )  0(k  1)の入力を上式に入れる
y (0)  D, x(1)  
y (1)  Cx(1)  C , x(2)  
y (2)  Cx(2)  C , x(3)   2
y (n)  Cx(n)  C  n 1
H 0 : D, H i : C  i 1(i  1)をマルコフパラメータという
H iはiステップ遅れて出力に影響を与える重みとなっている
パルス伝達関数とマルコフパラメータ（３）
 x(k  1)  x(k )  u (k )

 y (k )  Cx(k )  Du (k )
x(0)をゼロとして両辺をＺ変換すると
 zX ( z )  X ( z )  U ( z )

Y ( z )  CX ( z )  DU (k )
Y ( z )   C ( zI   ) 1   D U ( z )
Y ( z)
)
 C ( zI   ) 1   D : H ( z（パルス伝達関数）
U ( z)
2
1
2
1
1
1
C ( zI   )   D  C     C     C 2   
z
z
z
1
1
2
2
1
1
1
 H 0  H1    H 2    H 3   
z
z
z
1
 D 
z
0
確率過程（１）
水の上の花粉の軌跡は最初の位置が同じでも、‘熱的ノイズ’によってその軌跡は
非常に異なるものとなる。このような現象を数学的に取り扱いたい。
tは時間で、iは何回目の観測かを示すと考える。
これを数学的に、p(t ,  )のtを固定したとき、を
px (t )
p (t , 1 )
確率パラメータとしてもつ確率変数と考える。
p x (t , 1 )
p x (t ,  2 )
py
p (t ,  2 )
t
px
を固定して、p(t ,  )をtのみの関数と考えたとき、
見本過程という。（実際に観測されるのは見本過程
の一つである！）
確率過程（２）
P ( )
数学的には、確率過程とは時間tを変数に持つ確率変数
dP( )
の集合（族）と定義される。
Y (t ,  ), t  T  (Tは非負の整数とか、実数空間をとる）
   d
上記の設定だけで平均値を定義してみよう。だたし、 の軌道が生成される
確率を P( )により定義する。また、 の取りうる全体集合を形式的に とする。
P( )


E Y (t ) :  Y (t ,  )dP( )    Y (t ,  )
d  :  Y (t ,  ) p( )d  



 

このとき、いったい全体集合はどのような集合で, どうやって全体で積分を
実行したらいいだろうか？
確率過程（３）
これに対して答えを与えるのが次のDoob-Dynki nの補題である。
[ 補題]( Doob  Dynkin)
関数X , Y :   R nを考える。 X が確率変数であるとする。
この時、 Yが確率変数であるための必要十分条件は、 R n上
での確率変数g : R n  R nが存在して、以下の関係を満たす
ことである。
Y ( )  g ( X ( ))
確率過程（４）
この補題を使って x  X ( )  R nとし、   X 1 ( x)とした場合には次のようになる。
（ただし、 X 1は逆関数ではなくて集合関数として一般には定義）


Y ( )dP( )   n g ( X ( ))dP( X 1 ( x))   n Y ( x)dP ( x)
R
R
ただし、
Y ( x) : g ( x)


P X 1
1
dx
dP ( x) : dP( X ( x)) 
 x

つまり、 での確率分布関数から xでの確率分布関数( 確率密度関数）
が誘導される。
従って、 の分布関数を考えなくとも、 xの分布関数を考えて
期待値計算を行って良いことになる。 ( つまり、 Y ( x)の実現値が
観測されていると考える）
確率過程（５）
ここで、上記の結果を次の期待値計算に適用する。
E Y (t ) :  Y (t ,  )dP( )   Yt ( ) dP( )


xt  X t(  ) , Yt(  )  g ( X t ( ))とし、  =X t1 ( xt )として考え、前と同様に変形する。
E Y (t ) :  Y (t ,  )dP( )   n g ( X t ( )) dP( X t 1 ( xt ))   n Yt ( xt )dP ( xt )

R
R
ここで、
Yt ( xt ) : g ( xt )

1
dPt ( xt ) : dP ( X t ( xt ))
とし、 y (t ) : yt  Yt ( xt ), xt =Yt 1 ( yt )  Yt 1 ( y (t ))と考えると
E Y (t )   n y (t )dPt (Y 1 ( yt ))   n y (t )dPt ( y (t ))
R
R
となって、 y (t )の分布関数さえ分かれば計算可能となる。
ただし、 nの大きさは考える問題によって異なる。たくさんの時刻の結合（同時）確率分布を考える
と nは非常に大きくなる！
確率過程（６）
例
今、 x(t ,  )の平均値を計算したいとする。標本信号の考え方では
E{x(t )}   x(t ,  )dP ( )

としなければならないが、標本信号がどのような確率で生成され
るかは分らない。これに代わって
E{x(t )}   n x(t )dP ( x(t ))
R
と計算できることを示している。
確率過程（７）
また、例えば時刻tと sの相関を計算するために、
 x(t ,  ) 
xa (t , s,  ) : 

 x ( s,  ) 
として
E{xa (t ) xaT ( s )}   xa (t , s,  ) xaT (t , s,  ) dP( )

を計算しなければならなくなるが、これは
 x(t ) 
T
x
(
t
,
s
)
x
(
t
,
s
)
dP
(
x
(
t
,
s
)),
x
(
t
,
s
)
:

a
a
a
a
 x( s ) 
R2 n


E{xa (t ) xaT ( s )}  
として計算できることを示している。もちろん、 xa (t , s )の分布
を知る必要があるが。
確率過程（８）
一般に確率過程はy( t 1 ,  ) , , y( t n ,  ) の同時確率によってその性質が決定される。この同時確率が時間の推移に
関して不変であるとき、その確率過程は定常であるという。
[モーメント ]次式で定義される期待値を n次のモーメントという。
E{ y (t1 ,  ) y (t2 ,  )
y (tn ,  )}
特に、１次、２次のモーメントを平均値、相関関数と呼び、次のように表現する
平均値　：  y (t )  E{ y (t ,  )}
自己相関関数：  yy (t1 , t2 )  E{ y (t1 ,  ) y (t2 ,  )}
[ 共分散、分散] 共分散、分散を次式で定義する
自己共分散：  yy (t1 , t2 )  cov[ y (t1 ), y (t2 )] : E{( y (t1 )   y (t1 ))( y (t 2 )   y (t 2 ))}
分散　：  y 2 (t ) : cov[ y (t ), y (t )]
確率過程が定常であるとき、相関関数は時刻の差のみの関数となる、つまり   t2  t1とすると
 yy (t1 , t2 )   yy ( )   yy ( )
なぜなら、
E{ y (t1 ), y (t2 )}  E{ y (t1 ) y (  t1 )}  E{ y (t '1   ) y (t '1 )}
確率過程（８）
数学的（集合的）平均と時間平均
前のスライドで考えた平均値はいろいろな標本があり得るのに
たいして、それを集合的に考えてその期待値を計算していた。
その期待値を数学的期待値、または集合的期待値という。
しかし、現実の実験や観測では一つの時間関数が観測できるだけで、
同じ実験や観測を何度もすることができな場合も多い。そのような
一つの標本関数を用いて、集合的な性質が計算可能である場合を
エルゴード過程と呼ぶ。
例えば
1
T  2T
E{x(t )}  lim

T
T
x(t   )d
自己共分散とスペクトル密度関数
定常過程の場合は自己共分散は時間差のみの関数であった
xx ( )
これをフーリエ変換したものを考える

 xx ( j )   xx ( )e  j d

これはパワースペクトルと呼ばれる。また、逆フーリエ変換では
xx ( ) 
1
2



 xx ( j )e j d 
従って各時刻の xの分散xx (0)は
xx (0) 
1
2



 xx ( j )1d 
であるので、パワースペクトルの周波数領域での積分で与えられる。
パワースペクトルからの伝達関数の推定
相互相関関数とフィルター出力のパワースペクトル
相互相関関数は次式で定義される
xy ( ) : E{x(t ) y (t   )}
相互相関関数のフーリエ変換には次の関係がある

 xy ( j ) :  xy ( )e  j d

 xy ( j )   yx* ( j )
相互相関関数と自己相関関数の関係には次の関係がある .
ただし、 xから yへの伝達関数を H ( s )でその重み関数を h(t )とする .

xy (t )   h(t )xx (t   )d


 yy (t )   h(t ) yx (t   ) d

実際





1
T  2T
h(t )xx (t   )d   h(t ) lim
1
2T


T
T


T
T
x(t ')  h(t )x(t ' t   )d dt ' 

x(t ') x(t ' t   )dt 'd
1
2T

T
T
x(t ') y (t ' t )dt '  xy (t )
従って

 yy ( j )  F{ h(t ) yx (t   )d }  F{h}F{ yx }  H ( j ) yx ( j )

 H ( j ) xy ( j )  H ( j ) H * ( j ) xx* ( j ) | H ( j ) |2  xx ( j )
*
今yのスペクトル密度関数が分かっているとすると、パワースペクトルが
全周波数に渡って１のノイズ( 白色ノイズ）によって駆動されているとすると
 xx ( j )  1
であるので、
 yy ( j ) | H ( j ) |2
と考えることができる。また、 yの分散 yy (0)は
 yy (0) 
1
2



 yy ( j ) d 
1
2



| H ( j ) |2 d  || H ( s) ||2 2
となり、 H ( s )の H 2ノルムの２乗（の定数倍）となることが分かる
伊藤の確率微分方程式‘超入門’
 x(t  t )  x (t )  f c ( x(t ))t  L(t ) (t )  o(t ), x(t )  R n

p
 y (t  t )  y (t )  hc ( x(t )) t  V (t )  (t )  o( t ), y (t )  R
が任意の小さな正の tについて成り立つ時
dx(t )  f c ( x(t ), t )dt  L(t )d  (t )

dy (t )  hc ( x(t ), t ) dt  V (t ) d (t )
と表現する。ただし、  (t ), (t )は独立なブラウン運動で、
それぞれ対角な拡散行列Q(t ), R (t )を持つ。
 d    (t ), d  (t )d  T (t )  Q (t )dt , Q(t )  diag (q1 (t ), , qn (t ))

 d   (t ), d (t )d T (t )  R (t )dt , R (t )  diag (r (t ), , r (t ))
1
p

 (t )～N (0, Q (t ) t ),  (t )～N (0, R (t ) t )


o ( t )
lim
0
 t 0
t

伊藤の公式（１）
x(t )が伊藤の確率微分方程式を満たす時
dx(t )  f c ( x(t ), t )dt  L(t )d  (t )
y  f ( x)とすると yの微分方程式は次式となる。
f
1 T 2 f
dy 
dx  dx 2 dx
x
2
 x
伊藤の公式（１）
ただし、式を展開した後次の関係を利用。
dt 2  0, dtd i  0, d  i d  j  0(i  j )
 2
d i  qi dt
（ 関数が２次の係数を持つ場合、確率的要素が確定的成分に！）

44
伊藤の公式（２）

f
1 T 2 f
f
1  2 f
dy  dx  dx 2 dx  dx  Tr  2 dxdxT 
x
2
 x
x
2  x


f
1  2 f
 ( f c ( x(t ), t )dt  L(t )d  (t ))  Tr  2 L(t )QLT (t )  dt
x
2  x

Dynkinの公式( 伊藤の公式の応用）
 f

d
1  2 f
E{ f ( x(t ))}  E  f c ( x(t ), t )  Tr  2 L(t )QLT (t ) 
dt
2  x
 x

確率変数の収束（１）
確率変数の列X n (n  1, )とひとつの確率変数X (定数でも良い)を考える。
[確率収束]
任意の   0に対して、
lim{P(| X n  X |  )}  0(lim{P(| X n  X |  )}  1)
n 
n 
P
を満たすとき、 X nは X に確率収束するといい、pl i mX n  X ( X n  X )のように表現する。
n 
[概収束（確率１での収束）]
P(lim X n  X )  1
n 
a.s.
を満たすとき、 X nは X に概収束する（確率１で収束する）といい、 a . s . lim X n  X ( X n  X )のように記述する。
n 
確率収束の場合, X n(  ) を nに関する見本列( サンプルパス）であるとすると、 Xと X nの値が よ
り離れる確率はnが増加すると限りなく小さくなる。しかし、ひとつのサンプルパスについて
X n ( )が X に収束するかどうかは分らない。
確率変数の収束（２）
Lim
a.s.lim
l.i.m.
plim
[自乗平均収束]
E{|| X n ||2 }  , E{|| X ||2 }  で、
lim E{|| X n  X ||2 }  0
n 
を満たすとき、 X nはXに自乗平均収束するといい、l . i . m. X n  X のように表現する。
n 
[分布関数としての収束（法則収束）]
Fn , FをそれぞれX n , X の分布関数とする。このとき
連続な点全てで lim Fn ( )  F ( )(あるいは弱収束する )
n 
を満たすとき、 X nはXに分布関数として収束するといい、Li mX n  X のように記述する。
n 
確率収束しても概収束しない例
t ()
t
t
  [0,1]
lim P{| X t |  }  0
t 
  


p l im X t  0
t 
lim sup X t ( )  1となり、どのサンプルパスも０に収束しない！
t 
推定量の性質
真のパラメータ pとその推定値pˆの誤差を e : p  pˆとする
[ 不偏性]
E{e}  0のとき、 pˆを pの不偏推定量という
[ 一致性]
eNを N 個のデータから推定した際の誤差であるとする。このとき
p lim eN  0
N 
となるとき、 pˆ Nを pの（弱い）一致推定量という。
参考文献
1.
山北：システム制御特論テキスト
http://www.ac.ctrl.titech.ac.jp/~yamakita/text.html
2.
3.
4.
Ｂ．エクセンダール（谷口節男訳）：確率微分方
程式（シュプリンガー・ジャパン）（１９９９）
相良ら：システム同定（コロナ社）（１９９５）
足立修一：ユーザのためのシステム同定理論（
コロナ社）（１９９３）

Download Report