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.

Regresión lineal de dominio de frecuencia

Este ejemplo muestra cómo utilizar la transformada de Fourier discreta para construir un modelo de regresión lineal para una serie de tiempo. La serie temporal utilizada en este ejemplo es el número mensual de muertes accidentales en los Estados Unidos de 1973 a 1979. Los datos se publican en Brockwell y Davis (2006). La fuente original es el Consejo Nacional de seguridad de los Estados.

Introduzca los datos. Copie la matriz en el espacio de trabajo de MATLAB®.exdata

exdata = [         9007        7750        8162        7717        7792        7836         8106        6981        7306        7461        6957        6892         8928        8038        8124        7776        7726        7791         9137        8422        7870        7925        8106        8129        10017        8714        9387        8634        8890        9115        10826        9512        9556        8945        9299        9434        11317       10120       10093       10078       10625       10484        10744        9823        9620        9179        9302        9827         9713        8743        8285        8037        8314        9110         9938        9129        8433        8488        8850        9070         9161        8710        8160        7874        8265        8633         8927        8680        8034        8647        8796        9240];

es una matriz de 12 por 6.exdata Cada columna contiene 12 meses de datos.exdata La primera fila de cada columna contiene el número de muertes accidentales de Estados Unidos para enero del año correspondiente. La última fila de cada columna contiene el número de muertes accidentales de Estados Unidos para diciembre del año correspondiente.

Remodele la matriz de datos en una serie de tiempo de 72 por 1 y trace los datos para los años 1973 a 1978.

ts = reshape(exdata,72,1); years = linspace(1973,1979,72);  plot(years,ts,'o-','MarkerFaceColor','auto') xlabel('Year') ylabel('Number of Accidental Deaths')

Una inspección visual de los datos indica que el número de muertes accidentales varía de manera periódica. El período de la oscilación parece ser de aproximadamente 1 año (12 meses). La naturaleza periódica de los datos sugiere que un modelo adecuado puede ser

<math display="block">
<mrow>
<mi>X</mi>
<mo stretchy="false">(</mo>
<mi>n</mi>
<mo stretchy="false">)</mo>
<mo>=</mo>
<mi>μ</mi>
<mo>+</mo>
<munder>
<mrow>
<mo></mo>
</mrow>
<mrow>
<mi>k</mi>
</mrow>
</munder>
<mrow>
<mo>(</mo>
<msub>
<mrow>
<mi>A</mi>
</mrow>
<mrow>
<mi>k</mi>
</mrow>
</msub>
<mi mathvariant="normal">cos</mi>
<mrow>
<mfrac>
<mrow>
<mrow>
<mn>2</mn>
<mi>π</mi>
<mi>k</mi>
<mi>n</mi>
</mrow>
</mrow>
<mrow>
<mrow>
<mi>N</mi>
</mrow>
</mrow>
</mfrac>
</mrow>
<mo>+</mo>
<msub>
<mrow>
<mi>B</mi>
</mrow>
<mrow>
<mi>k</mi>
</mrow>
</msub>
<mi mathvariant="normal">cos</mi>
<mrow>
<mfrac>
<mrow>
<mrow>
<mn>2</mn>
<mi>π</mi>
<mi>k</mi>
<mi>n</mi>
</mrow>
</mrow>
<mrow>
<mrow>
<mi>N</mi>
</mrow>
</mrow>
</mfrac>
</mrow>
<mo>)</mo>
</mrow>
<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>
</mrow>
</math>
es la media general,
<math display="block">
<mrow>
<mi>N</mi>
</mrow>
</math>
es la longitud de la serie temporal y
<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 de variables aleatorias Gaussianas independientes y de distribución idéntica (IID) con media cero y cierta varianza. El término de ruido aditivo representa la aleatoriedad inherente a los datos. Los parámetros del modelo son la media general y las amplitudes de las senos cosenos y Sines. El modelo es lineal en los parámetros.

Para construir un modelo de regresión lineal en el dominio de tiempo, debe especificar qué frecuencias utilizar para las senos cosenos y Sines, formar la matriz de diseño y resolver las ecuaciones normales para obtener las estimaciones de mínimos cuadrados de los parámetros del modelo. En este caso, es más fácil utilizar la transformada discreta de Fourier para detectar las periodicidades, retener sólo un subconjunto de los coeficientes de Fourier e invertir la transformación para obtener la serie temporal ajustada.

