信息处理 - 有人可以详细说明如何将正弦模型的 AIC 应用于特定数据吗？ - 吾爱随笔录

有人可以详细说明如何将正弦模型的 AIC 应用于特定数据吗？

信息处理自回归模型

2022-01-25 16:53:47

注意： 这是另一个用户一直试图（未成功）提出的问题。因为从本质上来说，问的多个问题要么被我删除（因为它们彼此重复），要么被用户删除，所以我试图以对用户有所帮助的方式提出问题，并且是非常适合该网站。

假设我有一个数据测量向量， $x(n)$ ，我想为该数据拟合一个模型：

x (n) = \sum_{p = 1}^{P} a_{p} \cos (2 π f_{p} + ϕ_{p}) + ϵ (n)

$x(n) = \sum_{p=1}^{P} a_p \cos(2\pi f_p + \phi_p) + \epsilon(n)$ 其中每个

P

$P$ sinusoids 参数化为

(a_{p}, f_{p}, ϕ_{p})

$(a_p, f_p, \phi_p)$ 和

ϵ

$\epsilon$ 是具有零均值和方差的 iid 高斯白噪声过程

σ^{2}

$\sigma^2$ .

我的问题是：我们如何使用AIC（或BIC或MDL）来估计模型顺序， $P$ ?

维基百科条目只是说：

A I C = 2 P - 2 \ln (L)

$AIC = 2P - 2\ln(L)$ 在哪里

L

$L$ 是似然函数。

我如何计算 $L$ 给定 $x(n)$ ?

1个回答

多色调估计问题的 PDF（似然函数）为：

L (x; a, f, ϕ, σ^{2}) = \frac{1}{(2 π σ^{2})^{\frac{N}{2}}} \exp (- \frac{1}{2 σ^{2}} \sum_{n = 0}^{N - 1} (x (n) - \sum_{p = 1}^{P} a_{p} \cos (2 π f_{p} + ϕ_{p}))^{2})

${\cal L}({\bf x}; {\bf a}, {\bf f},{\bf \phi}, {\sigma}^2) = \frac{1}{(2\pi \sigma^2)^{\frac{N}{2}}} \exp\left( - \frac{1}{2\sigma^2} \sum_{n=0}^{N-1} ( x(n) - \sum_{p=1}^{P} {a}_p \cos(2\pi {f}_p + {\phi}_p) )^2 \right)$

Djuric 的论文将 AIC 和 MDL 值引用为：

在此处输入图像描述

最后的代码尝试：

生成 5 音噪声信号。
对给定值的幅度、频率和相位进行 ML 估计 $p$ .
找到每个值的（对数）似然值 $p$ .
从 Djuric 计算 AIC 和 MDL 值。

情节为 $p = 1 \ldots 60$ 是：

在此处输入图像描述

通常，AIC会高估模型阶数。在这种特殊情况下，AIC 得到 7，MDL 得到 5（正确值）。绘制每条曲线的最小值（如oAIC 曲线+上的 a 和 MDL 曲线上的 a）。

Djuric, PM，“高斯白噪声中正弦曲线的模型选择规则”，信号处理，IEEE Transactions on（第 44 卷，第 7 期），第 1744 - 1751 页。

scilab代码仅在下面

P = 5;
f = [ 0.1 0.12 0.231 0.333 0.478634];
phi = [0.329847 5.90324 3.09834 4.907983 1.984];
a = [ 1 1 2 2 3];
sigma = 1;

N = 100;

x = cos(2*%pi*[0:N-1]'*f + ones(N,1)* phi )*a';

xn = x + sigma*rand(N,1,'normal');


function [L,LL] = likelihood(data,ahat,fhat,phihat,sigmahat)
    N = length(data);
    datahat = cos(2*%pi*[0:N-1]'*fhat + ones(N,1)* phihat )*ahat';
    L = 1/(2*%pi*sigmahat*sigmahat)^(N/2) * exp(-1/(2*sigmahat*sigmahat)*sum((data - datahat).^2));
    LL = - log(2*%pi*sigmahat*sigmahat)*N/2 - -1/(2*sigmahat*sigmahat)*sum((data - datahat).^2);
endfunction

function [a,m] = aic(p,L,N)
    a = 3*p - log(L);
    m = 3*p/2*log(N) - log(L);
endfunction

function [ahat,fhat,phihat,regr] = sin_estimate(data)
    PAD = 10;
    N = length(data);
    DATA = fft([data; zeros((PAD-1)*N,1)],-1);
    [mx,ix] = max(abs(DATA(1:N*PAD/2,1)));
    ahat = mx/(N/2);
    fhat = (ix-1)/N/PAD;
    phihat = atan(imag(DATA(ix)),real(DATA(ix)));

    regr = data - cos(2*%pi*[0:N-1]'*fhat + ones(N,1)* phihat )*ahat
endfunction

pValues = [1:60];
likeValues = [];
loglikeValues = [];
aicValues = [];
mdlValues = [];
for p=pValues,
    signal = xn;
    a_est = zeros(1,p);
    f_est = zeros(1,p);
    phi_est = zeros(1,p);
    for pp=1:p
        [a_est(pp),f_est(pp),phi_est(pp),signal] =  sin_estimate(signal);   
    end

    [l,ll] = likelihood(xn,a_est,f_est,phi_est,sigma);

    likeValues = [likeValues l];
    loglikeValues = [loglikeValues ll];

    [aicv,mdlv] = aic(p,likeValues(p));
    aicValues = [aicValues; aicv]
    mdlValues = [mdlValues; mdlv]
end

clf;
plot(aicValues);
plot(mdlValues,'r');
title('AIC and MDL plots');

其它你可能感兴趣的问题

上一篇基于 ACF 结果求基频的算法下一篇当 FFT 大小小于信号持续时间时，幅度如何缩放？