当前位置：首页 > web >正文

基于多通道同步分析的智能听诊系统应用程序

web 2025/8/26 17:23:23

基于多通道同步分析的智能听诊系统应用程序

1. 项目概述

1.1 项目背景

听诊是医疗诊断中的重要手段，传统听诊方法依赖医生的经验和听觉判断，存在主观性强、难以量化分析等局限性。随着数字信号处理技术和人工智能的发展，智能听诊系统成为医疗科技领域的研究热点。本项目旨在开发一个能够处理16路时间同步音频数据的智能听诊系统，通过对多通道音频数据的同步采集、处理和分析，生成深度医疗报告，辅助医生进行更精准的诊断。

1.2 系统目标

实现16路音频数据的同步采集和存储
开发多通道音频信号处理和分析算法
构建基于机器学习的病理声音分类模型
生成结构化的深度分析报告
提供直观的用户界面和可视化功能

1.3 技术架构

系统采用分层架构设计，包括数据采集层、信号处理层、特征提取层、模型分析层和应用展示层。使用Python作为主要开发语言，结合多种科学计算和机器学习库实现核心功能。

2. 系统设计与实现

2.1 系统架构设计

# 系统主要模块定义
class IntelligentStethoscopeSystem:def __init__(self):self.data_acquisition = DataAcquisitionModule()self.signal_processing = SignalProcessingModule()self.feature_extraction = FeatureExtractionModule()self.ml_models = MLModelsModule()self.report_generation = ReportGenerationModule()self.visualization = VisualizationModule()def run_pipeline(self, audio_data=None):"""运行完整的处理流程"""if audio_data is None:audio_data = self.data_acquisition.capture_audio()processed_data = self.signal_processing.process(audio_data)features = self.feature_extraction.extract(processed_data)analysis_results = self.ml_models.analyze(features)report = self.report_generation.generate(analysis_results)return report, analysis_results

2.2 数据采集模块

import numpy as np
import sounddevice as sd
import threading
import time
from datetime import datetimeclass DataAcquisitionModule:def __init__(self, sample_rate=44100, channels=16, chunk_size=1024):self.sample_rate = sample_rateself.channels = channelsself.chunk_size = chunk_sizeself.is_recording = Falseself.audio_data = np.array([])self.stream = Noneself.lock = threading.Lock()def _callback(self, indata, frames, time, status):"""音频流回调函数"""if status:print(f"Stream status: {status}")with self.lock:if self.audio_data.size == 0:self.audio_data = indata.copy()else:self.audio_data = np.vstack((self.audio_data, indata))def capture_audio(self, duration=10):"""采集指定时长的音频数据"""print(f"开始采集音频，时长: {duration}秒")try:# 设置音频流参数self.audio_data = np.array([])self.is_recording = True# 创建输入流self.stream = sd.InputStream(samplerate=self.sample_rate,channels=self.channels,blocksize=self.chunk_size,callback=self._callback)# 开始录制self.stream.start()time.sleep(duration)self.stream.stop()self.stream.close()self.is_recording = Falseprint("音频采集完成")return self.audio_dataexcept Exception as e:print(f"音频采集错误: {e}")self.is_recording = Falseif self.stream:self.stream.close()return Nonedef save_audio(self, filename=None):"""保存音频数据到文件"""if filename is None:timestamp = datetime.now().strftime("%Y%m%d_%H%M%S")filename = f"audio_data_{timestamp}.npy"with self.lock:np.save(filename, self.audio_data)print(f"音频数据已保存到: {filename}")return filenamedef load_audio(self, filename):"""从文件加载音频数据"""try:self.audio_data = np.load(filename)print(f"音频数据已从 {filename} 加载")return self.audio_dataexcept Exception as e:print(f"加载音频数据错误: {e}")return None

2.3 信号处理模块

import numpy as np
from scipy import signal
import pywtclass SignalProcessingModule:def __init__(self):self.sample_rate = 44100  # 默认采样率def set_sample_rate(self, sample_rate):self.sample_rate = sample_ratedef bandpass_filter(self, data, lowcut=20, highcut=1000, order=4):"""带通滤波器"""nyquist = 0.5 * self.sample_ratelow = lowcut / nyquisthigh = highcut / nyquistb, a = signal.butter(order, [low, high], btype='band')filtered_data = signal.filtfilt(b, a, data, axis=0)return filtered_datadef notch_filter(self, data, freq=50, q=30):"""陷波滤波器，去除电源干扰"""b, a = signal.iirnotch(freq, q, self.sample_rate)filtered_data = signal.filtfilt(b, a, data, axis=0)return filtered_datadef wavelet_denoise(self, data, wavelet='db4', level=3):"""小波去噪"""# 对每个通道进行处理denoised_data = np.zeros_like(data)for ch in range(data.shape[1]):coeffs = pywt.wavedec(data[:, ch], wavelet, level=level)# 计算阈值sigma = np.median(np.abs(coeffs[-level])) / 0.6745uthresh = sigma * np.sqrt(2 * np.log(len(data[:, ch])))# 应用软阈值coeffs = [pywt.threshold(c, uthresh, mode='soft') for c in coeffs]denoised_data[:, ch] = pywt.waverec(coeffs, wavelet)return denoised_datadef normalize(self, data):"""归一化处理"""max_vals = np.max(np.abs(data), axis=0)normalized_data = data / max_valsreturn normalized_datadef synchronize_channels(self, data, ref_channel=0):"""通道同步，基于互相关计算延迟"""synchronized_data = np.zeros_like(data)synchronized_data[:, ref_channel] = data[:, ref_channel]for ch in range(data.shape[1]):if ch == ref_channel:continue# 计算互相关correlation = signal.correlate(data[:, ref_channel], data[:, ch], mode='full')lags = signal.correlation_lags(len(data[:, ref_channel]), len(data[:, ch]), mode='full')delay = lags[np.argmax(correlation)]# 调整延迟if delay > 0:synchronized_data[delay:, ch] = data[:-delay, ch]elif delay < 0:synchronized_data[:delay, ch] = data[-delay:, ch]else:synchronized_data[:, ch] = data[:, ch]return synchronized_datadef segment_audio(self, data, segment_duration=1.0):"""分割音频为固定时长的片段"""segment_length = int(segment_duration * self.sample_rate)num_segments = data.shape[0] // segment_lengthsegments = []for i in range(num_segments):start = i * segment_lengthend = start + segment_lengthsegment = data[start:end, :]segments.append(segment)return segmentsdef process(self, audio_data, apply_filters=True, sync_channels=True):"""完整的信号处理流程"""print("开始信号处理...")processed_data = audio_data.copy()# 通道同步if sync_channels:processed_data = self.synchronize_channels(processed_data)# 应用滤波器if apply_filters:processed_data = self.notch_filter(processed_data, freq=50)  # 去除50Hz电源干扰processed_data = self.bandpass_filter(processed_data, lowcut=20, highcut=1500)  # 心音频率范围processed_data = self.wavelet_denoise(processed_data)  # 小波去噪# 归一化processed_data = self.normalize(processed_data)print("信号处理完成")return processed_data

