Jakie jest najlepsze przybliżenie IIR (filtr AR) pierwszego rzędu do filtra średniej ruchomej (filtr FIR)?

Załóżmy następujący filtr IIR pierwszego rzędu:

$$y[n] = \alpha x[n] + (1 - \alpha) y[n - 1]$$

Jak mogę wybrać parametr IIR przybliża tak dobrze, jak to możliwe, FIR, który jest średnią arytmetyczną ostatnich próbek:$\alpha$$k$

$$z[n] = \frac{1}{k}x[n] + \frac{1}{k}x[n-1] + \ldots + \frac{1}{k}x[n-k+1]$$

Gdzie , co oznacza, że ​​dane wejściowe dla IIR mogą być dłuższe niż a mimo to chciałbym uzyskać najlepsze przybliżenie średniej z ostatnich danych wejściowych.$n \in [k, \infty)$$k$$k$

Wiem, że IIR ma nieskończoną odpowiedź impulsową, dlatego szukam najlepszego przybliżenia. Byłbym szczęśliwy dla rozwiązania analitycznego, czy to dla funkcji kosztu czy .L ${L}_{2}$${L}_{1}$

Jak można rozwiązać problemy z optymalizacją, biorąc pod uwagę tylko IIR pierwszego rzędu.

Dzięki.

Nie ma rozwiązania analitycznego dla będącego skalarem (tak myślę). Oto skrypt, który daje a dla danego K . Jeśli potrzebujesz go online, możesz zbudować LUT. Skrypt znajduje rozwiązanie, które minimalizuje$\alpha$$\alpha$$K$

$$\int_{0}^{\pi} dw \quad \left|H_1(jw) - H_2(jw)\right|^2$$

gdzie to odpowiedź częstotliwościowa FIR, a H 2 to odpowiedź częstotliwościowa IIR.$H_1$$H_2$

Nie określiłeś żadnego zakresu dla K. Ale chcę tylko wyjaśnić, że następujący system jest równoważny z twoim średnim filtrem i ma tę samą złożoność obliczeniową i twoje pierwsze zamówienie IIR!

$H(z) = \frac{1}{K} \frac{1 - z^{-K}}{1-z^{-1}}$

function a = find_a(K)

w = 0.0001:0.001:pi;
as = [-1:0.001:-0.001  0.001:0.001:1];

E = zeros(size(as));
for idx=1:length(as)
    fJ = J(w,as(idx),K);
    E(idx) = sum(fJ);
end

[Emin, indx] = min(E)
a = as(indx)

function f = J(w,a,K)
    num = 2*(2-a)*(1-cos(w*K)) + 2*(cos(w*(K-1)) - cos(w)) - 2*(1-a)*(cos(w)-cos(w*(K+1)));
    den = (2-a)^2 + 1 + (1-a)^2 + 2*(1-a)*cos(2*w) - 2*(2-a)^2*cos(w);
    f = -(a/K)*num./den;
    f = f+(1/K^2)*(1-cos(w*K))./(1-cos(w))+a^2./(1+(1-a)^2-2*(1-a)*cos(w));
end

end

Miło dyskutuje się na ten temat w przetwarzaniu sygnałów osadzonych za pomocą architektury mikro sygnałów , mniej więcej między stronami 63 i 69 . Na stronie 63 zawiera wyprowadzenie dokładnego filtru rekurencyjnej średniej ruchomej (który niaren podał w swojej odpowiedzi ),

$$H(z) = { 1 \over{N} } { 1 - z^{-N} \over { 1 - z^{-1} } }.$$

Dla wygody w odniesieniu do poniższej dyskusji odpowiada to równaniu różnicy:

$$y_n = y_{n-1} + {1\over{N} }(x_n - x_{n-N}).$$

Przybliżenie, które wprowadza filtr do określonej przez ciebie postaci, wymaga przyjęcia, że , ponieważ (i cytuję z str. 68 ) „ y n - 1 jest średnią z x n próbek”. To przybliżenie pozwala nam uprościć poprzednie równanie różnicy w następujący sposób:$x_{n-N} \approx y_{n-1}$$y_{n-1}$$x_n$

