1. 說明

?用傅立葉變換預(yù)測(cè)時(shí)序數(shù)據(jù)吩谦，原理是把時(shí)域數(shù)據(jù)轉(zhuǎn)換到頻域静袖，再轉(zhuǎn)換回來．python的numpy和scipy里面都有現(xiàn)成的轉(zhuǎn)換工具fft()和ifft()，但使用時(shí)會(huì)遇到一個(gè)問題：比如25天的數(shù)據(jù)轉(zhuǎn)到頻域再轉(zhuǎn)回時(shí)域，還是25天徐绑，雖然擬合了數(shù)據(jù)，但沒法直接預(yù)測(cè)未來莫辨，本篇介紹用它實(shí)現(xiàn)預(yù)測(cè)的方法．

2. 傅立葉變換

(1) 相關(guān)知識(shí)

?之前寫過關(guān)于傅立葉變換原理的文檔傲茄，這次就不再重復(fù)了，具體請(qǐng)見：http://www.reibang.com/p/9e786be6dccb
?本篇只從程序的角度看如何使用它． 經(jīng)過FFT轉(zhuǎn)換的數(shù)據(jù)和轉(zhuǎn)換前長(zhǎng)度一致沮榜，每個(gè)數(shù)據(jù)分為實(shí)部和虛部?jī)刹糠峙陶ィ僭O(shè)時(shí)序時(shí)數(shù)長(zhǎng)度為N（Ｎ最好是2的整數(shù)次冪，這樣算起來更快）蟆融，用fft()轉(zhuǎn)換后：下標(biāo)為0和 N /2的兩個(gè)復(fù)數(shù)的虛數(shù)部分為0草巡，下標(biāo)為i和 N - i 的兩個(gè)復(fù)數(shù)共輒，也就是其虛部數(shù)值相同型酥、符號(hào)相反山憨。再用ifft()從頻域轉(zhuǎn)回時(shí)域之后，出現(xiàn)了由誤差引起的很小的虛部弥喉，用np.real()取其實(shí)部即可．
?由于一半是另一半的共軛郁竟，因此只需要關(guān)心一半數(shù)據(jù)．fft轉(zhuǎn)換后下標(biāo)為0的實(shí)數(shù)表示時(shí)域信號(hào)中的直流成分（不隨時(shí)間變化），下標(biāo)為i的復(fù)數(shù) a + bj由境，其中a表示余弦成分棚亩，b表示其正弦成分．

(2) 示例功能

?數(shù)據(jù)是航空乘客數(shù)據(jù)＂AirPassengers.csv＂，可以從CSDN下載虏杰，其中包括從1949-1960年蔑舞，每月旅客的數(shù)量，程序預(yù)測(cè)未來幾年的旅客數(shù)據(jù)．

?如圖所示嘹屯，數(shù)據(jù)為非平穩(wěn)數(shù)據(jù)攻询，其趨勢(shì)向上，且波動(dòng)加俱州弟，為將其變?yōu)槠椒€(wěn)數(shù)據(jù)钧栖， 先對(duì)其做了對(duì)數(shù)和差分處理．

(3) 示例代碼

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

# 函數(shù)功能：將頻域數(shù)據(jù)轉(zhuǎn)換成時(shí)序數(shù)據(jù)
# bins為頻域數(shù)據(jù)，n設(shè)置使用前多少個(gè)頻域數(shù)據(jù)婆翔，loop設(shè)置生成數(shù)據(jù)的長(zhǎng)度

def fft_combine(bins, n, loops=1):
 length = int(len(bins) * loops)
 data = np.zeros(length)
 index = loops * np.arange(0, length, 1.0) / length * (2 * np.pi)
 for k, p in enumerate(bins[:n]):
   if k != 0 : p *= 2 # 除去直流成分之外, 其余的系數(shù)都 * 2
   data += np.real(p) * np.cos(k*index) # 余弦成分的系數(shù)為實(shí)數(shù)部分
   data -= np.imag(p) * np.sin(k*index) # 正弦成分的系數(shù)為負(fù)的虛數(shù)部分
 return index, data

if __name__ == '__main__':
 data = pd.read_csv('AirPassengers.csv')
 ts = data['Passengers']

# 平穩(wěn)化
 ts_log = np.log(ts)
 ts_diff = ts_log.diff(1)
 ts_diff = ts_diff.dropna()
 print(fy[:10]) # 顯示前10個(gè)頻域數(shù)據(jù)
 fy = np.fft.fft(ts_diff)
 conv1 = np.real(np.fft.ifft(fy)) # 逆變換
 index, conv2 = fft_combine(fy / len(ts_diff), int(len(fy)/2-1), 1.3) # 只關(guān)心一半數(shù)據(jù)
 plt.plot(ts_diff)
 plt.plot(conv1 - 0.5) # 為看清楚拯杠，將顯示區(qū)域下拉0.5
 plt.plot(conv2 - 1)
