資料平滑處理-均值|中值|Savitzky-Golay濾波器

清風紫雪發表於2022-03-03
Go

均值濾波器

均值濾波器是一種使用頻次較高的線性濾波器。它的實現原理很簡單，就是指定一個長度大小為奇數的視窗，使用視窗中所有資料的平均值來替換中間位置的值，然後平移該視窗，平移步長為 1，繼續重複上述操作，直至滑動到時序資料的末尾，如此一來，對時序資料的過濾操作就結束了。均值濾波器的思路簡單，計算速度快，但是它容易被視窗中的極值點或者峰值所左右，不能很好地保留序列的邊緣資訊，在去噪的同時也對資料訊號的細節特徵產生了一定的破壞，不能很好地去除噪聲點，這極大地影響了模型的預測精度。
均值濾波的公式

 

其中，m表示視窗大小，xi表示視窗中的第 i個資料，x表示視窗中所有資料的均值。
python程式碼實現均值濾波
# 1. 均值濾波函式
def moving_average(data, window=5):
    size = window - 1
    arr_value = list(data.values)
    fill_left = arr_value[0]
    for i in range(size):
        arr_value.insert(0, fill_left)
    dat = pd.Series(arr_value)
    dat_roll = dat.rolling(window).mean()
    return dat_roll.dropna().reset_index(drop=True)

具體應用到資料上：高頻資料未完全擬合，其他資料也與源資料擬合程度一般

 

 中值濾波器

中值濾波是一種非線性的濾波演算法，它是將指定長度大小為奇數的視窗中的所有資料按從小到大的順序進行排列，並將排好序的資料的中值取代視窗中間的值。中值濾波克服了均值濾波所存在的問題，對視窗中的極端值不敏感，從而可以有效保留區域中的邊緣資訊，並且能有效抑制椒鹽噪聲和脈衝噪聲，避免細節特徵的丟失。但是，在面對均勻分佈的高斯噪聲時，它表現很差。
中值濾波的公式見式

 

 具體的python實現程式碼如下：

# 2. 中值濾波函式
def median_filter(data: pd.Series, window=5):
    return pd.Series(scipy.signal.medfilt(data, window))

應用到資料上：中值濾波將高頻擬合，其他資料也與源資料擬合程度一般

 

 SG濾波器

 對曲線進行平滑處理，通過Savitzky-Golay 濾波器，可以在scipy庫裡直接呼叫，不需要再定義函式。

python程式碼實現：

from scipy.signal import savgol_filter
# 3. Savitzky-Golay濾波函式
newans = savgol_filter(data, 5, 3, mode= 'nearest')

plt.plot(index,data,label='源網路流量',color='r',linestyle='-',marker='*')
plt.plot(index,newans,label='SG濾波網路流量',color='b')#新增linestyle設定線條型別
plt.legend()
plt.show()

# 備註：
data：代表曲線點座標（x,y）中的y值陣列
window_length：視窗長度，該值需為正奇整數。例如：此處取值5
k值：polyorder為對視窗內的資料點進行k階多項式擬合，k的值需要小於window_length。例如：此處取值3
mode：確定了要應用濾波器的填充訊號的擴充套件型別。（This determines the type of extension to use for the padded signal to which the filter is applied. ）

python原理手寫程式碼實現：

# 3. Savitzky-Golay濾波函式
"""
  data - list格式的1×n緯資料
  window_size - 擬合的視窗大小
  rank - 擬合多項式階次
  ndata - 修正後的值
"""
def savgol(data: list, window_size: int, rank: int):
    m = int((window_size - 1) / 2)
    odata = data[:]
    # 處理邊緣資料，首尾增加m個首尾項
    for i in range(m):
        odata.insert(0, odata[0])
        odata.insert(len(odata), odata[len(odata)-1])
    # 建立X矩陣
    x = create_x(m, rank)
    # 計算加權係數矩陣B
    b = (x * (x.T * x).I) * x.T
    a0 = b[m]
    a0 = a0.T
    # 計算平滑修正後的值
    ndata = []
    for i in range(len(data)):
        y = [odata[i + j] for j in range(window_size)]
        y1 = np.mat(y) * a0
        y1 = float(y1)
        ndata.append(y1)
    return ndata

"""
  建立係數矩陣X
  size - 2×size+1 = window_size
  rank - 擬合多項式階次
  x - 建立的係數矩陣
"""
def create_x(size, rank):
    x = []
    for i in range(2 * size + 1):
        m = i - size
        row = [m**j for j in range(rank)]
        x.append(row)
    x = np.mat(x)
    return x

newans = savgol(list(data), 5, 3)

plt.plot(index,data,label='源網路流量',color='r',linestyle='-',marker='*')
plt.plot(index,newans,label='SG濾波網路流量',color='b')#新增linestyle設定線條型別
plt.legend()
plt.show()

應用到資料上效果：此圖更加接近源曲線，並且將高頻部分進行平滑處理，擬合低頻部分

 

 

 

相關文章