$$\begin{array}\\ y_{n} = y_{n-1} + {1\over{N}}(x_n - y_{n-1}) \\ y_{n} = y_{n-1} - {1\over{N}} y_{n-1} + {1\over{N}}x_n \\ y_{n} = (1-{1\over{N}})y_{n-1} + {1\over{N}}x_n. \end{array}$$

Ustawienie , dochodzimy do twojej pierwotnej postaci,yn=αxn+(1-α)yn-1, co pokazuje, że pożądany współczynnik (w odniesieniu do tego przybliżenia) wynosi dokładnie$\alpha = {1\over{N}}$$y_{n} = \alpha x_n + (1-\alpha)y_{n-1}$ (gdzieNjest liczbą próbek). $1\over{N}$$N$

Czy to przybliżenie jest w pewnym sensie „najlepsze”? Z pewnością jest elegancki. Oto jak porównuje się odpowiedź amplitudowa [przy 44,1 kHz] dla N = 3, a gdy N wzrasta do 10 (przybliżenie na niebiesko):

---

Jak sugeruje odpowiedź Petera , zbliżenie filtra FIR za pomocą filtra rekurencyjnego może być problematyczne pod normą dotyczącą najmniejszych kwadratów. Obszerną dyskusję na temat ogólnego rozwiązania tego problemu można znaleźć w pracy dyplomowej JOS, Techniki projektowania filtrów cyfrowych i identyfikacji systemu za pomocą aplikacji na skrzypce . Opowiada się za użyciem Normy Hankela, ale w przypadkach, w których reakcja fazowa nie ma znaczenia, obejmuje on także Metodę Kopca, która może w tym przypadku działać dobrze (i stosuje normę ). Szeroki przegląd technik zawartych w pracy można znaleźć tutaj . Mogą dać inne interesujące przybliżenia.$L^2$

