全髋关节置换手术的信号处理
全髋关节置换手术的信号处理
在医学音频信号研讨会的框架内，作为我们硕士学位的一部分，我们被要求解决一个与髋关节置换手术和从中记录的音频信号有关的研究问题。
在这篇文章中，详细描述了我们使用的数据、命中提取和预处理、交叉验证以及拟合最佳分类器以区分手术的首次命中和晚期命中之间的差异。
Taken from Google photos
全髋关节置换 (THR) 手术
在全髋关节置换术（也称为全髋关节置换术）中，受损的骨骼和软骨被移除并替换为假体组件。将受损的股骨头移除并替换为放置在股骨中空中心的金属杆。股骨柄可以胶合或压入骨中。
文献中描述了识别和准备受损软骨表面以及实现髋臼组件稳定固定的困难。因此，THR 由高级外科医生执行。
Total Hip Replacement (Taken from google photos)
我们的数据
我们有 31 种不同的 THR 手术的录音机，每个音频文件都包含手术的不同部分，并且充满了整个手术过程中的热门歌曲。
我们的研究
我们的自动命中查找器算法基于 Librosa 库。
数据，sr = librosa.load（文件名，sr=sample_rate）
D = np.abs(librosa.stft(数据，n_fft=1024，hop_length=512，window="hann"))
onset_env = librosa.onset.onset_strength(y=data, sr=sr, aggregate=np.max, S=D)
onset_frames = librosa.onset.onset_detect(y=data, onset_envelope=onset_env, sr=sr, normalize=True)
加载 .wav 文件，提取具有 1024 个样本长度的 n_fft 和 50% 重叠的跳跃长度的 STFT。从 STFT 的绝对值中，我们可以通过其最大幅度检测 .wav 文件中命中的峰值，我们使用起始（峰值）后 0.2 ms 的窗口长度来提取命中。
然后，我们手动验证了我们的命中查找器算法，通过在时域中绘制信号、具有峰值位置的平滑信号以及在 x 轴上按时间对齐的信号频谱图
Example of the automatic hits extraction algorithm
在组装包含 6051 个命中、按手术索引和记录编号的数据库后，我们希望提取每个命中的信号特征。在我们的研究中处理音频数据时，我们使用了 Librosa 库。对每个命中记录执行特征提取阶段，如下所示：
-
过零率—— 时域特征，即信号通过幅度轴的值 0（X 轴）的次数。我们在 64000Hz 的采样频率下使用 (sample size / 7) 个样本的窗口大小。返回总计 8个特点 每个命中记录。
-
节奏—— 时域特征，它是以每分钟节拍 (BPM) 为单位的音频的总体估计节奏，节奏是给定片段的速度或节奏，直接来自平均节拍持续时间。退货 1 个功能 每个命中记录。
-
日志 s— 频域特征，使用梅尔谱图。 Mel谱图是在频率尺度上通过非线性变换到Mel尺度的信号谱图，我们对音频数据的每帧使用1024个样本的窗口大小，其中512个样本重叠，并且在每帧上应用汉明窗，然后我们使用 DFT 将信号变换到频域。 对数 在 DFT 中应用函数将频谱图转换为 Mel **** 频谱图。频谱图从功率转换为分贝，我们在转换后的梅尔频谱图中选择每个频带的最大值。退货 128 项功能 每个命中记录。
-
多晶硅 — 频域特征，即梅尔频率倒谱系数是信号的紧凑频谱表示。我们使用了 20 个 mfcc 频段，并从每个频段中选择了最大值。退货 20 项功能 对于每个命中记录。
-
Mfcc 增量 — 频域特征，它是 mfcc 数据沿 mfcc 频带的一阶导数的局部估计。为了计算每个导数，我们使用了 mfcc 数据中的 3 个波段。退货 20 项功能 对于每个命中记录。
-
Mfcc 第二增量 — 频域特征，它是 mfcc 数据沿 mfcc 频带的二阶导数的局部估计。为了计算每个导数，我们使用了 mfcc 数据中的 3 个波段。退货 20 项功能 对于每个命中记录。
-
色度 — 频域特征，它使用短期傅里叶变换 (STFT) 来计算色度特征。我们选择了 1024 个样本的窗口大小，其中 512 个样本重叠。音频的色度值代表十二个不同音级的强度。退货 12 项功能 对于每个命中记录。
-
快速傅里叶变换 — 频域特征，即快速傅里叶变换数据。退货 6400 个功能 每个命中记录。
我们总共有 209 个信号特征和 6400 个 FFT 特征，每个命中记录总共有 6609 个特征。
我们必须先标记数据，然后才能对其进行分类。我们总共有 6051 个命中样本，每个命中有 6609 个特征列。之后，我们将只使用每次手术的第一次和最后一次拉刀击打，总共 1452 次击打。为了为我们的研究问题提供最佳答案，我们将第一个拉刀的前 100% 命中标记为 0（早期命中），最后一个拉刀的最后 50% 命中标记为 1。（后期命中）。标记阶段后我们有 912 次点击，其中 366 次后期点击和 546 次早期点击（比率为 40.1%）。
Samples and features downsampling
为了确保模型接收到的数据的质量，我们联系了 Mann-Whitney 测试，以了解早期命中与晚期命中中每列特征的分布。 Mann-Whitney 检验的原假设是组的分布是相同的，所以对于一个好的特征，我们想拒绝检验的原假设。上述测试以 95% 的置信度评估了所有特征。我们删除了所有无法使用测试拒绝零假设的特征。测试中有 3688 个特征在早期和晚期命中值之间没有显着差异，因此我们将它们从数据中删除。
此外，我们删除了与另一个特征相关性大于 90% 且方差小于 1% 的所有特征。我们删除了 2370 个相关列和 0 个方差小于 1% 的特征。
最后，数据有 554 列（特征）和 912 行（命中）。
拟合最佳分类模型
二元分类问题是我们研究的重点。为了防止训练、验证和测试之间的数据泄漏，我们将数据按手术 id 而不是命中 id 划分，并选择了四个手术作为测试数据。共记录134次命中，其中晚命中60次，早命中74次（比例为44.77%）。
对于剩余的手术，我们使用了 7-fold 交叉验证，这意味着每个折叠由 4 个手术组成（除了 1 个具有 3 个手术的手术），其余的 23/24 手术用于训练。
在我们的研究中，我们选择比较以下分类模型：
-
逻辑回归 — 惩罚弹性网络，0.5 l1 比率，正则化强度的倒数等于 150，类权重平衡
-
随机森林分类器—— 具有 1000 个估计器，最大深度为 7，每次拆分最少 3 个样本，创建新节点的最小杂质为 0.05，并且类权重平衡。
-
XGboost 分类器 — 学习率 0.5，最小损失减少 1，最大深度 12，基于梯度的采样方法，类权重平衡。
-
梯度提升分类器- 具有指数损失，0.2 学习率，5 个估计器，每个估计器只能使用 20% 的数据，每次拆分最少 3 个样本，最大特征日志以 2 为底，创建新节点的最小杂质为 0.05，最大深度为 4。
结果
对于 4 个选定模型中的每一个，我们计算了交叉验证中每个折叠的训练 AUC 和验证 AUC，并计算了所有模型的均值和标准差。
Example of the cross-validation for the random forest
我们假设相关方差是由从我们的数据中提取的“嘈杂”（手术室中不同机器的背景噪声和手术期间的谈话）命中引起的。
Confusion matrix of the test set after using the random forest classifier
Test scores of the fitted models
Histogram of the most significant feature by the random forest classifier for each of the classes
结论和未来工作
我们的研究问题是区分第一个拉刀的迟到和最后一个拉刀的迟到之间的区别。为了获得任何差异，我们必须构建一个数据管道，帮助我们理解和探索我们拥有的独特数据集。我们创建了一个分类模型，可以通过 80% AUC 区分早期命中和晚期命中，并且从它的学习中，我们从每个模型中提取了最有贡献的特征，并发现一些模型依赖于相同的特征或特征组（例如， 1130 和 1140 Hz）来提供预测。表 1 显示 40% 的特征在 1100 Hz 到 1300 Hz 的频率范围内，这些幅度的值越高，就越有可能是后期命中。例如，特征 mfcc delta 2 有助于 xgboost 模型和逻辑回归模型，其中较高的值表示早期命中的可能性较大。
对于未来的工作，我们建议在学习阶段之前将信号从噪声中清除。
所有相关代码都可以在 GitHub .
这项工作的成员——
布赖恩·里克什彭 — 领英
Guy Choen — 领英
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明