2.4 特征提取模块

import numpy as np
from scipy import stats
from scipy.fft import fft, fftfreq
import librosa
import pywtclass FeatureExtractionModule:def __init__(self, sample_rate=44100):self.sample_rate = sample_ratedef set_sample_rate(self, sample_rate):self.sample_rate = sample_ratedef time_domain_features(self, data):"""时域特征提取"""features = {}# 基本统计特征features['mean'] = np.mean(data, axis=0)features['std'] = np.std(data, axis=0)features['max'] = np.max(data, axis=0)features['min'] = np.min(data, axis=0)features['range'] = features['max'] - features['min']features['rms'] = np.sqrt(np.mean(data**2, axis=0))# 高阶统计特征features['skewness'] = stats.skew(data, axis=0)features['kurtosis'] = stats.kurtosis(data, axis=0)# 过零率features['zero_crossing_rate'] = np.mean(np.diff(np.sign(data)) != 0, axis=0)return featuresdef frequency_domain_features(self, data):"""频域特征提取"""features = {}# 计算FFTn = len(data)fft_vals = np.abs(fft(data, axis=0))[:n//2]freqs = fftfreq(n, 1/self.sample_rate)[:n//2]# 频谱特征features['spectral_centroid'] = np.sum(freqs.reshape(-1, 1) * fft_vals, axis=0) / np.sum(fft_vals, axis=0)features['spectral_bandwidth'] = np.sqrt(np.sum((freqs.reshape(-1, 1) - features['spectral_centroid'])**2 * fft_vals, axis=0) / np.sum(fft_vals, axis=0))features['spectral_rolloff'] = self._calculate_rolloff(fft_vals, freqs)# 频带能量bands = [(20, 50), (50, 100), (100, 200), (200, 400), (400, 800), (800, 1600)]for i, (low, high) in enumerate(bands):band_mask = (freqs >= low) & (freqs <= high)features[f'energy_{low}_{high}'] = np.sum(fft_vals[band_mask, :], axis=0)# 频谱熵spectral_power = fft_vals / np.sum(fft_vals, axis=0)spectral_entropy = -np.sum(spectral_power * np.log2(spectral_power + 1e-10), axis=0)features['spectral_entropy'] = spectral_entropyreturn featuresdef _calculate_rolloff(self, spectrum, freqs, percentile=85):"""计算频谱滚降点"""rolloff_points = np.zeros(spectrum.shape[1])for i in range(spectrum.shape[1]):total_energy = np.sum(spectrum[:, i])cumulative_energy = 0for j in range(len(spectrum)):cumulative_energy += spectrum[j, i]if cumulative_energy >= (percentile / 100) * total_energy:rolloff_points[i] = freqs[j]breakreturn rolloff_pointsdef wavelet_features(self, data, wavelet='db4', level=5):"""小波变换特征提取"""features = {}for ch in range(data.shape[1]):# 小波分解coeffs = pywt.wavedec(data[:, ch], wavelet, level=level)# 计算各层小波系数的能量for i, coeff in enumerate(coeffs):features[f'wavelet_energy_l{i}_ch{ch}'] = np.sum(coeff**2)# 计算各层小波系数的熵for i, coeff in enumerate(coeffs):energy = np.sum(coeff**2)if energy > 0:norm_coeff = coeff**2 / energyentropy = -np.sum(norm_coeff * np.log2(norm_coeff + 1e-10))features[f'wavelet_entropy_l{i}_ch{ch}'] = entropyreturn featuresdef mfcc_features(self, data, n_mfcc=13):"""MFCC特征提取"""mfcc_features = []for ch in range(data.shape[1]):mfcc = librosa.feature.mfcc(y=data[:, ch], sr=self.sample_rate, n_mfcc=n_mfcc)mfcc_mean = np.mean(mfcc, axis=1)mfcc_std = np.std(mfcc, axis=1)for i in range(n_mfcc):mfcc_features.append(mfcc_mean[i])mfcc_features.append(mfcc_std[i])return np.array(mfcc_features)def extract(self, processed_data, include_wavelet=True, include_mfcc=True):"""提取所有特征"""print("开始特征提取...")features = {}# 时域特征time_features = self.time_domain_features(processed_data)features.update(time_features)# 频域特征freq_features = self.frequency_domain_features(processed_data)features.update(freq_features)# 小波特征if include_wavelet:wavelet_features = self.wavelet_features(processed_data)features.update(wavelet_features)# MFCC特征if include_mfcc:for ch in range(processed_data.shape[1]):mfcc_feat = self.mfcc_features(processed_data[:, ch].reshape(-1, 1))for i, val in enumerate(mfcc_feat):features[f'mfcc_{i}_ch{ch}'] = valprint("特征提取完成")return features

2.5 机器学习模型模块

