MOTIVATION

论文中的辨识是认为系统的阶数已知的辨识,而且 MCMC 的抽样的变量是特征多项式矩阵,我想使用特征值来抽样,并同时考虑系统阶数未知的情况。

采用RJMMC的方法

细致平稳条件

α(x,y)=min⁡{1,R}π(x)p(x,y)=π(y)p(y,x)\alpha \left( x,y \right) =\min \left\{ 1,R \right\} \\ \pi \left( x \right) p\left( x,y \right) =\pi \left( y \right) p\left( y,x \right)

跨维度游动

k <-> k+1

  1. Birth move

    log⁡R=log⁡P(xk+1)P(y∣xk+1)P(xk)P(y∣xk)+log⁡P(dead)P(birth)+log⁡1k+1q(u)+log⁡∣J1∣\log R=\log \frac{P\left( x_{k+1} \right) P\left( y|x_{k+1} \right)}{P\left( x_k \right) P\left( y|x_k \right)}+\log \frac{P\left( dead \right)}{P\left( birth \right)}+\log \frac{\small{\frac{1}{k+1}}}{q\left( u \right)}+\log \left| J_1 \right|

  2. Death move

    log⁡R=log⁡P(xk)P(y∣xk)P(xk+1)P(y∣xk+1)+log⁡P(birth)P(dead)+log⁡q(u)1k+1−log⁡∣J1∣\log R=\log \frac{P\left( x_k \right) P\left( y|x_k \right)}{P\left( x_{k+1} \right) P\left( y|x_{k+1} \right)}+\log \frac{P\left( birth \right)}{P\left( dead \right)}+\log \frac{\small{q\left( u \right)}}{\frac{1}{k+1}}-\log \left| J_1 \right|

k <-> k+2

  1. Birth move

    log⁡R=log⁡P(xk+2)P(y∣xk+2)P(xk)P(y∣xk)+log⁡P(dead)P(birth)+log⁡1m+1q(θ)q(r)+log⁡∣J2∣\log R=\log \frac{P\left( x_{k+2} \right) P\left( y|x_{k+2} \right)}{P\left( x_k \right) P\left( y|x_k \right)}+\log \frac{P\left( dead \right)}{P\left( birth \right)}+\log \frac{\small{\frac{1}{m+1}}}{q\left( \theta \right) q\left( r \right)}+\log \left| J_2 \right|

  2. Death move

    log⁡R=log⁡P(xk)P(y∣xk)P(xk+2)P(y∣xk+2)+log⁡P(birth)P(dead)+log⁡q(θ)q(r)1m+1−log⁡∣J2∣\log R=\log \frac{P\left( x_k \right) P\left( y|x_k \right)}{P\left( x_{k+2} \right) P\left( y|x_{k+2} \right)}+\log \frac{P\left( birth \right)}{P\left( dead \right)}+\log \frac{\small{q\left( \theta \right) q\left( r \right)}}{\frac{1}{m+1}}-\log \left| J_2 \right|

其中m其中 m 是当前系统中复共轭极点对的数量 J1J_1 和 J2J_2 是雅可比矩阵行列式 ∣J1∣=∣∏i=1k(λi−ui)∣\left| J_1 \right|=\left| \prod_{\boldsymbol{i}=1}^{\boldsymbol{k}}{\left( \lambda _i-u_i \right)} \right| ∣J2∣=∣∏i=1k(λk+1−λi)(λk+2−λi)∣∣λk+1−λk+1∣\left| J_2 \right|=\left| \prod_{\boldsymbol{i}=1}^{\boldsymbol{k}}{\left( \lambda _{k+1}-\lambda _i \right) \left( \lambda _{k+2}-\lambda _i \right)} \right|\left| \lambda _{k+1}-\lambda _{k+1} \right|

同维度游动

根的类型转换

先确定映射关系,令

a=r1+r22、b=r1−r22r1=a+b、r2=a−ba=\frac{r_1+r_2}{2}\text{、}b=\frac{r_1-r_2}{2} \\ r_1=a+b\text{、}r_2=a-b

再确定雅可比矩阵行列式

JC→R=∣det⁡(∂a∂r1∂a∂r2∂b∂r1∂b∂r2)∣=12JR→C=JC→R−1=2J_{C\rightarrow R}=\left| \det \left( \begin{matrix} \frac{\partial a}{\partial r_1}& \frac{\partial a}{\partial r_2}\\ \frac{\partial b}{\partial r_1}& \frac{\partial b}{\partial r_2}\\ \end{matrix} \right) \right|=\frac{1}{2} \\ J_{R\rightarrow C}={J_{C\rightarrow R}}^{-1}=2

  1. 实根 <-> 复共轭根对 在nrn_r个实根中任选两个实根r1r_1和r2r_2,转换为复共轭根对a±jba\pm jb,其中a=r1+r22a=\frac{r_1+r_2}{2},b=r1−r22b=\frac{r_1-r_2}{2}。

q(x∣x′)=pm×1C(nr,2)q\left( x|x' \right) =p_m\times \frac{1}{C\left( n_r,2 \right)}

此时的接受率写作:Merge:α(x′,x)=min⁡{1,π(x)π(x′)⋅1q(x∣x′)⋅12}Merge\text{:}\alpha \left( x',x \right) =\min \left\{ 1,\frac{\pi \left( x \right)}{\pi \left( x' \right)}\cdot \frac{1}{q\left( x|x' \right)}\cdot \frac{1}{2} \right\}

  1. 复共轭根对 <-> 实根 在ncn_c个复共轭根对中任选一个复共轭根对a±jba\pm jb,转换为两个实根r1=a+br_1=a+b和r2=a−br_2=a-b。

q(x′∣x)=pc×1ncq\left( x'|x \right) =p_c\times \frac{1}{n_c}

此时的接受率写作:Split:α(x,x′)=min⁡{1,π(x′)π(x)⋅1q(x′∣x)⋅2}Split\text{:}\alpha \left( x,x' \right) =\min \left\{ 1,\frac{\pi \left( x' \right)}{\pi \left( x \right)}\cdot \frac{1}{q\left( x'|x \right)}\cdot 2 \right\}

同类型根的调整

  1. 实根调整 选择一个实根rr,通过添加噪声ϵ∼N(0,σ2)\epsilon \sim N\left( 0,\sigma ^2 \right)调整该实根的位置为r′=r+ϵr' = r + \epsilon。

α(x,x′)=min⁡{1,π(x′)π(x)}\alpha \left( x,x' \right) =\min \left\{ 1,\frac{\pi \left( x' \right)}{\pi \left( x \right)} \right\}

  1. 复共轭根对调整 选择一个复共轭根对a±jba\pm jb,通过添加噪声ϵa∼N(0,σa2)\epsilon _a \sim N\left( 0,\sigma _a^2 \right)和ϵb∼N(0,σb2)\epsilon _b \sim N\left( 0,\sigma _b^2 \right)调整该复共轭根对的位置为a′=a+ϵaa' = a + \epsilon _a和b′=b+ϵbb' = b + \epsilon _b。

α(x,x′)=min⁡{1,π(x′)π(x)}\alpha \left( x,x' \right) =\min \left\{ 1,\frac{\pi \left( x' \right)}{\pi \left( x \right)} \right\}