plt.show()

(4) 運(yùn)行結(jié)果

[ 1.34992672+0. j -0.09526905-0.14569535j -0.03664114-0.12007802j
 -0.2670005 +0.24512406j -0.10075074+0.0314084 j -0.26409417+0.04197159j
 0.14411338+0.18703009j 0.07467991+0.05367644j -0.26663142+0.15324939j
 0.03248223+0.14130114j]

(5). 示例分析

?輸出的是fft轉(zhuǎn)換后的數(shù)據(jù)，只顯示了前十個(gè)啃奴，形式為復(fù)數(shù)．復(fù)數(shù)模(絕對(duì)值)的兩倍為對(duì)應(yīng)頻率的余弦波的振幅潭陪；復(fù)數(shù)的輻角表示對(duì)應(yīng)頻率的余弦波的相位。第0個(gè)元素表示直流分量，虛部為0．在數(shù)據(jù)中的位置標(biāo)記了頻率大小依溯，值標(biāo)記了振幅大欣涎帷．
?圖中顯示的三條曲線分別為原始數(shù)據(jù)，做了fft以及ifft逆變換后的數(shù)據(jù)黎炉，以及fft后自己實(shí)現(xiàn)算法還原并預(yù)測(cè)了未來的數(shù)據(jù)枝秤，從圖中可見，基本擬合了原始曲線慷嗜，預(yù)測(cè)曲線看起來也比較合理．
?上述方法可實(shí)現(xiàn)用傅里葉變換預(yù)測(cè)時(shí)序數(shù)據(jù)．與ARMA算法相比淀弹，它沒有明顯衰減，更適合長(zhǎng)時(shí)間的預(yù)測(cè)．
?對(duì)于隨時(shí)間變化的波形庆械，比如語音數(shù)據(jù)薇溃，一般使用加窗后做傅立葉變量的方法擬合數(shù)據(jù)．

3. 小波變換

(1) 相關(guān)知識(shí)

?有了傅立葉變換，為什么還用小波呢缭乘？上面提到沐序，如果波型隨時(shí)間變化，就需要對(duì)波型加窗分段后再處理忿峻，而且有時(shí)候需要大窗口，有時(shí)候需要小窗口辕羽，處理起來就更加麻煩．于是引入了更靈活的小波．
?傅立葉變換的基是正余弦函數(shù)逛尚，而小波的基是各種形狀小波，也就是說它把整個(gè)波形看成是多個(gè)位置和寬度不同的小波的疊加．小波有兩個(gè)變量：尺度a和平移量t刁愿，尺度控制小伸的伸縮绰寞，平移量控制小波的平移，它不需將數(shù)據(jù)切分成段铣口，就可以處理時(shí)變數(shù)據(jù)．尤其對(duì)突變信號(hào)處理得更好．
?下圖是幾種常見的小波．

?離散小波變換滤钱，Discrete Wavelet Transformatio (dwt)，可以說是小波變換中最簡(jiǎn)單的一種脑题。這里使用Python調(diào)用pywt庫(kù)實(shí)現(xiàn)最簡(jiǎn)單的功能．
?經(jīng)過變換之后的返回值：cA:Approximation（近似）, cD:Detail（細(xì)節(jié)）件缸，其中近似cA是周期性有規(guī)律的部分，可以被模擬和預(yù)測(cè)叔遂，而cD可看做是噪聲他炊。換言之，用此方法可以拆分周期性數(shù)據(jù)已艰，和其上的擾動(dòng)數(shù)據(jù)痊末。

(2) 示例功能

?示例使用的仍然是乘客數(shù)據(jù)，下面代碼是將細(xì)節(jié)D設(shè)為0哩掺，然后還原凿叠。

(3) 示例代碼

import pywt
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

if __name__ == '__main__':
 data = pd.read_csv('AirPassengers.csv')
 ts = data['Passengers']
 # 平穩(wěn)化
 ts_log = np.log(ts)
 ts_diff = ts_log.diff(1)
 ts_diff = ts_diff.dropna()

 cA,cD = pywt.dwt(ts_diff, 'db2')
 cD = np.zeros(len(cD))
 new_data = pywt.idwt(cA, cD, 'db2')

 plt.plot(ts_diff)
 plt.plot(new_data - 0.5)

 plt.show()

(4) 運(yùn)行結(jié)果

(5) 示例分析

?可以看到，用小波擬合的效果也還可以，一般可以使用小波擬合cA盒件，使用ARMA擬合cD部分蹬碧，兩種方法配合使用．