HHT如何进行去噪去趋势

Question

我写了个HHT变化去噪去趋势的matlab代码，但是我随机写了个信号，并不能达到去噪效果，哪里出了问题呢。是根据归一化互相关性选取IMF的。

% 生成复杂信号，包括正弦波、方波和高斯噪声
Fs = 1000; % 采样率
t = 0:1/Fs:1; % 时间向量
signal = sin(2*pi*5*t) + square(2*pi*10*t) + 0.5*randn(size(t)); % 信号
%signal=XTrain(170,:);


% 进行 EMD 分解，假设结果存储在 imfs 矩阵中
imfs = emd(signal);

% 计算每个 IMF 与原始信号的交叉相关性
N = size(imfs, 2); % IMF 的数量
ncc_values = zeros(1, N);

figure;
subplot(N+1, 1, 1);
plot(signal, 'k');
title('原始复杂信号');
for i = 1:N
    subplot(N+1, 1, i+1);
    plot(imfs(:, i), 'b');
    title(['IMF ', num2str(i)]);
end

for i = 1:N
    % 计算滞后为 0 时的交叉相关性 Rxyi_0
    Rxyi_0 = sum(signal .* conj(imfs(:, i))');
    
    % 计算分母 denominator，用于归一化
    denominator = sqrt(sum(abs(signal).^2) * sum(abs(imfs(:, i)).^2));

    % 计算归一化滞后交叉相关性 ncc
    ncc = Rxyi_0 / denominator;
    ncc_values(i) = ncc;
end

% 计算阈值 TH
max_ncc = max(ncc_values);
TH = max_ncc / (10 * max_ncc - 3);

% 选择最相关的 IMF（假设需要大于阈值 TH）
relevant_imf_indices = find(ncc_values > TH);

% 重构信号，只使用最相关的 IMF
reconstructed_signal = sum(imfs(:,2:6), 2);

% 绘制原始信号和重构后的信号
figure;
subplot(2,1,1);
plot(signal);
title('Original Signal');

subplot(2,1,2);
plot(real(reconstructed_signal)); % 使用 real() 取实部
title('Reconstructed Signal (Selected IMF)');

公式是如下两张图片

Answer 1

【以下回答由 GPT 生成】

使用Hilbert-Huang变换（HHT）对信号进行去噪和去趋势处理是一个常用的方法，下面是对你提供的代码进行分析和改进的步骤：

1. 生成复杂信号部分：

Fs = 1000; #采样率
t = 0:1/Fs:1; #时间向量
signal = sin(2*pi*5*t) + square(2*pi*10*t) + 0.5*randn(size(t)); #信号

问题：生成信号的过程没有问题，包括正弦波、方波和高斯噪声。

1. EMD分解部分：

imfs = emd(signal);

问题：EMD分解是HHT的第一步，将信号分解为一组本征模态函数（Intrinsic Mode Functions, IMF）。此步骤没有问题。

1. 计算每个IMF与原始信号的相关性部分：

for i = 1:N
    #计算滞后为 0 时的交叉相关性 Rxyi_0
    Rxyi_0 = sum(signal .* conj(imfs(:, i))');

    #计算分母 denominator，用于归一化
    denominator = sqrt(sum(abs(signal).^2) * sum(abs(imfs(:, i)).^2));

    #计算归一化滞后交叉相关性 ncc
    ncc = Rxyi_0 / denominator;
    ncc_values(i) = ncc;
end

问题：计算相关性的公式没有问题，但代码中使用了复数的共轭运算（conj），这可能会导致错误的相关性计算结果。在这里，我们只需要对信号和IMF按元素进行乘法，不需要复数运算。因此，我们需要将conj(imfs(:, i))更改为imfs(:, i)。

4.计算阈值和选择IMF部分：

max_ncc = max(ncc_values);
TH = max_ncc / (10 * max_ncc - 3);

relevant_imf_indices = find(ncc_values > TH);

问题：计算阈值和选择IMF的方法没有问题。

5.重构信号部分：

reconstructed_signal = sum(imfs(:, 2:6), 2);

问题：在重构信号时，你只使用了IMF 2到6。你可以将此部分更改为使用具有最高相关性的IMF，即将2:6更改为relevant_imf_indices。

综上所述，下面是改进后的代码：

# 生成复杂信号，包括正弦波、方波和高斯噪声
Fs = 1000; # 采样率
t = 0:1/Fs:1; # 时间向量
signal = sin(2*pi*5*t) + square(2*pi*10*t) + 0.5*randn(size(t)); # 信号

# 进行 EMD 分解，假设结果存储在 imfs 矩阵中
imfs = emd(signal);

# 计算每个 IMF 与原始信号的交叉相关性
N = size(imfs, 2); # IMF 的数量
ncc_values = zeros(1, N);

figure;
subplot(N+1, 1, 1);
plot(signal, 'k');
title('原始复杂信号');
for i = 1:N
    subplot(N+1, 1, i+1);
    plot(imfs(:, i), 'b');
    title(['IMF ', num2str(i)]);
end

for i = 1:N
    # 计算滞后为 0 时的交叉相关性 Rxyi_0
    Rxyi_0 = sum(signal .* imfs(:, i)');

    # 计算分母 denominator，用于归一化
    denominator = sqrt(sum(abs(signal).^2) * sum(abs(imfs(:, i)).^2));

    # 计算归一化滞后交叉相关性 ncc
    ncc = Rxyi_0 / denominator;
    ncc_values(i) = ncc;
end

# 计算阈值 TH
max_ncc = max(ncc_values);
TH = max_ncc / (10 * max_ncc - 3);

# 选择最相关的 IMF（假设需要大于阈值 TH）
relevant_imf_indices = find(ncc_values > TH);

# 重构信号，只使用最相关的 IMF
reconstructed_signal = sum(imfs(:, relevant_imf_indices), 2);

# 绘制原始信号和重构后的信号
figure;
subplot(2, 1, 1);
plot(signal, 'k');
title('原始信号');
subplot(2, 1, 2);
plot(reconstructed_signal, 'r');
title('重构信号');

希望这些改进可以帮助你解决信号去噪和去趋势处理的问题。如果还有其他问题，请随时提问！

如果你已经解决了该问题, 非常希望你能够分享一下解决方案, 写成博客, 将相关链接放在评论区, 以帮助更多的人 ^-^