基于经验模态分解或局部均值分解的信号特征提取

基于经验模态分解（EMD）或局部均值分解（LMD）的信号特征提取。

通常将能量密度与平均周期的乘积 称为 K值，或者更专业地称为 加权平均频率（Weighted Average Frequency） 或 边际谱的一阶矩（First Moment of the Marginal Spectrum）。这个参数是Hilbert-Huang变换谱分析中的一个核心特征量。

下面我将详细解释这个概念、其物理意义、计算方法以及应用。

1. 核心概念解析

首先，我们明确两个基本量：

• 能量密度 (Energy Density, E_i): 对于EMD/LMD分解得到的第i个本征模态函数（IMF）或乘积函数（PF），其能量通常定义为该分量信号幅值的平方和或方差。能量密度则描述了该分量单位时间（或单位样本）内所包含的能量大小。
  E_i = sum( IMF_i(t)^2 ) / N （N为数据点数）

• 平均周期 (Average Period, T_i): 对于单分量的振荡模式，其平均周期可以通过计算过零点数（Zero-Crossings） 或极值点数（Extrema） 来估计。
  T_i = (N - 1) / (Number of Zero Crossings) * 采样间隔

K值 (K_i) 就是这两个量的乘积：
K_i = E_i * T_i

2. 物理意义：为什么这个乘积很重要？

K_i 的物理意义非常深刻，它可以被理解为 第i个IMF分量所代表的振荡模式对原信号“平均能量”的贡献。

• 从量纲分析看: E_i 的量纲是 [振幅]²，T_i 的量纲是 [时间]。K_i 的量纲是 [振幅]² * [时间]。这与功率谱密度（Power Spectral Density, PSD） 图的积分面积（即总功率）的量纲是一致的。K_i 描绘了在Hilbert谱或边际谱中，不同频率分量所携带的“能量”在频率轴上的分布情况。

• 从信息聚合看: K_i 将一个IMF分量的两个最基本属性（强度E_i和尺度T_i）融合成了一个单一的特征值。这使得它成为一个非常强大和简洁的信号特征提取工具。相比于单独使用能量或频率，K_i 能更好地区分不同物理过程产生的模态。

3. 计算方法与步骤

计算每个IMF分量的 K_i 值的标准流程：

1. 信号分解:

   • 使用EMD或LMD算法将原始信号 x(t) 分解为一系列IMF/PF分量 c_i(t) 和一个残差 r_n(t)。
   • x(t) = sum_{i=1}^{n} c_i(t) + r_n(t)

2. 计算每个IMF的能量密度 E_i:

   • E_i = (1/N) * sum_{t=1}^{N} [c_i(t)]^2
   • 这里用方差作为能量的一种度量。

3. 计算每个IMF的平均周期 T_i:

   • 过零点法: 计算该IMF分量穿过零点的次数 N_z。平均周期为：
     T_i = (N - 1) / (N_z * Δt) （Δt是采样间隔，单位是秒/点）
   • 极值点法: 计算该IMF分量的局部极大值点个数 N_e。平均周期为：
     T_i = (N - 1) / (N_e * Δt)
   • 过零点法通常更常用。

4. 计算每个IMF的K值 K_i:

   • K_i = E_i * T_i

5. (可选) 构建特征向量:

   • 将所有IMF分量的 K_i 值按顺序排列，可以形成一个特征向量 [K_1, K_2, ..., K_n]，该向量可以高度概括原始信号的频率-能量结构，用于后续的机器学习分类或故障诊断。

4. MATLAB

计算EMD分解后各IMF能量和K值

