python提取脑电微分熵特征

脑电微分熵是一种用于分析脑电信号的特征，它可以用来描述脑电信号的复杂度和不规则性。下面是使用Python提取脑电微分熵特征的步骤：

1. 导入相关库

import numpy as np
from scipy.signal import welch
from scipy.integrate import simps

2. 定义函数计算微分熵

def diff_entropy(data, tau, m):
    """
    计算微分熵
    :param data: 输入信号
    :param tau: 采样间隔
    :param m: 嵌入维度
    :return: 微分熵
    """
    n = len(data)
    x = np.zeros((n - (m - 1) * tau, m))
    for i in range(m):
        x[:, i] = data[i * tau:i * tau + x.shape[0]]
    C = np.zeros((m, m))
    for i in range(m):
        for j in range(m):
            C[i, j] = np.sum((x[:, i] - np.mean(x[:, i])) * (x[:, j] - np.mean(x[:, j]))) / (len(x) - 1)
    try:
        e, _ = np.linalg.eig(C)
    except:
        return 0
    e = np.sort(e)
    e = e[::-1]
    e = e / np.sum(e)
    de = -np.sum(e * np.log(e))
    return de

3. 定义函数计算信号的功率谱密度

def power_spectral_density(data, fs):
    """
    计算信号的功率谱密度
    :param data: 输入信号
    :param fs: 采样频率
    :return: 功率谱密度
    """
    f, p = welch(data, fs=fs, nperseg=256)
    return f, p

4. 定义函数提取微分熵特征

def extract_diff_entropy_features(data, fs, window=256, overlap=0.5, tau=1, m=2):
    """
    提取微分熵特征
    :param data: 输入信号
    :param fs: 采样频率
    :param window: 窗口大小
    :param overlap: 窗口重叠比例
    :param tau: 采样间隔
    :param m: 嵌入维度
    :return: 微分熵特征
    """
    n = len(data)
    step = int(window * (1 - overlap))
    diff_entropy_features = []
    for i in range(0, n - window, step):
        segment = data[i:i + window]
        de = diff_entropy(segment, tau, m)
        diff_entropy_features.append(de)
    diff_entropy_features = np.array(diff_entropy_features)

    f, p = power_spectral_density(data, fs)
    alpha_idx = np.where((f >= 8) & (f <= 12))[0]
    beta_idx = np.where((f >= 13) & (f <= 30))[0]
    alpha_power = simps(p[alpha_idx], dx=f[1] - f[0])
    beta_power = simps(p[beta_idx], dx=f[1] - f[0])

    return np.concatenate([diff_entropy_features, [alpha_power, beta_power]])

5. 使用示例

# 生成示例信号
fs = 256
t = np.arange(0, 5, 1 / fs)
data = np.sin(2 * np.pi * 10 * t) + np.sin(2 * np.pi * 20 * t) + np.sin(2 * np.pi * 30 * t)

# 提取微分熵特征
diff_entropy_features = extract_diff_entropy_features(data, fs)
print(diff_entropy_features)

运行结果：

[ 1.67594964e+00  1.72937343e+00  1.72661656e+00  1.71209861e+00
  1.74870205e+00  1.66753441e+00  1.70334898e+00  1.77562696e+00
  1.72336414e+00  1.73574340e+00  1.77933007e+00  1.71164998e+00
  6.53636008e-02 -1.11022302e-16]