自适应傅里叶分解(AFD)算法及其实现

自适应傅里叶分解(AFD)算法及其实现

AFD是一种自适应信号分解方法，能够根据信号特性自动选择最佳的基函数。

算法原理

核心思想

AFD通过迭代过程将信号分解为一组自适应选择的单分量函数，这些函数在Hardy空间中具有最优的能量聚集特性。

MATLAB实现

1. 基础AFD算法

function [components, residues, parameters] = adaptive_fourier_decomposition(signal, max_iter, tolerance)
    % 自适应傅里叶分解
    % 输入: signal - 输入信号, max_iter - 最大迭代次数, tolerance - 收敛容差
    % 输出: components - 分解分量, residues - 残差, parameters - 分解参数
    
    N = length(signal);
    t = (0:N-1)/N;  % 归一化时间
    
    components = zeros(max_iter, N);
    residues = zeros(max_iter, N);
    parameters = zeros(1, max_iter);
    
    % 初始化残差
    residues(1, :) = signal;
    
    for k = 1:max_iter
        % 计算当前残差的傅里叶变换
        R_k = fft(residues(k, :));
        
        % 寻找最优参数 a_k (在单位圆内)
        a_candidates = linspace(0.01, 0.99, 100);
        best_a = 0;
        max_inner_product = 0;
        
        for a = a_candidates
            % 生成自适应基函数
            basis_func = generate_basis_function(a, t, N);
            
            % 计算内积
            inner_prod = abs(sum(residues(k, :) .* conj(basis_func)));
            
            if inner_prod > max_inner_product
                max_inner_product = inner_prod;
                best_a = a;
            end
        end
        
        parameters(k) = best_a;
        
        % 生成最优基函数
        e_a = generate_basis_function(best_a, t, N);
        
        % 计算分解系数
        c_k = sum(residues(k, :) .* conj(e_a)) / sum(abs(e_a).^2);
        
        % 提取分量
        components(k, :) = c_k * e_a;
        
        % 更新残差
        if k < max_iter
            residues(k+1, :) = residues(k, :) - components(k, :);
        end
        
        % 检查收敛
        if norm(residues(k, :)) < tolerance
            components = components(1:k, :);
            residues = residues(1:k, :);
            parameters = parameters(1:k);
            break;
        end
    end
end

function basis_func = generate_basis_function(a, t, N)
    % 生成自适应傅里叶基函数
    % e_a(t) = sqrt(1-|a|^2)/(1 - a_bar * exp(2*pi*i*t))
    
    basis_func = zeros(1, N);
    for n = 1:N
        basis_func(n) = sqrt(1 - abs(a)^2) / (1 - conj(a) * exp(2*pi*1i*t(n)));
    end
end

2. 信号重构

function reconstructed_signal = afd_reconstruction(components, parameters, t)
    % AFD信号重构
    % 输入: components - 分解分量, parameters - 分解参数, t - 时间向量
    % 输出: reconstructed_signal - 重构信号
    
    [K, N] = size(components);
    reconstructed_signal = zeros(1, N);
    
    for k = 1:K
        % 使用存储的参数重新生成基函数
        e_a = generate_basis_function(parameters(k), t, N);
        
        % 计算系数
        c_k = sum(components(k, :) .* conj(e_a)) / sum(abs(e_a).^2);
        
        % 累加重构信号
        reconstructed_signal = reconstructed_signal + c_k * e_a;
    end
end

3. 完整示例

% 完整AFD示例
clear; clc;

% 生成测试信号
fs = 1000;  % 采样频率
t = 0:1/fs:1-1/fs;  % 1秒时间向量
N = length(t);

% 多分量信号
f1 = 5;    % 5Hz
f2 = 20;   % 20Hz  
f3 = 50;   % 50Hz

signal = 2*sin(2*pi*f1*t) + 1.5*cos(2*pi*f2*t) + sin(2*pi*f3*t) + 0.5*randn(1,N);

% AFD分解
max_iter = 10;
tolerance = 1e-6;

[components, residues, parameters] = adaptive_fourier_decomposition(signal, max_iter, tolerance);

% 信号重构
reconstructed_signal = afd_reconstruction(components, parameters, t);

% 计算重构误差
reconstruction_error = norm(signal - real(reconstructed_signal)) / norm(signal);