function [K_values, E_values, T_values] = calculate_K_EMD(signal, fs)
    % calculate_K_EMD 计算信号EMD分解后各IMF的K值
    %  输入:
    %   - signal: 输入信号 (向量)
    %   - fs: 采样频率 (Hz)
    %  输出:
    %   - K_values: 各IMF的K值 (向量)
    %   - E_values: 各IMF的能量 (向量)
    %   - T_values: 各IMF的平均周期 (向量)

    % 1. 执行EMD分解
    imf = emd(signal, 'Interpolation', 'pchip'); % 需要安装Signal Processing Toolbox
    
    [num_imf, N] = size(imf); % num_imf: IMF个数, N: 信号长度
    dt = 1/fs; % 采样间隔

    % 初始化输出
    K_values = zeros(1, num_imf);
    E_values = zeros(1, num_imf);
    T_values = zeros(1, num_imf);

    for i = 1:num_imf
        current_imf = imf(i, :);
        
        % 2. 计算能量密度 E_i
        E_i = sum(current_imf.^2) / N;
        E_values(i) = E_i;
        
        % 3. 计算平均周期 T_i (过零点法)
        % 找到过零点的索引 (从正到负或负到正都算)
        zero_crossings = find(diff(sign(current_imf)) ~= 0);
        num_zero_cross = length(zero_crossings);
        
        % 避免除零错误
        if num_zero_cross > 1
            % 平均周期 = 总时间 / 过零点间隔数
            % 总时间 = (N-1)*dt
            % 过零点间隔数 = num_zero_cross
            T_i = (N-1)*dt / num_zero_cross;
        else
            T_i = 0; % 如果没有或只有一个过零点，周期视为无穷大或无效，设为0
        end
        T_values(i) = T_i;
        
        % 4. 计算 K_i
        K_values(i) = E_i * T_i;
    end
    
    % 可视化结果
    figure;
    subplot(2, 2, 1);
    plot(signal);
    title('原始信号');
    xlabel('采样点');
    
    subplot(2, 2, 2);
    plot(imf');
    title('EMD分解结果 (IMFs)');
    xlabel('采样点');
    
    subplot(2, 2, 3);
    stem(E_values);
    title('各IMF能量密度 E_i');
    xlabel('IMF 序号');
    
    subplot(2, 2, 4);
    bar(K_values);
    title('各IMF K 值 (E_i * T_i)');
    xlabel('IMF 序号');
    ylabel('K_i');
end

% 使用示例:
% 生成一个测试信号
fs = 1000; % 采样率1kHz
t = 0:1/fs:2;
signal = sin(2*pi*10*t) + 0.5*sin(2*pi*50*t) + 0.2*randn(size(t)); % 10Hz + 50Hz + 噪声

% 调用函数计算K值
[K_vals, E_vals, T_vals] = calculate_K_EMD(signal, fs);

5. 应用场景

K_i 值最大的用途是作为信号的特征向量，广泛应用于状态监测和模式识别：

1. 机械故障诊断:

   • 轴承故障: 不同位置（内圈、外圈、滚动体）的故障会激发不同频率的共振，产生不同K_i值的IMF，从而可以区分故障类型。
   • 齿轮箱故障: 齿根裂纹、断齿、磨损等故障会导致调制现象，产生特定的边频带IMF分量，其K_i值具有明显的区分度。

2. 生物医学信号处理:

   • EEG/ECG信号分析: 不同节律的脑电波（δ, θ, α, β, γ）或心电图的P波、QRS波群、T波对应着不同的振荡模式，其K_i值可以作为疾病诊断的特征。

3. 结构健康监测:

   • 桥梁、建筑物的振动信号在损伤前后，其模态频率和能量会发生变化，K_i值可以有效地捕捉这种变化。

参考代码 EMD或者LMD的能量密度与其平均周期的乘积K值 www.youwenfan.com/contentcnh/55174.html

总结

特性	描述
定义	K_i = E_i * T_i，即IMF分量的能量密度与其平均周期的乘积。
物理意义	表征了该振荡模式对原信号平均能量的贡献，是Hilbert谱分析中的核心特征量。
优势	将强度（能量） 和尺度（周期） 信息融合为一个简洁、有物理意义的特征值。
主要用途	特征提取，用于模式识别、故障诊断和状态分类。构建特征向量 [K_1, K_2, ..., K_n] 来描述信号。
计算方法	1. EMD/LMD分解 → 2. 计算各IMF能量 E_i → 3. 计算各IMF平均周期 T_i → 4. 相乘得到 K_i。

简而言之，K值是一个能够同时反映信号特定振荡模式强度和尺度的强大特征指标，是利用EMD/LMD进行信号分析时不可或缺的工具。