import numpy as np
import joblib
from sklearn.ensemble import RandomForestClassifier, IsolationForest
from sklearn.svm import SVC
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score, classification_report, confusion_matrix
import tensorflow as tf
from tensorflow import keras
from tensorflow.keras import layersclass MLModelsModule:def __init__(self):self.models = {}self.scalers = {}self.history = {}def train_random_forest(self, X, y, model_name='random_forest', test_size=0.2):"""训练随机森林分类器"""print(f"训练随机森林模型: {model_name}")# 数据预处理X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=test_size, random_state=42)# 标准化scaler = StandardScaler()X_train_scaled = scaler.fit_transform(X_train)X_test_scaled = scaler.transform(X_test)self.scalers[model_name] = scaler# 训练模型model = RandomForestClassifier(n_estimators=100, random_state=42)model.fit(X_train_scaled, y_train)# 评估模型y_pred = model.predict(X_test_scaled)accuracy = accuracy_score(y_test, y_pred)print(f"模型准确率: {accuracy:.4f}")print("分类报告:")print(classification_report(y_test, y_pred))self.models[model_name] = modelreturn model, accuracydef train_svm(self, X, y, model_name='svm', test_size=0.2):"""训练支持向量机分类器"""print(f"训练SVM模型: {model_name}")# 数据预处理X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=test_size, random_state=42)# 标准化scaler = StandardScaler()X_train_scaled = scaler.fit_transform(X_train)X_test_scaled = scaler.transform(X_test)self.scalers[model_name] = scaler# 训练模型model = SVC(kernel='rbf', random_state=42, probability=True)model.fit(X_train_scaled, y_train)# 评估模型y_pred = model.predict(X_test_scaled)accuracy = accuracy_score(y_test, y_pred)print(f"模型准确率: {accuracy:.4f}")print("分类报告:")print(classification_report(y_test, y_pred))self.models[model_name] = modelreturn model, accuracydef build_cnn_model(self, input_shape, num_classes):"""构建CNN模型"""model = keras.Sequential([layers.Conv1D(32, kernel_size=3, activation='relu', input_shape=input_shape),layers.MaxPooling1D(pool_size=2),layers.Conv1D(64, kernel_size=3, activation='relu'),layers.MaxPooling1D(pool_size=2),layers.Conv1D(128, kernel_size=3, activation='relu'),layers.MaxPooling1D(pool_size=2),layers.Flatten(),layers.Dense(128, activation='relu'),layers.Dropout(0.5),layers.Dense(num_classes, activation='softmax')])model.compile(optimizer='adam',loss='sparse_categorical_crossentropy',metrics=['accuracy'])return modeldef train_cnn(self, X, y, model_name='cnn', test_size=0.2, epochs=50, batch_size=32):"""训练CNN模型"""print(f"训练CNN模型: {model_name}")# 数据预处理X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=test_size, random_state=42)# 确保输入是3D格式 (samples, timesteps, channels)if len(X_train.shape) == 2:X_train = np.expand_dims(X_train, axis=2)X_test = np.expand_dims(X_test, axis=2)# 构建模型input_shape = (X_train.shape[1], X_train.shape[2])num_classes = len(np.unique(y))model = self.build_cnn_model(input_shape, num_classes)# 训练模型history = model.fit(X_train, y_train,batch_size=batch_size,epochs=epochs,validation_data=(X_test, y_test),verbose=1)# 评估模型test_loss, test_accuracy = model.evaluate(X_test, y_test, verbose=0)print(f"测试准确率: {test_accuracy:.4f}")self.models[model_name] = modelself.history[model_name] = historyreturn model, test_accuracydef train_anomaly_detection(self, X, model_name='anomaly_detection', contamination=0.1):"""训练异常检测模型"""print(f"训练异常检测模型: {model_name}")# 标准化scaler = StandardScaler()X_scaled = scaler.fit_transform(X)self.scalers[model_name] = scaler# 训练隔离森林模型model = IsolationForest(contamination=contamination, random_state=42)model.fit(X_scaled)self.models[model_name] = modelreturn modeldef predict(self, X, model_name):"""使用指定模型进行预测"""if model_name not in self.models:raise ValueError(f"模型 {model_name} 不存在")model = self.models[model_name]scaler = self.scalers.get(model_name, None)if scaler:X = scaler.transform(X)# 根据模型类型进行预测if isinstance(model, IsolationForest):return model.predict(X)else:return model.predict_proba(X) if hasattr(model, 'predict_proba') else model.predict(X)def save_model(self, model_name, filename):"""保存模型到文件"""if model_name not in self.models:raise ValueError(f"模型 {model_name} 不存在")joblib.dump({'model': self.models[model_name],'scaler': self.scalers.get(model_name, None)}, filename)print(f"模型已保存到: {filename}")def load_model(self, model_name, filename):"""从文件加载模型"""data = joblib.load(filename)self.models[model_name] = data['model']if data['scaler']:self.scalers[model_name] = data['scaler']print(f"模型已从 {filename} 加载")def analyze(self, features):"""使用所有模型进行分析"""print("开始机器学习分析...")# 将特征字典转换为特征向量feature_vector = self._features_to_vector(features)feature_vector = feature_vector.reshape(1, -1)results = {}# 使用每个模型进行预测for model_name in self.models:try:prediction = self.predict(feature_vector, model_name)results[model_name] = predictionexcept Exception as e:print(f"模型 {model_name} 预测错误: {e}")results[model_name] = Noneprint("机器学习分析完成")return resultsdef _features_to_vector(self, features):"""将特征字典转换为特征向量"""vector = []for key in sorted(features.keys()):value = features[key]if np.isscalar(value):vector.append(value)else:vector.extend(value.flatten())return np.array(vector)

2.6 报告生成模块

