Esta página aún no se ha traducido para esta versión. Puede ver la versión más reciente de esta página en inglés.

Sesgo y variabilidad en el periodograma

Este ejemplo muestra cómo reducir el sesgo y la variabilidad en el periodograma. El uso de una ventana puede reducir el sesgo en el periodograma, y el uso de ventanas con promedio puede reducir la variabilidad.

Utilice procesos autoregresivos estacionarios (AR) de amplio sentido para mostrar los efectos del sesgo y la variabilidad en el periodograma. Los procesos de realidad aumentada presentan un modelo conveniente porque sus PSD tienen expresiones de forma cerrada. Cree un modelo AR(2) del siguiente formulario:

<math display="block">
<mrow>
<mi>y</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo stretchy="false">)</mo>
<mo>-</mo>
<mn>0</mn>
<mo>.</mo>
<mn>7</mn>
<mn>5</mn>
<mi>y</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
<mo stretchy="false">)</mo>
<mo>+</mo>
<mn>0</mn>
<mo>.</mo>
<mn>5</mn>
<mi>y</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>2</mn>
<mo stretchy="false">)</mo>
<mo>=</mo>
<mi>ε</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo stretchy="false">)</mo>
<mo>,</mo>
</mrow>
</math>

Dónde

<math display="block">
<mrow>
<mi>ε</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo stretchy="false">)</mo>
</mrow>
</math>
es una secuencia de ruido blanco media cero con alguna varianza especificada. En este ejemplo, suponga que la varianza y el período de muestreo son 1. Para simular el proceso AR(2) anterior, cree un filtro de todo polo (IIR). Vea la respuesta de magnitud del filtro.

B2 = 1; A2 = [1 -0.75 0.5]; fvtool(B2,A2)

Este proceso es paso de banda. El rango dinámico del PSD es de aproximadamente 14,5 dB, como se puede determinar con el código siguiente.

[H2,W2] = freqz(B2,A2,1e3,1); dr2 = max(20*log10(abs(H2)))-min(20*log10(abs(H2)))
dr2 = 14.4984 

Al examinar la colocación de los polos, verá que este proceso AR(2) es estable. Los dos polos están dentro del círculo de la unidad.

fvtool(B2,A2,'Analysis','polezero')

A continuación, cree un proceso AR(4) descrito por la siguiente ecuación:

<math display="block">
<mrow>
<mi>y</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo stretchy="false">)</mo>
<mo>-</mo>
<mn>2</mn>
<mo>.</mo>
<mn>7</mn>
<mn>6</mn>
<mn>0</mn>
<mn>7</mn>
<mi>y</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
<mo stretchy="false">)</mo>
<mo>+</mo>
<mn>3</mn>
<mo>.</mo>
<mn>8</mn>
<mn>1</mn>
<mn>0</mn>
<mn>6</mn>
<mi>y</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>2</mn>
<mo stretchy="false">)</mo>
<mo>-</mo>
<mn>2</mn>
<mo>.</mo>
<mn>6</mn>
<mn>5</mn>
<mn>3</mn>
<mn>5</mn>
<mi>y</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>3</mn>
<mo stretchy="false">)</mo>
<mo>+</mo>
<mn>0</mn>
<mo>.</mo>
<mn>9</mn>
<mn>2</mn>
<mn>3</mn>
<mn>8</mn>
<mi>y</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>4</mn>
<mo stretchy="false">)</mo>
<mo>=</mo>
<mi>ε</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo stretchy="false">)</mo>
<mo>.</mo>
</mrow>
</math>

Utilice el código siguiente para ver la respuesta de magnitud de este sistema IIR.

B4 = 1; A4 = [1 -2.7607 3.8106 -2.6535 0.9238]; fvtool(B4,A4)

Examinando la colocación de los polos, puede ver que este proceso AR(4) también es estable. Los cuatro polos están dentro del círculo de la unidad.

fvtool(B4,A4,'Analysis','polezero')

El rango dinámico de este PSD es de aproximadamente 65 dB, mucho mayor que el modelo AR(2).

[H4,W4] = freqz(B4,A4,1e3,1); dr4 = max(20*log10(abs(H4)))-min(20*log10(abs(H4)))
dr4 = 64.6213 

Para simular realizaciones de estos procesos AR( ), utilice y .prandnfilter Establezca el generador de números aleatorios en la configuración predeterminada para producir resultados repetibles. Traza las realizaciones.

rng default x = randn(1e3,1); y2 = filter(B2,A2,x); y4 = filter(B4,A4,x);  subplot(2,1,1) plot(y2) title('AR(2) Process') xlabel('Time')  subplot(2,1,2) plot(y4) title('AR(4) Process') xlabel('Time')

Calcular y trazar los periodogramas de las realizaciones AR(2) y AR(4). Compare los resultados con el verdadero PSD. Tenga en cuenta que convierte las frecuencias en milihercios para el trazado.periodograma

Fs = 1; NFFT = length(y2);  subplot(2,1,1) periodogram(y2,rectwin(NFFT),NFFT,Fs) hold on plot(1000*W2,20*log10(abs(H2)),'r','linewidth',2) title('AR(2) PSD and Periodogram')  subplot(2,1,2) periodogram(y4,rectwin(NFFT),NFFT,Fs) hold on plot(1000*W4,20*log10(abs(H4)),'r','linewidth',2) title('AR(4) PSD and Periodogram') text(350,20,'\downarrow Bias')

En el caso del proceso AR(2), la estimación del periodograma sigue la forma del verdadero PSD, pero presenta una variabilidad considerable. Esto se debe a los bajos grados de libertad. Las desviaciones negativas pronunciadas (en dB) en el periodograma se explican tomando el registro de una variable aleatoria chi-cuadrada con dos grados de libertad.

En el caso del proceso AR(4), el periodograma sigue la forma del verdadero PSD a bajas frecuencias, pero se desvía del PSD en las frecuencias altas. Este es el efecto de la convolución con el kernel de Fejer. El gran rango dinámico del proceso AR(4) en comparación con el proceso AR(2) es lo que hace que el sesgo sea más pronunciado.

Mitigue el sesgo demostrado en el proceso AR(4) mediante una cónico o una ventana. En este ejemplo, utilice una ventana Hamming para reducir la realización de AR(4) antes de obtener el periodograma.

figure periodogram(y4,hamming(length(y4)),NFFT,Fs) hold on plot(1000*W4,20*log10(abs(H4)),'r','linewidth',2) title('AR(4) PSD and Periodogram with Hamming Window') legend('Periodogram','AR(4) PSD')

Tenga en cuenta que la estimación del periodograma ahora sigue el verdadero PSD AR(4) en todo el rango de frecuencia Nyquist. Las estimaciones de periodograma sólo tienen dos grados de libertad, por lo que el uso de una ventana no reduce la variabilidad del periodograma, pero sí aborda el sesgo.

En la estimación espectral no paramétrica, dos métodos para aumentar los grados de libertad y reducir la variabilidad del periodograma son el promedio del segmento superpuesto de Welch y la estimación espectral multicómayor.

Obtenga una estimación multitaper de la serie temporal AR(4) utilizando un producto de medio ancho de banda de tiempo de 3.5. Trazar el resultado.

NW = 3.5;  figure pmtm(y4,NW,NFFT,Fs) hold on plot(1000*W4,20*log10(abs(H4)),'r','linewidth',2) legend('Multitaper Estimate','AR(4) PSD')

El método multitaper produce una estimación de PSD con una variabilidad significativamente menor que el periodograma. Dado que el método multitaper también utiliza ventanas, verá que también se aborda el sesgo del periodograma.

Consulte también

|