OK, spróbujmy wyprowadzić to, co najlepsze: tak, że współczynnik zx[n-m]toα(1-

$$\begin{array}{lcl} y[n] &=& \alpha x[n] + (1 - \alpha) y[n - 1] \\ &=& \alpha x[n] + (1 - \alpha) \alpha x[n-1] + (1 - \alpha)^2 y[n - 2]\\ &=& \alpha x[n] + (1 - \alpha) \alpha x[n-1] + (1 - \alpha)^2 \alpha x[n-2] + (1 - \alpha)^3 y[n - 3]\\ \end{array}$$ $x[n-m]$ .$\alpha(1-\alpha)^m$

Najlepsze przybliżenie średnich kwadratów zminimalizuje:

$$\begin{array}{lcl} J(\alpha) &=& \sum_{m=0}^{k-1} (\alpha(1-\alpha)^m - \frac{1}{k})^2 + \sum_{m=k}^\infty \alpha^2(1-\alpha)^{2m}\\ &=& \sum_{m=0}^{k-1} \left(\alpha^2(1-\alpha)^{2m} - \frac{2}{k}\alpha(1-\alpha)^m + \frac{1}{k^2}\right) + \alpha^2 (1-\alpha)^{2k} \sum_{m=0}^\infty (1-\alpha)^{2m} \\ &=& \alpha^2\frac{1- (1-\alpha)^{2k}}{1 - (1-\alpha)^2} + \frac{2\alpha}{k} \frac{1 - (1 - \alpha)^k}{1 - (1 - \alpha)} + \frac{\alpha^2(1-\alpha)^{2k}}{1 - (1 - \alpha)^2}+ \frac{1}{k}\\ &=& \frac{\alpha^2}{1 - (1 - \alpha)^2} + \frac{2}{k} (1-(1-\alpha)^k) + \frac{1}{k}\\ &=& \frac{\alpha^2}{2\alpha - \alpha^2 }+ \frac{2}{k} (1-(1-\alpha)^k)+ \frac{1}{k}\\ &=& \frac{\alpha}{2 - \alpha }+ \frac{2}{k} (1-(1-\alpha)^k)+ \frac{1}{k}\\ \end{array}$$ because the FIR coefficients are zero for $m > k - 1$.

Next step is to take derivatives and equate to zero.

---

Looking at a plot of the derived $J$ for $K = 1000$ and $\alpha$ from 0 to 1, it looks like the problem (as I've set it up) is ill-posed, because the best answer is $\alpha = 0$.

---

I think there's a mistake here. The way it should be according to my calculations is:

$$\begin{array}{lcl} J(\alpha) &=& \sum_{m=0}^{k-1} (\alpha(1-\alpha)^m - \frac{1}{k})^2 + \sum_{m=k}^\infty \alpha^2(1-\alpha)^{2m} \\ &=& \sum_{m=0}^{k-1} \left(\alpha^2(1-\alpha)^{2m} - \frac{2}{k}\alpha(1-\alpha)^m + \frac{1}{k^2}\right) + \alpha^2 (1-\alpha)^{2k} \sum_{m=0}^\infty (1-\alpha)^{2m} \\ &=& \alpha^2\frac{1- (1-\alpha)^{2k}}{1 - (1-\alpha)^2} - \frac{2\alpha}{k} \frac{1 - (1 - \alpha)^k}{1 - (1 - \alpha)} + \frac{1}{k} + \frac{\alpha^2(1-\alpha)^{2k}}{1 - (1 - \alpha)^2} \end{array}$$

Simplifying it according to Mathematica yields:

$$J(\alpha) = \frac{\alpha}{2 - \alpha} + \frac{2 {(1 - \alpha)}^{k} -1}{k}$$

Using the following code on MATLAB yields something equivalent though different:

syms a k;

expr1 = (a ^ 2) * ((1 - ((1 - a) ^ (2 * k))) / (1 - ((1 - a) ^ 2)));
expr2 = ((2 * a) / k) * ((1 - ((1 - a) ^ (k))) / (1 - (1 - a)));
expr3 = (1 / k);
expr4 = ((a ^ 2) * ((1 - a) ^ (2 * k))) / (1 - ((1 - a) ^ (2)));

simpExpr = simplify(expr1 - expr2 + expr3 + expr4);

$$J(\alpha) = \frac{-2}{\alpha - 2} - \frac{k - 2{(1 - \alpha)}^{k} + 1}{k}$$

Anyhow, those functions do have minimum.

---

So let's assume that we really only care about the approximation over the support (length) of the FIR filter. In that case, the optimization problem is just:

$$J_2(\alpha) = \sum_{m=0}^{k-1} (\alpha(1-\alpha)^m - \frac{1}{k})^2$$

Plotting $J_2(\alpha)$ for various values of $K$ versus $\alpha$ results in the date in the plots and table below.

For $K$ = 8. $\alpha_{\tt min}$ = 0.1533333
For $K$ = 16. $\alpha_{\tt min}$ = 0.08
For $K$ = 24. $\alpha_{\tt min}$ = 0.0533333
For $K$ = 32. $\alpha_{\tt min}$ = 0.04
For $K$ = 40. $\alpha_{\tt min}$ = 0.0333333
For $K$ = 48. $\alpha_{\tt min}$ = 0.0266667
For $K$ = 56. $\alpha_{\tt min}$ = 0.0233333
For $K$ = 64. $\alpha_{\tt min}$ = 0.02
For $K$ = 72. $\alpha_{\tt min}$ = 0.0166667

The red dashed lines are $1/K$ and the green lines are $\alpha_{\tt min}$, the value of $\alpha$ that minimizes $J_2(\alpha)$ (chosen from $\tt alpha = [0:.01:1]/3;$).

Based in experimental tests with k in range (2 to 100) the best fit (sum squared error) gives a relation of alfa = 1/k^0.865 being k number of samples for MovAvg filter

I stumbled upon this old question and I would like to share my solution. As mentioned in other answers, there is no analytical solution, but the function to be minimized behaves nicely and the optimal value of $\alpha$ can be found easily with a few Newton iterations. There is also a formula to check the optimality of the result.

The impulse response of the length $N$ FIR moving average filter is given by

$$h_{FIR}[n]=\frac{1}{N}(u[n]-u[n-N])\tag{1}$$

where $u[n]$ is the unit step function. The first order IIR filter

$$y[n]=\alpha x[n]+(1-\alpha)y[n-1]\tag{2}$$

has the impulse response

$$h_{IIR}[n]=\alpha(1-\alpha)^nu[n]\tag{3}$$

The goal is now to minimize the squared error

$$\epsilon=\sum_{n=0}^{\infty}\left(h_{FIR}[n]-h_{IIR}[n]\right)^2\tag{4}$$

Using $(1)$ and $(3)$, the error can be written as

$$\begin{align}\epsilon(\alpha)&=\sum_{n=0}^{N-1}\left(\alpha(1-\alpha)^n-\frac{1}{N}\right)^2+\sum_{n=N}^{\infty}\alpha^2(1-\alpha)^{2n}\\&=\alpha^2\sum_{n=0}^{\infty}(1-\alpha)^{2n}-\frac{2\alpha}{N}\sum_{n=0}^{N-1}(1-\alpha)^n+\sum_{n=0}^{N-1}\frac{1}{N^2}\\&=\frac{\alpha^2}{1-(1-\alpha)^2}-\frac{2\alpha}{N}\frac{1-(1-\alpha)^N}{1-(1-\alpha)}+\frac{1}{N}\\&=\frac{\alpha}{2-\alpha}-\frac{2}{N}\left(1-(1-\alpha)^N\right)+\frac{1}{N},\qquad 0<\alpha<2\tag{5}\end{align}$$

This expression is very similar to the one given in this answer, but it's not identical. The restriction on $\alpha$ in $(5)$ makes sure that the infinite sum converges, and it is identical to the stability condition for the IIR filter given by $(2)$.

Setting the derivative of $(5)$ to zero results in

$$(1-\alpha)^{N-1}(2-\alpha)^2=1\tag{6}$$

Note that the optimal $\alpha$ must be in the interval $(0,1]$ because larger values of $\alpha$ result in an alternating impulse response $(3)$, which cannot approximate the constant impulse repsonse of the FIR moving average filter.

Taking the square root of $(6)$ and introducing $\beta=1-\alpha$, we obtain

$$\beta^{(N+1)/2}+\beta^{(N-1)/2}-1=0\tag{7}$$

This equation cannot be solved analytically for $\beta$, but it can be solved for $N$:

$$N=-2\frac{\log(1+\beta)}{\log(\beta)},\qquad \beta\neq 0\tag{8}$$

Equation $(8)$ can be used to double-check a numerical solution of $(7)$; it must return the specified value of $N$.

Equation $(7)$ can be solved with a few lines of (Matlab/Octave) code:

N = 50;     % desired filter length of FIR moving average filter

if ( N == 1 )    % no iteration for trivial case
    b = 0;
else
    % Newton iteration
    b = 1;       % starting value
    Nit = 7;
    n = (N+1)/2;
    for k = 1:Nit,
        f = b^n + b^(n-1) -1;
        fp = n*b^(n-1) + (n-1)*b^(n-2);
        b = b - f/fp;
    end

    % check result
    N0 = -2*log(1+b)/log(b) + 1     % must equal N
end

a = 1 - b;

Below is a table with the optimal values of $\alpha$ for a range of filter lengths $N$:

   N     alpha

   1   1.0000e+00
   2   5.3443e-01
   3   3.8197e-01
   4   2.9839e-01
   5   2.4512e-01
   6   2.0809e-01
   7   1.8083e-01
   8   1.5990e-01
   9   1.4333e-01
  10   1.2987e-01
  20   6.7023e-02
  30   4.5175e-02
  40   3.4071e-02
  50   2.7349e-02
  60   2.2842e-02
  70   1.9611e-02
  80   1.7180e-02
  90   1.5286e-02
 100   1.3768e-02
 200   6.9076e-03 
 300   4.6103e-03
 400   3.4597e-03
 500   2.7688e-03
 600   2.3078e-03
 700   1.9785e-03
 800   1.7314e-03
 900   1.5391e-03
1000   1.3853e-03