如何从加速度计数据中提取冻结指数进行FoG检测？

Question

我正在做一个项目，检测帕金森病患者的冻结步态发作情况，根据这个 纸张 和其他的，提取冻结指数是"冻结 "频段（3-8Hz）的功率除以运动频段（0.5-3Hz）的功率。"将是一个很好的功能。

下面是解释。从原始信号中提取的一个标准特征是冻结指数(FI)，定义为所谓的冻结频段和运动频段(分别为3-8赫兹和0.5-3赫兹)所包含的功率之比。这个功能很方便，因为它只需要FFT计算。

但我无法在Python中实现它。

我有一个这样的数据框架。然后，我做了这样的事情，从传感器的时间序列数据中提取特征。

win_size=200
step_size=40
for col in df.columns:
    if col != 'Label':
        df[col + '_avg'] = df[col].rolling(win_size).mean()[step_size - 1::step_size]

现在，我想提取冰冻指数，我该怎么做呢？ 我需要有人给我解释一下，因为我并不完全理解。

Answer 1

我发现这个 文章 是非常有用的。本文中有一个函数可以计算特定频段内信号的平均功率。而由于冻结指数是所谓的冻结频段和运动频段（分别为3-8赫兹和0.5-3赫兹）所包含的功率之间的比值，所以我们可以用这个函数得到每个频段的功率，然后再除以它们。

这就是这个函数。

def bandpower(data, sf, band, window_sec=None, relative=False):
    """Compute the average power of the signal x in a specific frequency band.

    Parameters
    ----------
    data : 1d-array
        Input signal in the time-domain.
    sf : float
        Sampling frequency of the data.
    band : list
        Lower and upper frequencies of the band of interest.
    window_sec : float
        Length of each window in seconds.
        If None, window_sec = (1 / min(band)) * 2
    relative : boolean
        If True, return the relative power (= divided by the total power of the signal).
        If False (default), return the absolute power.

    Return
    ------
    bp : float
        Absolute or relative band power.
    """
    from scipy.signal import welch
    from scipy.integrate import simps
    band = np.asarray(band)
    low, high = band

    # Define window length
    if window_sec is not None:
        nperseg = window_sec * sf
    else:
        nperseg = (2 / low) * sf

    # Compute the modified periodogram (Welch)
    freqs, psd = welch(data, sf, nperseg=nperseg)

    # Frequency resolution
    freq_res = freqs[1] - freqs[0]

    # Find closest indices of band in frequency vector
    idx_band = np.logical_and(freqs >= low, freqs <= high)

    # Integral approximation of the spectrum using Simpson's rule.
    bp = simps(psd[idx_band], dx=freq_res)

    if relative:
        bp /= simps(psd, dx=freq_res)
    return bp

然后，我创建了这个简单的函数来返回FI。

def freeze_index(data, sf, band1, band2, win_sec=None, relative=False):
    fi = bandpower(data, sf, band1, win_sec, relative) / bandpower(data, sf, band2, win_sec, relative)
    return fi

这是我在滚动窗口函数中的调用方法。

for col in df.columns:
    if col != 'Label':
        df[col + '_fi'] = df[col].rolling(win_size).apply(freeze_index, args=(fs, [3, 8], [0.5, 3], win_size,))[step_size - 1::step_size]

希望这是正确的解决方案，希望对你有所帮助。