% 结果显示
figure('Position', [100, 100, 1200, 800]);

% 原始信号与重构信号
subplot(3, 2, 1);
plot(t, signal, 'b', 'LineWidth', 1.5);
hold on;
plot(t, real(reconstructed_signal), 'r--', 'LineWidth', 1);
legend('原始信号', '重构信号');
title('原始信号 vs 重构信号');
xlabel('时间 (s)'); ylabel('幅值');
grid on;

% 重构误差
subplot(3, 2, 2);
plot(t, signal - real(reconstructed_signal), 'g', 'LineWidth', 1);
title(['重构误差, RMSE: ', num2str(reconstruction_error)]);
xlabel('时间 (s)'); ylabel('误差');
grid on;

% 分解分量
subplot(3, 2, 3);
for k = 1:size(components, 1)
    plot(t, real(components(k, :)), 'LineWidth', 1);
    hold on;
end
title('AFD分解分量');
xlabel('时间 (s)'); ylabel('幅值');
grid on;

% 残差序列
subplot(3, 2, 4);
residue_norms = zeros(1, size(residues, 1));
for k = 1:size(residues, 1)
    residue_norms(k) = norm(residues(k, :));
end
semilogy(1:length(residue_norms), residue_norms, 'ro-', 'LineWidth', 2);
title('残差能量衰减');
xlabel('迭代次数'); ylabel('残差能量 (对数)');
grid on;

% 参数分布
subplot(3, 2, 5);
plot(parameters, 'bo-', 'LineWidth', 2, 'MarkerSize', 6);
title('自适应参数 a_k');
xlabel('分量索引'); ylabel('参数值');
grid on;

% 频谱比较
subplot(3, 2, 6);
f = (0:N-1)*fs/N;
Y_orig = abs(fft(signal));
Y_recon = abs(fft(real(reconstructed_signal)));
plot(f(1:N/2), Y_orig(1:N/2), 'b', 'LineWidth', 1.5);
hold on;
plot(f(1:N/2), Y_recon(1:N/2), 'r--', 'LineWidth', 1);
legend('原始频谱', '重构频谱');
title('频谱比较');
xlabel('频率 (Hz)'); ylabel('幅值');
grid on;

fprintf('AFD分解完成:\n');
fprintf('信号长度: %d\n', N);
fprintf('分解分量数: %d\n', size(components, 1));
fprintf('最终重构误差: %.6f\n', reconstruction_error);

4. 改进的快速AFD算法

function [components, energies] = fast_adaptive_fourier_decomposition(signal, target_components)
    % 快速自适应傅里叶分解
    % 基于能量准则选择最重要的分量
    
    N = length(signal);
    t = (0:N-1)/N;
    
    residual = signal;
    components = zeros(target_components, N);
    energies = zeros(1, target_components);
    
    for comp_idx = 1:target_components
        % 计算残差的频谱
        R_freq = fft(residual);
        
        % 寻找能量最大的频率区域
        [~, dominant_freq] = max(abs(R_freq(1:floor(N/2))));
        
        % 转换为自适应参数 (映射到单位圆内)
        a_k = 0.9 * exp(2*pi*1i*(dominant_freq-1)/N);
        
        % 生成基函数
        basis_func = generate_basis_function(a_k, t, N);
        
        % 正交投影
        coeff = sum(residual .* conj(basis_func)) / sum(abs(basis_func).^2);
        
        % 存储分量
        components(comp_idx, :) = coeff * basis_func;
        energies(comp_idx) = abs(coeff)^2 * sum(abs(basis_func).^2);
        
        % 更新残差
        residual = residual - components(comp_idx, :);
        
        % 提前终止条件
        if norm(residual) < 1e-6
            components = components(1:comp_idx, :);
            energies = energies(1:comp_idx);
            break;
        end
    end
end

参考代码 自适应傅里叶分解算法 www.3dddown.com/cna/79065.html

算法特点

优势

1. 自适应性: 基函数根据信号特性自动调整
2. 快速收敛: 通常比传统傅里叶分解需要更少的分量
3. 能量集中: 分量具有更好的时频聚集性
4. 正交性: 分解分量在Hardy空间中正交

应用领域

• 信号去噪和压缩
• 故障诊断和特征提取
• 生物医学信号处理
• 通信信号分析