1. 理解
1.1 Matlab 帮助:
a = arburg(x,p)返回与输入数组x的p阶模型相对应的归一化自回归(AR)参数。 如果x是一个向量,则输出数组a是一个行向量。 如果x是矩阵,则参数沿模型的第n行位于x的第n列。 a有p + 1列。 p必须小于x的元素(或行)数。
[a,e] = arburg(x,p)返回白噪声输入的估计方差e。
[a,e,rc] = arburg(x,p)返回rc中的反射系数。
1.2 自我的一些理解
AR(P) 模型作用: 使用前几个数据来预测后面的数据,p是阶数,对应p个历史数据
In an AR model of order p, the current output is a linear combination of the past p outputs plus a white noise input.
The weights on the p past outputs minimize the mean-square prediction error of the autoregression.
If y(n) is the current value of the output and x(n) is a zero mean white noise input, the AR(p) model is:
在阶数为p的AR模型中,当前输出是过去p输出的线性组合加上白噪声
p过去的输出上的权重使自回归的均方预测误差最小。
如果y(n)是输出的当前值并且x(n)是零平均白噪声输入,则AR(p)模型为:
y(n) = -(负号) p个累加a(k)*y(n-k) + 0均值的高斯白噪声x(n) [通常其他表达里使用 w(n)]
-----------------------------------
[ ar_coeffs,NoiseVariance] = arburg(data,order)
ar_coeffs 第一位不需要用,是用于归一化,后面是对应的系数
NoiseVariance 是方差
举例:
p=3: [a E]=arburg(x,3)
a =1.0000 -0.6982 -0.2626 0.0739
E =0.4567
p=12: [a E]=arburg(x,12)
a =1.0000 -0.6495 -0.3066 -0.0934 0.0987 0.4076 -0.1786 -0.0126 -0.0805 -0.0899 0.0382 0.1628 -0.2501
E =0.3237
1.3 网络
https://blog.csdn.net/weixin_43165881/article/details/106878784
1.4 python
Python 中 librosa 库中有 lpc 函数是使用的 burg 方法。
2使用:
向前向后预测
寻找一个案列,进行AR拟合,与原函数一致
这个AR(1)模型在当初陀螺仪卡尔曼滤波中用过,w(n) 是高斯白噪声,先要通过一组数据算出方差。
burg的具体的应用方法还需要再组织一下想法
2021-03-25 日更新
恰巧这日用python 实现了下, python种的librosa 库种lpc 函数是使用 burg法的线性预测 ,代码如下
# %% 0_0.import libs import numpy as np import librosa as lb import matplotlib.pyplot as plt #%% 0_1.def funcs def next_predict(sig,M_nexts,p): a= lb.lpc(sig,p) # 计算模型系数 len_sig = len(sig) len_sig_pred = M_nexts+ sig sig_predict = np.zeros(len_sig_pred) sig_predict[:len_sig]=sig for i in range(M_nexts): sig_predict[ len_sig + i ] = -np.dot(a[1:],sig_predict[len_sig + i :len_sig+ i-p:-1]) # + error sig_predict_out = sig_predict[len_sig:] return sig_predict_out def forward_predict(sig,M_forward,p): a= lb.lpc(sig,p) # 计算模型系数 len_sig = len(sig) len_sig_pred = M_forward + len_sig sig_predict = np.zeros(len_sig_pred) sig_predict[M_forward:]=sig for i in range(M_forward): sig_predict[ M_forward -1 - i ] = - np.dot(a[1:],sig_predict[ M_forward-i : M_forward-i+ p]) # + error sig_predict_out = sig_predict[:M_forward]
return sig_predict_out
# %% 1. 数据的ar-burg 预测系数 t = np.arange(0,10,0.01) sig1 = np.sin(2*np.pi*t) # sin 函数 sig2 = 0.2 *np.sin(10*np.pi*t) noise = np.random.randn(len(sig1)) sig = sig1 + sig2 + noise # %% 2. 求arburg 预测的数据 # ar,err = arburg(sig,n_Order) # librosa.lpc 使用了burg法计算 lpc系数 p=300 a= lb.lpc(sig,p) # %% 根据得到的系数进行向前向后预测 M_nexts = 200 # 向后预测的个数 M_befores = 200 # 向后预测的个数 # 设第 n个数的预测为,预测M个数 # sig_predict(n)= - a[1]*sig(n-1) -a[2]*sig(n-2) - a[3]*sig(n-3)... -a[p]*sig(n-p) len_sig = len(sig) len_sig_pred = M_nexts+ len_sig sig_predict = np.zeros(len_sig_pred) sig_predict[:len(sig)]=sig # 预测 for i in range(M_nexts): sig_predict[ len_sig + i ] = -np.dot(a[1:],sig_predict[len_sig + i :len_sig+ i-p:-1]) # + error # 如果我们知道err,计算出err的方差,也可以进行预测 plt.figure() plt.plot(sig_predict) plt.plot(sig ) sig_f_pred = forward_predict(sig,M_befores,p) sig_all = np.hstack([sig_f_pred,sig]) plt.figure() plt.plot(sig_all) plt.plot(sig_f_pred)
图1.向前预测
图2. 向后预测