【信号处理】MATLAB一维信号时序全特征分析（统计域/频域/时域）(附代码)

我们将探索时序特征分析，这是信号处理中的关键技术之一。信号在我们的日常生活中无处不在，从声音到图像，从传感器数据到金融市场的波动。理解信号的特征对于分类、检测和监测至关重要。因此，让我们深入了解这三类特征分析算法，它们分别基于统计域、谱域和时域。

基于统计域的时序特征分析

在统计域中，我们使用各种统计量来描述信号的整体特征。这些包括最大值、最小值、均值、中位数等。例如，最大值和最小值能够告诉我们信号的振幅范围，而均值和中位数则反映了信号的整体位置。此外，偏度和峰度等统计量能够描述信号分布的形状，揭示信号的非正态性。通过直方图和四分位距等进一步分析，我们可以深入了解信号的分布情况和波动性。

• 最大值(Maximum)、最小值(Minimum)、均值(Mean)、中位数(Median)：这些是描述信号整体振幅和位置的基本统计量，能够提供信号的整体特征。
• 偏度(Skewness)、峰度(Kurtosis)：描述信号分布形状的统计量，可以反映信号的非正态性。
• 直方图(Histogram)、四分位距(Interquartile Range)、绝对误差均值(Mean Absolute Deviation)、绝对误差中位数(Median Absolute Deviation)：提供了信号分布的更详细信息，有助于了解信号的波动情况。
• 均方根(Root Mean Square)、标准差(Standard Deviation)、方差(Variance)：用于衡量信号的离散程度和波动性。
• 经验分布函数百分位数(Empirical Distribution Function Percentile Count)、经验分布函数斜率(ECDF Slope)：提供了信号在不同百分位数处的分布情况，有助于发现信号中的异常值或特定趋势。

基于谱域的时序特征分析

谱域分析涉及将信号从时域转换为频域，以便我们能够观察信号在不同频率上的成分。这种分析的关键工具是快速傅里叶变换（FFT），它能够将信号转换为频谱表示。通过分析频率成分和能量分布，我们可以了解信号的频率特征。而小波变换则提供了对信号在不同尺度下变化的认识，例如小波绝对均值和小波标准差等特征能够描述信号的局部特性。

• 快速傅里叶变换(Fast Fourier Transform)：将信号从时域转换到频域，提供了信号在不同频率上的成分信息。
• 傅里叶变换平均系数(FFT Mean Coefficient)、小波变换(Wavelet Transform)：用于分析信号的频率成分和能量分布。
• 小波绝对均值(Wavelet Absolute Mean)、小波标准差(Wavelet Standard Deviation)、小波方差(Wavelet Variance)：描述小波变换后的信号特征，有助于了解信号在不同尺度下的变化。
• 谱距离(Spectral Distance)、频谱基频(Spectral Fundamental Frequency)、频谱最大频率(Spectral Maximum Frequency)、频谱中频(Spectral Median Frequency)、频谱最大峰值(Spectral Maximum Peaks)：提供了信号在频域上的更多详细信息，有助于分析信号的频率分布和频率成分。

基于时域的时序特征分析

时域分析涉及对信号在时间上的特征进行分析。自相关和质心等特征可以告诉我们信号在时间上的自相关性和集中程度。差分均值和差分绝对值均值等特征则描述了信号在不同时间点上的变化情况。此外，熵和波峰与波谷距离等特征可以揭示信号的周期性和波动性。

• 自相关(Autocorrelation)、质心(Centroid)：提供了信号在时间上的自相关性和集中程度。
• 差分均值(Mean Differences)、差分绝对值均值(Mean Absolute Differences)、差分中位数(Median Differences)、差分绝对值中位数(Median Absolute Differences)、差分绝对值之和(Sum of Absolute Differences)：描述信号在不同时间点上的变化情况。
• 熵(Entropy)、波峰与波谷距离(Peak to Peak Distance)、曲线覆盖面积(Area Under the Curve)、最大峰值个数(The Number of Maximum Peaks)、最小峰值个数(The Number of Minimum Peaks)、跨零率(Zero Crossing Rate)：提供了信号的更多时域特征，有助于发现信号的周期性、波动性和其他变化特征。

