核心流程首先从原始振动信号出发，通过多几何并行的特征提取架构同时进行四个维度的几何分析：在仿射几何维度，系统计算信号的仿射不变量矩、仿射曲率和仿射弧长等特征，这些特征在尺度变换、剪切变形和平移操作下保持不变，提高了模型对工况变化的鲁棒性；在复流形维度，算法通过希尔伯特变换将信号映射到复数域，提取幅度、相位、瞬时频率等复数特征，并利用复数小波变换获取时频域的复数系数，全面捕捉信号的复数域几何结构；在芬斯勒几何维度，系统计算方向相关的芬斯勒度量、芬斯勒曲率和Berwald度量，处理信号的非对称性和方向依赖性，为非线性动力学分析提供更灵活的几何框架；在子流形几何维度，算法通过相空间重构将一维信号嵌入高维相空间，将信号视为高维流形中的低维子流形，计算主曲率、高斯曲率、平均曲率、体积、表面积及拓扑不变量等几何特征。
随后，将这四种几何理论提取的特征进行维度对齐和融合，形成综合了60维几何特征的多维表征，最后使用梯度提升树分类器进行训练和优化，实现对轴承四种工作状态的精准智能诊断。整个系统不仅实现了卓越的分类性能，还通过特征重要性分析揭示了各几何理论对故障诊断的贡献度，为故障机理分析提供了理论依据。
详细算法步骤
第一步，数据采集与信号预处理，系统从轴承振动传感器获取原始的一维时序振动信号，对信号进行去除直流分量和标准化处理以消除量纲影响，然后将长序列分割成固定长度的重叠样本段，为后续的几何特征提取准备规整的数据格式。
第二步，仿射几何不变性特征提取，计算信号的仿射不变量矩如Hu矩等特征，这些特征在尺度变换、剪切变形和平移操作下保持不变，同时计算仿射曲率和仿射弧长参数化特征，增强模型对不同工况下信号变化的鲁棒性。
第三步，复流形结构特征分析，通过希尔伯特变换将实数信号转换为复数解析信号，提取幅度、相位和瞬时频率等复数域特征，利用复数小波变换获取信号的时频域复数系数，分析复平面上的几何结构和相位信息。
第四步，芬斯勒几何特征提取，计算方向相关的芬斯勒度量来捕捉信号的非对称特性，分析芬斯勒曲率和Berwald度量，处理信号在不同方向上的变化敏感性，提供更灵活的几何描述框架。
第五步，子流形几何特征分析，通过相空间重构技术将一维信号嵌入高维相空间，将信号轨迹视为高维流形中的低维子流形，计算主曲率、高斯曲率、平均曲率等几何不变量以及体积、表面积等几何测度。
第六步，多几何特征融合与对齐，将四种几何理论提取的特征向量进行维度统一和拼接处理，形成综合的多维特征表示，全面表征信号在不同几何视角下的特性。
第七步，智能分类器训练与优化，使用梯度提升树算法在融合后的多维特征上进行模型训练，通过集成多个弱学习器构建强分类器，优化模型参数以获得最佳分类性能。
第八步，故障诊断与性能评估，将新的振动信号输入训练好的模型进行故障类型预测，通过准确率、精确率、召回率、F1分数等指标全面评估系统性能，并生成混淆矩阵和特征重要性可视化图表。
效果分析
在训练集和测试集上均实现了完美的100%准确率和F1分数，对四种轴承工作状态（正常、滚珠故障、外圈故障、内圈故障）实现了完全无误差的分类识别，得益于系统融合了四种不同的高等几何理论：通过仿射几何提取尺度、剪切和平移变换下保持不变的鲁棒特征；通过复流形分析在复数域同时捕捉信号的幅度和相位信息；通过芬斯勒几何处理方向依赖的非对称特征；通过子流形几何分析高维流形中的低维嵌入结构，最终将这四种几何理论提取的共60维特征融合后输入梯度提升树分类器进行训练，形成了对轴承振动信号从时频域、复数域、仿射变换域到流形几何域的全方位多层次特征表征，从而在测试集248个样本上展现了完美的分类能力，验证了多几何融合策略在复杂故障诊断任务中的强大优势和理论创新价值。
class AffineInvariantFeatureExtractor: """仿射几何不变性特征提取 - 创新点：保持仿射变换下的特征不变性""" def __init__(self, n_moments=10, n_invariants=20): self.n_moments = n_moments self.n_invariants = n_invariants self.scaler = StandardScaler() def compute_affine_moments(self, signal): """计算仿射不变量矩""" n = len(signal) x = np.arange(n) # 中心矩 mu_00 = np.sum(signal) mu_10 = np.sum(x * signal) mu_01 = np.sum(signal**2) mu_20 = np.sum(x**2 * signal) mu_11 = np.sum(x * signal**2) mu_02 = np.sum(signal**3) # 归一化矩 if mu_00 > 0: eta_10 = mu_10 / mu_00 eta_01 = mu_01 / mu_00 eta_20 = mu_20 / mu_00 eta_11 = mu_11 / mu_00 eta_02 = mu_02 / mu_00 else: eta_10 = eta_01 = eta_20 = eta_11 = eta_02 = 0 # 仿射不变量矩 (Affine Moment Invariants) invariants = [] # Hu矩 (平移、旋转、尺度不变) nu_20 = eta_20 - eta_10**2 nu_11 = eta_11 - eta_10 * eta_01 nu_02 = eta_02 - eta_01**2 # 一阶不变量 I1 = nu_20 + nu_02 invariants.append(I1) # 二阶不变量 I2 = (nu_20 - nu_02)**2 + 4 * nu_11**2 invariants.append(I2) # 高阶不变量 if n > 3: # 三阶中心矩 mu_30 = np.sum((x - eta_10)**3 * signal) mu_21 = np.sum((x - eta_10)**2 * (signal - eta_01) * signal) mu_12 = np.sum((x - eta_10) * (signal - eta_01)**2 * signal) mu_03 = np.sum((signal - eta_01)**3 * signal) # 归一化三阶矩 gamma_30 = mu_30 / (mu_00**2.5) if mu_00 > 0 else 0 gamma_21 = mu_21 / (mu_00**2.5) if mu_00 > 0 else 0 gamma_12 = mu_12 / (mu_00**2.5) if mu_00 > 0 else 0 gamma_03 = mu_03 / (mu_00**2.5) if mu_00 > 0 else 0 # 三阶不变量 I3 = gamma_30**2 + gamma_03**2 + 3 * (gamma_21**2 + gamma_12**2) invariants.append(I3) # 计算信号曲率特征（仿射不变性） curvature_features = self._compute_curvature_invariants(signal) invariants.extend(curvature_features) # 计算仿射弧长参数化特征 arc_length_features = self._compute_affine_arc_length(signal) invariants.extend(arc_length_features) return np.array(invariants[:self.n_invariants]) def _compute_curvature_invariants(self, signal): """计算曲率仿射不变量""" n = len(signal) if n < 3: return [0, 0, 0] # 一阶差分 dx = np.gradient(signal) # 二阶差分 d2x = np.gradient(dx) # 计算曲率（欧氏几何） curvature = np.abs(d2x) / (1 + dx**2)**1.5 # 计算仿射曲率（仿射几何） # 仿射曲率公式: κ_affine = (x'y'' - x''y') / (x'^2 + y'^2)^(3/2) # 这里y=signal, x=时间索引 t = np.arange(n) dt = np.gradient(t) d2t = np.gradient(dt) numerator = dt * d2x - d2t * dx denominator = (dt**2 + dx**2)**1.5 affine_curvature = numerator / (denominator + 1e-10) # 曲率统计特征 features = [ np.mean(curvature), np.std(curvature), np.max(curvature), np.mean(np.abs(affine_curvature)), np.std(affine_curvature), np.sum(curvature**2), # 曲率能量 np.sum(affine_curvature**2) # 仿射曲率能量 ] return features def _compute_affine_arc_length(self, signal): """计算仿射弧长参数化特征""" n = len(signal) if n < 2: return [0, 0, 0] # 时间参数 t = np.arange(n) # 计算仿射弧长微分: ds = (x'y'' - x''y')^(1/3) dt dx = np.gradient(signal) d2x = np.gradient(dx) dt = np.gradient(t) d2t = np.gradient(dt) # 仿射弧长微分 det = dt * d2x - d2t * dx ds = np.abs(det)**(1/3) # 总仿射弧长 total_arc_length = np.sum(ds) # 弧长参数化特征 features = [ total_arc_length, np.mean(ds), np.std(ds), np.max(ds) / (np.min(ds) + 1e-10), # 弧长变化率 np.sum(ds**2) # 弧长能量 ] return features参考文章：
基于仿射不变性-复流形结构-芬斯勒度量-子流形嵌入四维几何融合的机械故障诊断方法（Python） - 哥廷根数学学派的文章 -
https://zhuanlan.zhihu.com/p/1986377999415792360
工学博士，担任《Mechanical System and Signal Processing》审稿专家，担任
《中国电机工程学报》优秀审稿专家，《控制与决策》，《系统工程与电子技术》，《电力系统保护与控制》，《宇航学报》等EI期刊审稿专家。
擅长领域：现代信号处理，机器学习，深度学习，数字孪生，时间序列分析，设备缺陷检测、设备异常检测、设备智能故障诊断与健康管理PHM等。