import json
from datetime import datetime
import matplotlib.pyplot as plt
import numpy as npclass ReportGenerationModule:def __init__(self):self.template = {"report_id": "","timestamp": "","patient_info": {},"acquisition_parameters": {},"signal_quality_metrics": {},"feature_summary": {},"analysis_results": {},"diagnostic_impressions": {},"recommendations": []}def generate(self, analysis_results, patient_info=None, acquisition_params=None):"""生成完整报告"""print("生成分析报告...")report = self.template.copy()report["report_id"] = f"RPT_{datetime.now().strftime('%Y%m%d_%H%M%S')}"report["timestamp"] = datetime.now().isoformat()# 添加患者信息if patient_info:report["patient_info"] = patient_infoelse:report["patient_info"] = {"patient_id": "未知","age": "未知","gender": "未知","recording_location": "未知"}# 添加采集参数if acquisition_params:report["acquisition_parameters"] = acquisition_paramselse:report["acquisition_parameters"] = {"sample_rate": 44100,"channels": 16,"duration": "未知"}# 添加分析结果report["analysis_results"] = analysis_results# 生成诊断印象report["diagnostic_impressions"] = self._generate_diagnostic_impressions(analysis_results)# 生成建议report["recommendations"] = self._generate_recommendations(analysis_results)print("报告生成完成")return reportdef _generate_diagnostic_impressions(self, analysis_results):"""生成诊断印象"""impressions = {}# 解析随机森林模型结果if 'random_forest' in analysis_results and analysis_results['random_forest'] is not None:rf_probs = analysis_results['random_forest'][0]impressions['cardiac_abnormality_probability'] = float(rf_probs[1]) if len(rf_probs) > 1 else float(rf_probs[0])# 解析SVM模型结果if 'svm' in analysis_results and analysis_results['svm'] is not None:svm_probs = analysis_results['svm'][0]impressions['respiratory_abnormality_probability'] = float(svm_probs[1]) if len(svm_probs) > 1 else float(svm_probs[0])# 解析异常检测结果if 'anomaly_detection' in analysis_results and analysis_results['anomaly_detection'] is not None:anomaly_score = analysis_results['anomaly_detection'][0]impressions['anomaly_detected'] = True if anomaly_score == -1 else Falsereturn impressionsdef _generate_recommendations(self, analysis_results):"""生成建议"""recommendations = []# 基于随机森林模型结果if 'random_forest' in analysis_results and analysis_results['random_forest'] is not None:rf_probs = analysis_results['random_forest'][0]cardiac_prob = rf_probs[1] if len(rf_probs) > 1 else rf_probs[0]if cardiac_prob > 0.7:recommendations.append("检测到心脏异常的高可能性，建议进行心电图和心脏超声检查")elif cardiac_prob > 0.4:recommendations.append("检测到心脏异常的中等可能性，建议进行定期监测和随访")else:recommendations.append("心脏功能正常，建议保持健康生活方式")# 基于SVM模型结果if 'svm' in analysis_results and analysis_results['svm'] is not None:svm_probs = analysis_results['svm'][0]respiratory_prob = svm_probs[1] if len(svm_probs) > 1 else svm_probs[0]if respiratory_prob > 0.7:recommendations.append("检测到呼吸异常的高可能性，建议进行胸部X光或CT检查")elif respiratory_prob > 0.4:recommendations.append("检测到呼吸异常的中等可能性，建议进行肺功能测试")# 基于异常检测结果if 'anomaly_detection' in analysis_results and analysis_results['anomaly_detection'] is not None:anomaly_score = analysis_results['anomaly_detection'][0]if anomaly_score == -1:recommendations.append("检测到异常信号模式，建议进行进一步临床评估")if not recommendations:recommendations.append("未检测到明显异常，建议定期健康检查")return recommendationsdef save_report(self, report, filename=None, format='json'):"""保存报告到文件"""if filename is None:timestamp = datetime.now().strftime("%Y%m%d_%H%M%S")filename = f"stethoscope_report_{timestamp}.{format}"if format == 'json':with open(filename, 'w') as f:json.dump(report, f, indent=4)elif format == 'txt':with open(filename, 'w') as f:f.write(self._format_text_report(report))else:raise ValueError("不支持的格式，请使用 'json' 或 'txt'")print(f"报告已保存到: {filename}")return filenamedef _format_text_report(self, report):"""将报告格式化为文本"""text = f"智能听诊系统分析报告\n"text += f"报告ID: {report['report_id']}\n"text += f"生成时间: {report['timestamp']}\n\n"text += "患者信息:\n"for key, value in report['patient_info'].items():text += f"  {key}: {value}\n"text += "\n采集参数:\n"for key, value in report['acquisition_parameters'].items():text += f"  {key}: {value}\n"text += "\n诊断印象:\n"for key, value in report['diagnostic_impressions'].items():text += f"  {key}: {value}\n"text += "\n建议:\n"for i, rec in enumerate(report['recommendations'], 1):text += f"  {i}. {rec}\n"return textdef generate_visual_report(self, report, audio_data, features, filename=None):"""生成可视化报告"""if filename is None:timestamp = datetime.now().strftime("%Y%m%d_%H%M%S")filename = f"visual_report_{timestamp}.png"fig, axes = plt.subplots(3, 2, figsize=(15, 12))fig.suptitle('智能听诊系统可视化报告', fontsize=16)# 绘制音频波形axes[0, 0].set_title('音频波形 (前4个通道)')for i in range(min(4, audio_data.shape[1])):axes[0, 0].plot(audio_data[:1000, i] + i * 0.5, label=f'通道 {i+1}')axes[0, 0].legend()axes[0, 0].set_xlabel('采样点')axes[0, 0].set_ylabel('振幅')# 绘制频谱axes[0, 1].set_title('频谱图 (通道1)')n = len(audio_data)fft_vals = np.abs(np.fft.fft(audio_data[:, 0]))[:n//2]freqs = np.fft.fftfreq(n, 1/44100)[:n//2]axes[0, 1].plot(freqs, fft_vals)axes[0, 1].set_xlabel('频率 (Hz)')axes[0, 1].set_ylabel('振幅')axes[0, 1].set_xlim(0, 2000)# 绘制特征分布axes[1, 0].set_title('时域特征分布')time_features = {k: v for k, v in features.items() if any(t in k for t in ['mean', 'std', 'rms'])}axes[1, 0].bar(range(len(time_features)), list(time_features.values())[:10])axes[1, 0].set_xticks(range(len(time_features))[:10])axes[1, 0].set_xticklabels(list(time_features.keys())[:10], rotation=45)axes[1, 0].set_ylabel('值')# 绘制频域特征axes[1, 1].set_title('频域特征分布')freq_features = {k: v for k, v in features.items() if any(t in k for t in ['energy', 'spectral'])}axes[1, 1].bar(range(len(freq_features)), list(freq_features.values())[:10])axes[1, 1].set_xticks(range(len(freq_features))[:10])axes[1, 1].set_xticklabels(list(freq_features.keys())[:10], rotation=45)axes[1, 1].set_ylabel('值')# 绘制诊断概率axes[2, 0].set_title('异常概率')diagnoses = []probabilities = []if 'cardiac_abnormality_probability' in report['diagnostic_impressions']:diagnoses.append('心脏异常')probabilities.append(report['diagnostic_impressions']['cardiac_abnormality_probability'])if 'respiratory_abnormality_probability' in report['diagnostic_impressions']:diagnoses.append('呼吸异常')probabilities.append(report['diagnostic_impressions']['respiratory_abnormality_probability'])if diagnoses:axes[2, 0].bar(diagnoses, probabilities)axes[2, 0].set_ylabel('概率')axes[2, 0].set_ylim(0, 1)# 添加文本报告axes[2, 1].set_title('报告摘要')text_content = f"报告ID: {report['report_id']}\n\n"text_content += "诊断印象:\n"for k, v in report['diagnostic_impressions'].items():text_content += f"{k}: {v:.3f}\n" if isinstance(v, float) else f"{k}: {v}\n"text_content += "\n建议:\n"for i, rec in enumerate(report['recommendations'], 1):text_content += f"{i}. {rec}\n"axes[2, 1].text(0.05, 0.95, text_content, transform=axes[2, 1].transAxes, verticalalignment='top', fontsize=10)axes[2, 1].axis('off')plt.tight_layout()plt.savefig(filename, dpi=300, bbox_inches='tight')plt.close()print(f"可视化报告已保存到: {filename}")return filename