这些特征分析算法在实际应用中通常用于提取信号的特征，从而进行信号分类、故障检测、状态监测等任务。通过对不同域的特征进行综合分析，可以更全面地理解信号的特性和行为。

代码已经整理好，如下：

%%  生成示例一维信号
% 参数设置
fs = 1000; % 采样率
t = 0:1/fs:1; % 时间向量
f_carrier = 10; % 载波频率
A = 0.05; % 基础信号振幅

% 生成多个频率分量的基础信号
f_carrier = [50, 10, 20,100,1,9,5]; % 多个载波频率
carrier_signal = zeros(size(t));
for i = 1:length(f_carrier)
    carrier_signal = carrier_signal + A * sin(2*pi*f_carrier(i)*t);
end

% 添加表面缺陷扰动和噪声
defect_amplitude = 0.8; % 缺陷振幅
defect_location = 0.5; % 缺陷位置
defect_signal = zeros(size(t));
defect_signal(abs(t - defect_location-0.01)<0.01) = defect_amplitude;
defect_signal(abs(t - defect_location+0.01)<0.01) = -defect_amplitude;

noise_amplitude =0.1; % 噪声振幅
noise_signal_with_defect = noise_amplitude * randn(size(t));
% 合并信号（有缺陷）
signal_with_defect = carrier_signal + defect_signal + noise_signal_with_defect;
% 添加无缺陷的噪声
noise_signal_without_defect = noise_amplitude * randn(size(t));
% 合并信号（无缺陷）
signal_without_defect = carrier_signal + noise_signal_without_defect;

% 绘制信号图像
figure;
subplot(2,1,1);
plot(signal_with_defect);
title('Signal with Defect ');
xlabel('Distance (mm)');
ylabel('Amplitude');
axis([0 1000 -1 1])
grid on;

subplot(2,1,2);
plot(signal_without_defect);
title('Signal without Defect');
xlabel('Distance (mm)');
ylabel('Amplitude');
axis([0 1000 -1 1])
grid on;

%% 三类处理方法

signal = signal_with_defect;

% 计算基于统计域的时序特征
maximum_value = max(signal);
minimum_value = min(signal);
mean_value = mean(signal);
median_value = median(signal);
skewness_value = skewness(signal);
kurtosis_value = kurtosis(signal);
histogram_data = hist(signal, 10); % 计算直方图，分为10个箱子
interquartile_range = iqr(signal);
mean_absolute_deviation = mad(signal);
median_absolute_deviation = mad(signal, 1);
root_mean_square = rms(signal);
standard_deviation = std(signal);
variance = var(signal);
percentile_count = prctile(signal, [25, 50, 75]); % 计算25、50、75百分位数
ecdf_slope = gradient(ecdf(signal), range(signal)/numel(signal)); % 计算经验分布函数斜率

% 输出基于统计域的时序特征值
disp('---基于统计域的时序特征---');
disp(['最大值(Maximum): ', num2str(maximum_value)]);
disp(['最小值(Minimum): ', num2str(minimum_value)]);
disp(['均值(Mean): ', num2str(mean_value)]);
disp(['中位数(Median): ', num2str(median_value)]);
disp(['偏度(Skewness): ', num2str(skewness_value)]);
disp(['峰度(Kurtosis): ', num2str(kurtosis_value)]);
disp(['四分位距(Interquartile Range): ', num2str(interquartile_range)]);
disp(['绝对误差均值(Mean Absolute Deviation): ', num2str(mean_absolute_deviation)]);
disp(['绝对误差中位数(Median Absolute Deviation): ', num2str(median_absolute_deviation)]);
disp(['均方根(Root Mean Square): ', num2str(root_mean_square)]);
disp(['标准差(Standard Deviation): ', num2str(standard_deviation)]);
disp(['方差(Variance): ', num2str(variance)]);
disp(['经验分布函数百分位数(Empirical Distribution Function Percentile Count): ', num2str(percentile_count)]);
disp(['经验分布函数斜率(ECDF Slope): ', num2str(ecdf_slope(1))]);

