使用Eigen和Rcpp的Kalman滤波器和平滑器的执行时间

问题描述

我有一个用c ++编写的多元时间序列的卡尔曼滤波器和干扰平滑器的实现。过滤器计算似然并存储向量K，标量v和s进行平滑：

for(int i = 0; i < N; i++) {
    m = At*m;
    P = At*P*At.transpose() + Q;
    
    for(int j = 0; j < D; j++) {
      if(Y(i,j)>0) {
        v = Y(i,j) - mu(j) - H.row(j)*m;
        K = P*H.row(j).transpose();
        s = H.row(j)*K + eta(j);
        m += K*v/s;
        P = P.selfadjointView<Lower>().rankUpdate(K,-1/s);
        ll += -0.5*(log(s)+pow(v,2)/s);
        K_arr.col(i*D+j) = K;
        s_arr(i,j) = s;
        v_arr(i,j) = v;
      }
    }
  }

平滑器是类似的双循环，其中L和Nm是与滤波器中状态协方差矩阵P大小相同的矩阵。计算观测噪声向量的梯度：

for(int i = N-1; i > -1; i--) {
  for(int j = D-1; j > -1; j--) {
    if(s_arr(i,j)>0) { 
      s = s_arr(i,j); v = v_arr(i,j); K = K_arr.col(i*D+j)/s;
      L = I - K*H.row(j);
      Nm = H.row(j).transpose()*H.row(j)/s + L.transpose()*Nm*L;
      r = H.row(j).transpose()*v/s + L.transpose()*r;
      ehat = eta(j)/s * (v - K.dot(r));
      var_ehat = pow( (eta(j)/s),2) * (s + K.dot(Nm * K));
      g_eta(j) += log(eta(j)) + (pow(ehat,2) + var_ehat)/eta(j);
    }
  }
  r = At.transpose()*r;
  Nm = At.transpose()*Nm*At;

}

让我感到困扰的是，当我从R调用Rcpp编译函数时，添加平滑器会使执行时间与仅过滤相比增加大约十倍。有什么想法导致如此长的运动时间增加吗？

解决方法

暂无找到可以解决该程序问题的有效方法，小编努力寻找整理中！

如果你已经找到好的解决方法，欢迎将解决方案带上本链接一起发送给小编。

小编邮箱:dio#foxmail.com (将#修改为@）