Realizar un análisis espectral de los datos para revelar qué frecuencias contribuyen significativamente a la variabilidad en los datos. Debido a que la media general de la señal es de aproximadamente 9.000 y es proporcional a la transformada de Fourier a una frecuencia de 0, reste la media antes del análisis espectral. Esto reduce el coeficiente de Fourier de gran magnitud a 0 de frecuencia y hace que las oscilaciones significativas sean más fáciles de detectar. Las frecuencias de la transformada de Fourier se espacian en un intervalo que es el recíproco de la longitud de la serie temporal, 1/72. Muestrear los datos mensualmente, la frecuencia más alta en el análisis espectral es 1 ciclo/2 meses. En este caso, es conveniente mirar el análisis espectral en términos de ciclos/año para escalar las frecuencias en consecuencia para la visualización.

tsdft = fft(ts-mean(ts)); freq = 0:1/72:1/2;  plot(freq.*12,abs(tsdft(1:length(ts)/2+1)),'o-', ...     'MarkerFaceColor','auto') xlabel('Cycles/Year') ylabel('Magnitude') ax = gca; ax.XTick = [1/6 1 2 3 4 5 6];

En base a las magnitudes, la frecuencia de 1 ciclo/12 meses es la oscilación más significativa en los datos. La magnitud en 1 ciclo/12 meses es más del doble de grande que cualquier otra magnitud. Sin embargo, el análisis espectral revela que también hay otros componentes periódicos en los datos. Por ejemplo, parece que hay componentes periódicos en armónicos (múltiplos enteros) de 1 ciclo/12 meses. También parece haber un componente periódico con un período de 1 ciclo/72 meses.

Basado en el análisis espectral de los datos, ajuste un modelo de regresión lineal simple usando un término coseno y seno con una frecuencia del componente más significativo: 1 ciclo/año (1 ciclo/12 meses).

Determine la bandeja de frecuencias en la transformada de Fourier discreta que corresponde a 1 ciclo/12 meses. Debido a que las frecuencias se espacian en 1/72 y la primera bin corresponde a la frecuencia 0, la ubicación correcta es 72/12 + 1. Esta es la ubicación de frecuencia de la frecuencia.positive También debe incluir la bandeja de frecuencias correspondiente a la frecuencia:-1 ciclo/12 meses.negative Con MATLAB Indexing, la ubicación de frecuencia de la frecuencia negativa es 72-72/12 + 1.

Cree un vector 72-by-1 de ceros. Rellene los elementos apropiados del vector con los coeficientes de Fourier correspondientes a una frecuencia positiva y negativa de 1 ciclo/12 meses. Invierta la transformada de Fourier y añada la media general para obtener un ajuste a los datos de muerte accidentales.

freqbin = 72/12; freqbins = [freqbin 72-freqbin]+1; tsfit = zeros(72,1); tsfit(freqbins) = tsdft(freqbins); tsfit = ifft(tsfit); mu = mean(ts); tsfit = mu+tsfit;

Trace los datos originales junto con la serie ajustada utilizando dos coeficientes de Fourier.

plot(years,ts,'o-','MarkerFaceColor','auto') xlabel('Year') ylabel('Number of Accidental Deaths') hold on plot(years,tsfit,'linewidth',2) legend('Data','Fitted Model') hold off

El modelo ajustado parece captar la naturaleza periódica general de los datos y apoya la conclusión inicial de que los datos oscilan con un ciclo de 1 año.

Para evaluar cuán adecuadamente la frecuencia única de 1 ciclo/12 meses representa la serie temporal observada, forme los residuos. Si los residuos se asemejan a una secuencia de ruido blanco, el modelo lineal simple con una frecuencia ha modelado adecuadamente la serie temporal.

Para evaluar los residuos, utilice la secuencia de autocorrelación con intervalos de confianza de 95% para un ruido blanco.

resid = ts-tsfit; [xc,lags] = xcorr(resid,50,'coeff');  stem(lags(51:end),xc(51:end),'filled') hold on lconf = -1.96*ones(51,1)/sqrt(72); uconf = 1.96*ones(51,1)/sqrt(72); plot(lags(51:end),lconf,'r') plot(lags(51:end),uconf,'r') xlabel('Lag') ylabel('Correlation Coefficient') title('Autocorrelation of Residuals') hold off

Los valores de autocorrelación caen fuera de los límites de confianza del 95% en un número de retrasos. No parece que los residuos sean ruido blanco. La conclusión es que el modelo lineal simple con un componente sinusoidal no tiene en cuenta todas las oscilaciones en el número de muertes accidentales. Esto se espera porque el análisis espectral reveló componentes periódicos adicionales además de la oscilación dominante. La creación de un modelo que incorpore términos periódicos adicionales indicados por el análisis espectral mejorará el ajuste y blanqueará los residuos.

Ajuste un modelo que consista en las tres magnitudes de coeficiente de Fourier más grandes. Debido a que tiene que retener los coeficientes de Fourier correspondientes a las frecuencias negativas y positivas, conserve los 6 índices más grandes.

