星期日, 16. 六月 2019 07:41下午
最近忙于項(xiàng)目，因此本期博客就簡(jiǎn)單梳理一下信號(hào)(筆者做的多數(shù)情況下是三相電流數(shù)據(jù))下的基本時(shí)域和頻域特征。
提取特征可以用tsfresh這個(gè)庫(kù)來(lái)實(shí)現(xiàn)，但是這里還是采用了numpy。
該部分的程序代碼筆者寫了一個(gè)可以直接調(diào)用的類，見筆者git

1 時(shí)域特征

假設(shè)共有m條數(shù)據(jù)，每一條數(shù)據(jù)長(zhǎng)度為n，第ｉ條數(shù)據(jù)第ｊ個(gè)數(shù)據(jù)點(diǎn)用\(z_{ij}\)表示，并且z為數(shù)組，不能是列表，否則以下一些程序會(huì)報(bào)錯(cuò)。以下主要從數(shù)學(xué)公式和python實(shí)現(xiàn)來(lái)敘述。

(1) 含量綱的時(shí)域特征

含量綱的時(shí)域特征，筆者簡(jiǎn)單整理出了十個(gè)，其中包括最大值(maximum)、最小值(minimum)、極差(range)、均值(mean)、中位數(shù)(media)、眾數(shù)(mode)、標(biāo)準(zhǔn)差(standard deviation)、均方根值(root mean square/rms)、均方值(mean square/ms)、k階中心/原點(diǎn)矩。
導(dǎo)入所需要的庫(kù)

import numpy as np

最大值

max_z = np.max(z, axis=1)

最小值

min_z = np.min(z, axis=1)   

極差

range_z = np.max(z, axis=1)-np.min(z, axis=1)

均值

mean_z = np.mean(z, axis=1)

中位數(shù)

將一組數(shù)從小到大排序，出現(xiàn)在中間的數(shù)(當(dāng)n為奇數(shù)時(shí))或者中間兩個(gè)數(shù)的均值(當(dāng)n為偶數(shù)時(shí))

media_z = np.median(z, axis=1)

眾數(shù)

一組數(shù)從大到小排序，出現(xiàn)次數(shù)最多的數(shù)(當(dāng)有多個(gè)數(shù)出現(xiàn)次數(shù)一樣，取最小的數(shù))

import scipy.stats
mode_z = scipy.stats.mode(z, axis=1)[0].reshape([-1])

標(biāo)準(zhǔn)差

std_z = np.std(z, axis=1)

均方根值

from sklearn.metrics import mean_squared_error
rms_z = [np.sqrt(mean_squared_error(zi, np.zeros(len(zi)))) for zi in z]
rms_z = np.array(rms_z)

均方值(二階原點(diǎn)矩)

from sklearn.metrics import mean_squared_error
ms_z = [mean_squared_error(zi, np.zeros(len(zi))) for zi in z]
ms_z = np.array(ms_z)

ｋ階中心矩/原點(diǎn)矩

k階中心矩

k階原點(diǎn)矩

def k_order_moment(z, k, is_center=True, is_origin=True):
    """
    Calculate k-order center moment and k-order origin moment of z
    :param z: array_like
    :param k: int
    :param is_center: bool; whether calculate k-order center moment
    :param is_origin: bool; whether calculate k-order origin moment
    :return: tuple; return k-order center moment and k-order origin moment
    """
    if (is_center is False) and (is_origin is False):
        raise ValueError("At least one of is_center and is_origin is True")
    if (type(k) is not int) or (k < 0):
        raise TypeError("k must be a integrate and more than 0")
    if type(z) is list:
        z = np.array(z)
        
    mean_z = np.mean(z, axis=1)
    if is_origin is False:
        k_center = np.mean([(z[i]-mean_z[i])**k for i in range(z.shape[0])], axis=1)
        return (k_center, None)
    if is_center is False:
        k_origin = np.mean([z[i]**k for i in range(z.shape[0])], axis=1)
        return (None, k_origin)
    if is_center and is_origin:
        k_center = np.mean([(z[i] - mean_z[i]) ** k for i in range(z.shape[0])], axis=1)
        k_origin = np.mean([z[i] ** k for i in range(z.shape[0])], axis=1)
        return (k_center, k_origin)

(2) 無(wú)量綱的時(shí)域特征

無(wú)量綱的時(shí)域筆者主要列舉了６個(gè)，分別為偏度(skewness)，峰度(kurtosis)，峰度因子(kurtosis factor)、波形因子(waveform factor)、脈沖因子(pulse factor)、裕度因子(margin factor)。

偏度(三階標(biāo)準(zhǔn)矩)

\(\mu,\sigma\)為總體均值和標(biāo)準(zhǔn)差，不是樣本均值和標(biāo)準(zhǔn)差!!!
偏度可通過(guò)下面兩種方法計(jì)算：
方法１：

方法２：

python和大多數(shù)軟件采用方法２求偏度

import pandas as pd
skew_z = pd.DataFrame(z.transpose()).skew().values

峰度(四階標(biāo)準(zhǔn)矩)

為總體均值和標(biāo)準(zhǔn)差，不是樣本均值和標(biāo)準(zhǔn)差!!!
峰度同偏度一樣也有兩種方法計(jì)算：
方法１：

方法２()：

python和大多數(shù)軟件采用方法２求峰度

import pandas as pd
kurt_z = pd.DataFrame(z.transpose()).kurt().values

峰度因子

kurt_factor_z = max_z/rms_z

波形因子

wave_factor_z = rms_z/mean_z

脈沖因子

pulse_factor_z = max_z/abs(mean_z)

裕度因子

margin_factor_z = max_z/ms_z