2.7 可视化模块

import matplotlib.pyplot as plt
import numpy as np
from matplotlib.gridspec import GridSpecclass VisualizationModule:def __init__(self):self.fig_size = (15, 10)def plot_waveforms(self, audio_data, channels=None, title="音频波形"):"""绘制音频波形图"""if channels is None:channels = range(min(4, audio_data.shape[1]))plt.figure(figsize=self.fig_size)for i, ch in enumerate(channels):plt.plot(audio_data[:2000, ch] + i * 0.5, label=f'通道 {ch+1}')plt.title(title)plt.xlabel('采样点')plt.ylabel('振幅')plt.legend()plt.grid(True)plt.show()def plot_spectrograms(self, audio_data, channels=None, sample_rate=44100):"""绘制频谱图"""if channels is None:channels = range(min(4, audio_data.shape[1]))fig, axes = plt.subplots(2, 2, figsize=self.fig_size)axes = axes.flatten()for i, ch in enumerate(channels):if i >= len(axes):breakn = len(audio_data[:, ch])fft_vals = np.abs(np.fft.fft(audio_data[:, ch]))[:n//2]freqs = np.fft.fftfreq(n, 1/sample_rate)[:n//2]axes[i].plot(freqs, fft_vals)axes[i].set_title(f'通道 {ch+1} 频谱')axes[i].set_xlabel('频率 (Hz)')axes[i].set_ylabel('振幅')axes[i].set_xlim(0, 2000)axes[i].grid(True)plt.tight_layout()plt.show()def plot_feature_importance(self, feature_names, importance_scores, top_n=20):"""绘制特征重要性图"""# 排序特征重要性indices = np.argsort(importance_scores)[-top_n:]plt.figure(figsize=(12, 8))plt.barh(range(top_n), importance_scores[indices])plt.yticks(range(top_n), [feature_names[i] for i in indices])plt.xlabel('特征重要性')plt.title('Top 特征重要性')plt.tight_layout()plt.show()def plot_confidence_scores(self, analysis_results):"""绘制置信度分数图"""model_names = []confidence_scores = []for model_name, result in analysis_results.items():if result is not None and hasattr(result, 'shape') and result.shape[0] > 0:if result.shape[1] > 1:  # 多分类问题confidence = np.max(result[0])else:  # 二分类或异常检测confidence = np.abs(result[0][0])model_names.append(model_name)confidence_scores.append(confidence)plt.figure(figsize=(10, 6))plt.bar(model_names, confidence_scores)plt.xlabel('模型')plt.ylabel('置信度')plt.title('模型预测置信度')plt.ylim(0, 1)plt.xticks(rotation=45)plt.tight_layout()plt.show()def create_dashboard(self, audio_data, features, analysis_results, report):"""创建综合仪表板"""fig = plt.figure(figsize=(20, 15))gs = GridSpec(3, 3, figure=fig)# 波形图ax1 = fig.add_subplot(gs[0, 0])for i in range(min(4, audio_data.shape[1])):ax1.plot(audio_data[:1000, i] + i * 0.5, label=f'通道 {i+1}')ax1.set_title('音频波形')ax1.legend()# 频谱图ax2 = fig.add_subplot(gs[0, 1])n = len(audio_data)fft_vals = np.abs(np.fft.fft(audio_data[:, 0]))[:n//2]freqs = np.fft.fftfreq(n, 1/44100)[:n//2]ax2.plot(freqs, fft_vals)ax2.set_title('通道1频谱')ax2.set_xlim(0, 2000)# 特征分布ax3 = fig.add_subplot(gs[0, 2])time_features = {k: v for k, v in features.items() if any(t in k for t in ['mean', 'std', 'rms'])}ax3.bar(range(len(time_features)), list(time_features.values())[:10])ax3.set_xticks(range(len(time_features))[:10])ax3.set_xticklabels(list(time_features.keys())[:10], rotation=45)ax3.set_title('时域特征')# 频域特征ax4 = fig.add_subplot(gs[1, 0])freq_features = {k: v for k, v in features.items() if any(t in k for t in ['energy', 'spectral'])}ax4.bar(range(len(freq_features)), list(freq_features.values())[:10])ax4.set_xticks(range(len(freq_features))[:10])ax4.set_xticklabels(list(freq_features.keys())[:10], rotation=45)ax4.set_title('频域特征')# 诊断概率ax5 = fig.add_subplot(gs[1, 1])diagnoses = []probabilities = []if 'cardiac_abnormality_probability' in report['diagnostic_impressions']:diagnoses.append('心脏异常')probabilities.append(report['diagnostic_impressions']['cardiac_abnormality_probability'])if 'respiratory_abnormality_probability' in report['diagnostic_impressions']:diagnoses.append('呼吸异常')probabilities.append(report['diagnostic_impressions']['respiratory_abnormality_probability'])if diagnoses:ax5.bar(diagnoses, probabilities)ax5.set_ylim(0, 1)ax5.set_title('异常概率')# 模型置信度ax6 = fig.add_subplot(gs[1, 2])model_names = []confidence_scores = []for model_name, result in analysis_results.items():if result is not None and hasattr(result, 'shape') and result.shape[0] > 0:if result.shape[1] > 1:  # 多分类问题confidence = np.max(result[0])else:  # 二分类或异常检测confidence = np.abs(result[0][0])model_names.append(model_name)confidence_scores.append(confidence)ax6.bar(model_names, confidence_scores)ax6.set_xticklabels(model_names, rotation=45)ax6.set_title('模型置信度')ax6.set_ylim(0, 1)# 报告文本ax7 = fig.add_subplot(gs[2, :])text_content = f"报告ID: {report['report_id']}\n\n"text_content += "诊断印象:\n"for k, v in report['diagnostic_impressions'].items():text_content += f"{k}: {v:.3f}\n" if isinstance(v, float) else f"{k}: {v}\n"text_content += "\n建议:\n"for i, rec in enumerate(report['recommendations'], 1):text_content += f"{i}. {rec}\n"ax7.text(0.05, 0.95, text_content, transform=ax7.transAxes, verticalalignment='top', fontsize=12)ax7.axis('off')ax7.set_title('报告摘要')plt.tight_layout()plt.show()

3. 系统集成与测试

3.1 主应用程序

import tkinter as tk
from tkinter import ttk, filedialog, messagebox
import threading
import matplotlib.pyplot as plt
from matplotlib.backends.backend_tkagg import FigureCanvasTkAgg
import numpy as np
from PIL import Image, ImageTkclass IntelligentStethoscopeApp:def __init__(self, root):self.root = rootself.root.title("基于多通道同步分析的智能听诊系统")self.root.geometry("1200x800")# 初始化系统self.system = IntelligentStethoscopeSystem()# 存储数据self.audio_data = Noneself.processed_data = Noneself.features = Noneself.analysis_results = Noneself.report = None# 创建界面self.create_widgets()def create_widgets(self):"""创建界面组件"""# 创建选项卡self.notebook = ttk.Notebook(self.root)self.notebook.pack(fill=tk.BOTH, expand=True, padx=10, pady=10)# 数据采集选项卡self.data_frame = ttk.Frame(self.notebook)self.notebook.add(self.data_frame, text="数据采集")self.create_data_acquisition_tab()# 信号处理选项卡self.processing_frame = ttk.Frame(self.notebook)self.notebook.add(self.processing_frame, text="信号处理")self.create_processing_tab()# 特征提取选项卡self.features_frame = ttk.Frame(self.notebook)self.notebook.add(self.features_frame, text="特征提取")self.create_features_tab()# 分析报告选项卡self.analysis_frame = ttk.Frame(self.notebook)self.notebook.add(self.analysis_frame, text="分析报告")self.create_analysis_tab()# 可视化选项卡self.visualization_frame = ttk.Frame(self.notebook)self.notebook.add(self.visualization_frame, text="可视化")self.create_visualization_tab()def create_data_acquisition_tab(self):"""创建数据采集选项卡"""# 控制框架control_frame = ttk.LabelFrame(self.data_frame, text="采集控制")control_frame.pack(fill=tk.X, padx=10, pady=5)ttk.Button(control_frame, text="开始采集", command=self.start_acquisition).pack(side=tk.LEFT, padx=5, pady=5)ttk.Button(control_frame, text="停止采集", command=self.stop_acquisition).pack(side=tk.LEFT, padx=5, pady=5)ttk.Button(control_frame, text="加载数据", command=self.load_audio_data).pack(side=tk.LEFT, padx=5, pady=5)ttk.Button(control_frame, text="保存数据", command=self.save_audio_data).pack(side=tk.LEFT, padx=5, pady=5)# 参数框架params_frame = ttk.LabelFrame(self.data_frame, text="采集参数")params_frame.pack(fill=tk.X, padx=10, pady=5)ttk.Label(params_frame, text="采样率:").grid(row=0, column=0, padx=5, pady=5, sticky=tk.W)self.sample_rate_var = tk.StringVar(value="44100")ttk.Entry(params_frame, textvariable=self.sample_rate_var, width=10).grid(row=0, column=1, padx=5, pady=5)ttk.Label(params_frame, text="通道数:").grid(row=0, column=2, padx=5, pady=5, sticky=tk.W)self.channels_var = tk.StringVar(value="16")ttk.Entry(params_frame, textvariable=self.channels_var, width=10).grid(row=0, column=3, padx=5, pady=5)ttk.Label(params_frame, text="采集时长(秒):").grid(row=0, column=4, padx=5, pady=5, sticky=tk.W)self.duration_var = tk.StringVar(value="10")ttk.Entry(params_frame, textvariable=self.duration_var, width=10).grid(row=0, column=5, padx=5, pady=5)# 状态框架status_frame = ttk.LabelFrame(self.data_frame, text="状态")status_frame.pack(fill=tk.X, padx=10, pady=5)self.status_var = tk.StringVar(value="就绪")ttk.Label(status_frame, textvariable=self.status_var).pack(padx=5, pady=5)# 波形显示框架waveform_frame = ttk.LabelFrame(self.data_frame, text="波形显示")waveform_frame.pack(fill=tk.BOTH, expand=True, padx=10, pady=5)self.waveform_fig, self.waveform_ax = plt.subplots(figsize=(10, 4))self.waveform_canvas = FigureCanvasTkAgg(self.waveform_fig, waveform_frame)self.waveform_canvas.get_tk_widget().pack(fill=tk.BOTH, expand=True)def create_processing_tab(self):"""创建信号处理选项卡"""# 控制框架control_frame = ttk.LabelFrame(self.processing_frame, text="处理控制")control_frame.pack(fill=tk.X, padx=10, pady=5)ttk.Button(control_frame, text="处理信号", command=self.process_signal).pack(side=tk.LEFT, padx=5, pady=5)# 参数框架params_frame = ttk.LabelFrame(self.processing_frame, text="处理参数")params_frame.pack(fill=tk.X, padx=10, pady=5)self.apply_filters_var = tk.BooleanVar(value=True)ttk.Checkbutton(params_frame, text="应用滤波器", variable=self.apply_filters_var).grid(row=0, column=0, padx=5, pady=5)self.sync_channels_var = tk.BooleanVar(value=True)ttk.Checkbutton(params_frame, text="同步通道", variable=self.sync_channels_var).grid(row=0, column=1, padx=5, pady=5)# 结果显示框架result_frame = ttk.LabelFrame(self.processing_frame, text="处理结果")result_frame.pack(fill=tk.BOTH, expand=True, padx=10, pady=5)# 创建左右分栏paned_window = ttk.PanedWindow(result_frame, orient=tk.HORIZONTAL)paned_window.pack(fill=tk.BOTH, expand=True)# 原始信号框架original_frame = ttk.Frame(paned_window)paned_window.add(original_frame, weight=1)ttk.Label(original_frame, text="原始信号").pack()self.original_fig, self.original_ax = plt.subplots(figsize=(6, 3))self.original_canvas = FigureCanvasTkAgg(self.original_fig, original_frame)self.original_canvas.get_tk_widget().pack(fill=tk.BOTH, expand=True)# 处理后的信号框架processed_frame = ttk.Frame(paned_window)paned_window.add(processed_frame, weight=1)ttk.Label(processed_frame, text="处理后的信号").pack()self.processed_fig, self.processed_ax = plt.subplots(figsize=(6, 3))self.processed_canvas = FigureCanvasTkAgg(self.processed_fig, processed_frame)self.processed_canvas.get_tk_widget().pack(fill=tk.BOTH, expand=True)# 其他选项卡的创建方法类似，限于篇幅不再详细展开def start_acquisition(self):"""开始采集音频"""try:sample_rate = int(self.sample_rate_var.get())channels = int(self.channels_var.get())duration = int(self.duration_var.get())self.system.data_acquisition.sample_rate = sample_rateself.system.data_acquisition.channels = channels# 在新线程中采集音频，避免界面冻结def acquisition_thread():self.status_var.set("采集中...")self.audio_data = self.system.data_acquisition.capture_audio(duration)self.status_var.set("采集完成")# 更新波形显示self.update_waveform_display()thread = threading.Thread(target=acquisition_thread)thread.daemon = Truethread.start()except ValueError:messagebox.showerror("错误", "请输入有效的参数值")def stop_acquisition(self):"""停止采集音频"""if self.system.data_acquisition.is_recording:self.system.data_acquisition.is_recording = Falseif self.system.data_acquisition.stream:self.system.data_acquisition.stream.stop()self.system.data_acquisition.stream.close()self.status_var.set("采集已停止")def load_audio_data(self):"""加载音频数据"""filename = filedialog.askopenfilename(filetypes=[("NumPy files", "*.npy")])if filename:self.audio_data = self.system.data_acquisition.load_audio(filename)if self.audio_data is not None:self.update_waveform_display()self.status_var.set(f"已加载: {filename}")def save_audio_data(self):"""保存音频数据"""if self.audio_data is not None:filename = filedialog.asksaveasfilename(defaultextension=".npy", filetypes=[("NumPy files", "*.npy")])if filename:self.system.data_acquisition.save_audio(filename)self.status_var.set(f"已保存: {filename}")else:messagebox.showwarning("警告", "没有可保存的数据")def update_waveform_display(self):"""更新波形显示"""if self.audio_data is not None:self.waveform_ax.clear()for i in range(min(4, self.audio_data.shape[1])):self.waveform_ax.plot(self.audio_data[:1000, i] + i * 0.5, label=f'通道 {i+1}')self.waveform_ax.legend()self.waveform_ax.set_xlabel('采样点')self.waveform_ax.set_ylabel('振幅')self.waveform_ax.set_title('音频波形')self.waveform_canvas.draw()def process_signal(self):"""处理信号"""if self.audio_data is not None:self.status_var.set("处理中...")# 在新线程中处理信号def processing_thread():self.processed_data = self.system.signal_processing.process(self.audio_data,apply_filters=self.apply_filters_var.get(),sync_channels=self.sync_channels_var.get())# 更新显示self.update_processing_display()self.status_var.set("处理完成")thread = threading.Thread(target=processing_thread)thread.daemon = Truethread.start()else:messagebox.showwarning("警告", "请先采集或加载音频数据")def update_processing_display(self):"""更新处理结果显示"""if self.audio_data is not None and self.processed_data is not None:# 更新原始信号显示self.original_ax.clear()self.original_ax.plot(self.audio_data[:1000, 0])self.original_ax.set_title('原始信号 (通道1)')self.original_ax.set_xlabel('采样点')self.original_ax.set_ylabel('振幅')self.original_canvas.draw()# 更新处理后的信号显示self.processed_ax.clear()self.processed_ax.plot(self.processed_data[:1000, 0])self.processed_ax.set_title('处理后的信号 (通道1)')self.processed_ax.set_xlabel('采样点')self.processed_ax.set_ylabel('振幅')self.processed_canvas.draw()# 其他方法类似，限于篇幅不再详细展开def main():"""主函数"""root = tk.Tk()app = IntelligentStethoscopeApp(root)root.mainloop()if __name__ == "__main__":main()

3.2 系统测试与验证

import unittest
import numpy as np
import tempfile
import osclass TestIntelligentStethoscopeSystem(unittest.TestCase):def setUp(self):"""设置测试环境"""self.system = IntelligentStethoscopeSystem()# 生成模拟音频数据self.sample_rate = 44100self.duration = 2  # 2秒self.t = np.linspace(0, self.duration, self.sample_rate * self.duration, endpoint=False)# 生成16通道的测试信号self.test_audio = np.zeros((len(self.t), 16))for i in range(16):freq = 100 + i * 10  # 每个通道频率不同self.test_audio[:, i] = 0.5 * np.sin(2 * np.pi * freq * self.t)# 添加一些噪声noise = 0.1 * np.random.normal(0, 1, len(self.t))self.test_audio[:, i] += noisedef test_data_acquisition(self):"""测试数据采集模块"""# 测试保存和加载功能with tempfile.NamedTemporaryFile(suffix='.npy', delete=False) as tmp_file:tmp_filename = tmp_file.nametry:# 保存数据self.system.data_acquisition.audio_data = self.test_audioself.system.data_acquisition.save_audio(tmp_filename)# 检查文件是否存在self.assertTrue(os.path.exists(tmp_filename))# 加载数据loaded_data = self.system.data_acquisition.load_audio(tmp_filename)# 检查数据是否一致self.assertIsNotNone(loaded_data)self.assertEqual(loaded_data.shape, self.test_audio.shape)np.testing.assert_array_almost_equal(loaded_data, self.test_audio)finally:# 清理临时文件if os.path.exists(tmp_filename):os.unlink(tmp_filename)def test_signal_processing(self):"""测试信号处理模块"""# 测试处理功能processed_data = self.system.signal_processing.process(self.test_audio)# 检查处理后的数据形状self.assertEqual(processed_data.shape, self.test_audio.shape)# 检查数据是否归一化for i in range(processed_data.shape[1]):self.assertLessEqual(np.max(np.abs(processed_data[:, i])), 1.0)def test_feature_extraction(self):"""测试特征提取模块"""# 首先处理信号processed_data = self.system.signal_processing.process(self.test_audio)# 提取特征features = self.system.feature_extraction.extract(processed_data)# 检查是否提取了特征self.assertIsNotNone(features)self.assertGreater(len(features), 0)# 检查特定特征是否存在self.assertIn('mean', features)self.assertIn('std', features)self.assertIn('spectral_centroid', features)def test_ml_analysis(self):"""测试机器学习分析"""# 首先处理信号和提取特征processed_data = self.system.signal_processing.process(self.test_audio)features = self.system.feature_extraction.extract(processed_data)# 训练一个简单的模型X = np.random.rand(100, 10)  # 模拟训练数据y = np.random.randint(0, 2, 100)  # 模拟标签model, accuracy = self.system.ml_models.train_random_forest(X, y)# 检查模型是否训练成功self.assertIsNotNone(model)self.assertGreaterEqual(accuracy, 0)  # 准确率应该在0-1之间self.assertLessEqual(accuracy, 1)def test_report_generation(self):"""测试报告生成"""# 模拟分析结果analysis_results = {'random_forest': np.array([[0.3, 0.7]]),  # 心脏异常概率70%'svm': np.array([[0.8, 0.2]]),  # 呼吸异常概率20%'anomaly_detection': np.array([1])  # 无异常}# 生成报告report = self.system.report_generation.generate(analysis_results)# 检查报告结构self.assertIn('report_id', report)self.assertIn('timestamp', report)self.assertIn('diagnostic_impressions', report)self.assertIn('recommendations', report)# 检查诊断印象impressions = report['diagnostic_impressions']self.assertIn('cardiac_abnormality_probability', impressions)self.assertIn('respiratory_abnormality_probability', impressions)self.assertIn('anomaly_detected', impressions)# 检查建议self.assertGreater(len(report['recommendations']), 0)if __name__ == '__main__':unittest.main()

4. 系统部署与优化

4.1 性能优化策略

import numba
from numba import jit, prange
import multiprocessing as mp
from concurrent.futures import ThreadPoolExecutorclass OptimizedSignalProcessing:def __init__(self):self.sample_rate = 44100@staticmethod@jit(nopython=True, parallel=True)def _apply_filter_numba(b, a, data):"""使用Numba加速滤波器应用"""filtered_data = np.zeros_like(data)n = len(b)for ch in prange(data.shape[1]):# 直接实现FIR滤波器for i in range(n-1, len(data)):filtered_data[i, ch] = 0for j in range(n):filtered_data[i, ch] += b[j] * data[i-j, ch]if j > 0 and i-j >= 0:filtered_data[i, ch] -= a[j] * filtered_data[i-j, ch]filtered_data[i, ch] /= a[0]return filtered_datadef parallel_process(self, audio_data, num_workers=None):"""并行处理多通道数据"""if num_workers is None:num_workers = mp.cpu_count()# 分割数据到不同通道channels = audio_data.shape[1]chunks = []for i in range(channels):chunks.append(audio_data[:, i])# 使用多进程并行处理with ThreadPoolExecutor(max_workers=num_workers) as executor:results = list(executor.map(self._process_channel, chunks))# 合并结果processed_data = np.column_stack(results)return processed_datadef _process_channel(self, channel_data):"""处理单个通道"""# 应用滤波器filtered = self.bandpass_filter_single(channel_data)# 归一化normalized = filtered / np.max(np.abs(filtered))return normalizeddef bandpass_filter_single(self, data, lowcut=20, highcut=1000, order=4):"""单通道带通滤波器"""nyquist = 0.5 * self.sample_ratelow = lowcut / nyquisthigh = highcut / nyquistb, a = signal.butter(order, [low, high], btype='band')filtered_data = signal.filtfilt(b, a, data)return filtered_data

4.2 部署配置

# config.yaml
system:sample_rate: 44100channels: 16chunk_size: 1024processing:apply_filters: truesync_channels: truelowcut: 20highcut: 1500wavelet: db4wavelet_level: 3feature_extraction:include_wavelet: trueinclude_mfcc: truen_mfcc: 13ml_models:random_forest:n_estimators: 100random_state: 42svm:kernel: rbfrandom_state: 42cnn:epochs: 50batch_size: 32anomaly_detection:contamination: 0.1report:default_format: jsoninclude_visualizations: true