tsfit2dft = zeros(72,1); [Y,I] = sort(abs(tsdft),'descend'); indices = I(1:6); tsfit2dft(indices) = tsdft(indices);

Demostrar que preservar sólo 6 de los 72 coeficientes de Fourier (3 frecuencias) retiene la mayor parte de la energía de la señal. En primer lugar, demostrar que retener todos los coeficientes de Fourier produce equivalencia de energía entre la señal original y la transformada de Fourier.

norm(1/sqrt(72)*tsdft,2)/norm(ts-mean(ts),2)
ans = 1.0000 

La relación es 1. Ahora, examine la relación de energía donde sólo 3 frecuencias se conservan.

norm(1/sqrt(72)*tsfit2dft,2)/norm(ts-mean(ts),2)
ans = 0.8991 

Casi el 90% de la energía se retiene. De manera equivalente, el 90% de la varianza de la serie temporal se contabiliza en 3 componentes de frecuencia.

Forme una estimación de los datos basados en 3 componentes de frecuencia. Compare los datos originales, el modelo con una frecuencia y el modelo con 3 frecuencias.

tsfit2 = mu+ifft(tsfit2dft,'symmetric');  plot(years,ts,'o-','markerfacecolor','auto') xlabel('Year') ylabel('Number of Accidental Deaths') hold on plot(years,tsfit,'linewidth',2) plot(years,tsfit2,'linewidth',2) legend('Data','1 Frequency','3 Frequencies') hold off

El uso de 3 frecuencias ha mejorado el ajuste a la señal original. Puede ver esto examinando la autocorrelación de los residuos del modelo de 3 frecuencias.

resid = ts-tsfit2; [xc,lags] = xcorr(resid,50,'coeff');  stem(lags(51:end),xc(51:end),'filled') hold on lconf = -1.96*ones(51,1)/sqrt(72); uconf = 1.96*ones(51,1)/sqrt(72); plot(lags(51:end),lconf,'r') plot(lags(51:end),uconf,'r') xlabel('Lag') ylabel('Correlation Coefficient') title('Autocorrelation of Residuals') hold off

El uso de 3 frecuencias ha resultado en residuos que se aproximan más estrechamente a un proceso de ruido blanco.

Demuestre que los valores de parámetro obtenidos de la transformada de Fourier equivalen a un modelo de regresión lineal de dominio de tiempo. Encuentre las estimaciones de mínimos cuadrados para la media general, las amplitudes del coseno y las amplitudes sinusoidales para las tres frecuencias formando la matriz de diseño y resolviendo las ecuaciones normales. Compare la serie de tiempo ajustada con la obtenida de la transformada de Fourier.

X = ones(72,7); X(:,2) = cos(2*pi/72*(0:71))'; X(:,3) = sin(2*pi/72*(0:71))'; X(:,4) = cos(2*pi*6/72*(0:71))'; X(:,5) = sin(2*pi*6/72*(0:71))'; X(:,6) = cos(2*pi*12/72*(0:71))'; X(:,7) = sin(2*pi*12/72*(0:71))'; beta = X\ts; tsfit_lm = X*beta; max(abs(tsfit_lm-tsfit2))
ans = 7.2760e-12 

Los dos métodos producen resultados idénticos. El valor absoluto máximo de la diferencia entre las dos formas de onda es del orden de 10-12. En este caso, el enfoque de dominio de frecuencia era más fácil que el enfoque de dominio de tiempo equivalente. Naturalmente, se utiliza un análisis espectral para inspeccionar visualmente las oscilaciones que están presentes en los datos. A partir de ese paso, es sencillo utilizar los coeficientes de Fourier para construir un modelo para la señal consistente en una suma de cosenos y senos.

Para obtener más detalles sobre el análisis espectral en series temporales y la equivalencia con la regresión del dominio del tiempo, vea (Shumway y Stoffer, 2006).

Mientras que el análisis espectral puede responder a qué componentes periódicos contribuyen significativamente a la variabilidad de los datos, no explica por qué esos componentes están presentes. Si examina estos datos de cerca, verá que los valores mínimos en el ciclo de 12 meses tienden a producirse en febrero, mientras que los valores máximos se producen en julio. Una explicación plausible de estos datos es que las personas son naturalmente más activas en verano que en invierno. Desafortunadamente, como resultado de este aumento de la actividad, hay una mayor probabilidad de la ocurrencia de accidentes mortales.

Referencias

Brockwell, Peter J., y Richard A. Davis. .Time Series: Theory and Methods Nueva York: Springer, 2006.

Shumway, Robert H., y David S. Stoffer. .Time Series Analysis and Its Applications with R Examples Nueva York: Springer, 2006.

Consulte también

| |