% 计算基于谱域的时序特征
f_signal = fft(signal); % 计算傅里叶变换
fft_mean_coefficient = mean(abs(f_signal)); % 计算傅里叶变换系数的均值
wavelet_transform = cwt(signal, 1:64, 'cmor'); % 进行连续小波变换 确保您的 MATLAB 版本支持 cwt 函数
wavelet_absolute_mean = mean(abs(wavelet_transform(:))); % 计算小波变换的绝对均值
wavelet_standard_deviation = std(wavelet_transform(:)); % 计算小波变换的标准差
wavelet_variance = var(wavelet_transform(:)); % 计算小波变换的方差

% 输出基于谱域的时序特征值
disp('---基于谱域的时序特征---');
disp(['傅里叶变换平均系数(FFT Mean Coefficient): ', num2str(fft_mean_coefficient)]);
disp(['小波绝对均值(Wavelet Absolute Mean): ', num2str(wavelet_absolute_mean)]);
disp(['小波标准差(Wavelet Standard Deviation): ', num2str(wavelet_standard_deviation)]);
disp(['小波方差(Wavelet Variance): ', num2str(wavelet_variance)]);

% 计算基于时域的时序特征
autocorrelation = xcorr(signal); % 计算自相关函数
centroid = sum(t .* signal) / sum(signal); % 计算质心
mean_differences = mean(diff(signal)); % 计算差分均值
mean_absolute_differences = mean(abs(diff(signal))); % 计算差分绝对值均值
median_differences = median(diff(signal)); % 计算差分中位数
median_absolute_differences = median(abs(diff(signal))); % 计算差分绝对值中位数
sum_absolute_differences = sum(abs(diff(signal))); % 计算差分绝对值之和
entropy_value = entropy(signal); % 计算信号的熵
peak_to_peak_distance = max(signal) - min(signal); % 计算波峰与波谷距离
area_under_curve = trapz(t, signal); % 计算曲线覆盖面积
num_max_peaks = length(findpeaks(signal)); % 计算最大峰值个数
num_min_peaks = length(findpeaks(-signal)); % 计算最小峰值个数
zero_crossing_rate = sum(abs(diff(sign(signal)))) / (2 * length(signal)); % 计算跨零率

% 输出基于时域的时序特征值
disp('---基于时域的时序特征---');
disp(['自相关(Autocorrelation)(前三个): ', num2str(autocorrelation(1:3))]);
disp(['质心(Centroid): ', num2str(centroid)]);
disp(['差分均值(Mean Differences): ', num2str(mean_differences)]);
disp(['差分绝对值均值(Mean Absolute Differences): ', num2str(mean_absolute_differences)]);
disp(['差分中位数(Median Differences): ', num2str(median_differences)]);
disp(['差分绝对值中位数(Median Absolute Differences): ', num2str(median_absolute_differences)]);
disp(['差分绝对值之和(Sum of Absolute Differences): ', num2str(sum_absolute_differences)]);
disp(['熵(Entropy): ', num2str(entropy_value)]);
disp(['波峰与波谷距离(Peak to Peak Distance): ', num2str(peak_to_peak_distance)]);
disp(['曲线覆盖面积(Area Under the Curve): ', num2str(area_under_curve)]);
disp(['最大峰值个数(The Number of Maximum Peaks): ', num2str(num_max_peaks)]);
disp(['最小峰值个数(The Number of Minimum Peaks): ', num2str(num_min_peaks)]);
disp(['跨零率(Zero Crossing Rate): ', num2str(zero_crossing_rate)]);

运行结果：