【机器学习】乳腺肿瘤超声图像二分类 利用小波变换和自相关系数

超声图像二分类 利用小波变换和自相关系数"/>

【机器学习】乳腺肿瘤超声图像二分类 利用小波变换和自相关系数

利用机器学习分类乳腺超声图像

利用小波变换和自相关系数来获取图像特征

小波变换

对图像做小波变换可以有效的提取出图片中的高频部分，比如边缘、形状等，具体细节可参考其他文章。

这里做完小波变换后的图片，我们将h1,v1,c1叠加求均值，得到图像的轮廓。

左边是原图，右边是经过小波变换叠加后的图，需要聚义的是经过小波变换后图像的大小会为原来的1/2。

自相关系数（纹理信息）

具体内容可以参考link

自相关系数在一定程度上可以反映图像的纹理信息，也可以用灰度共生矩阵来计算纹理特征，这个csdn有就不展开了。

得到的纹理信息可视化是这样的。

将这两个信息融合就能进行简单的分类了。

实验结果

五折交叉验证结果

输出报告

精确度有0.8已结很不错了。

最后要说的是图像不要加预处理，加降噪或者图像均衡都会导致准确率下降，因为预处理可能会导致纹理信息等丢失。

代码

import datetimestarttime = datetime.datetime.now()import numpy as np
import os
import cv2
import sklearn.model_selection as ms
from sklearn.model_selection import train_test_split
from sklearn.metrics import confusion_matrix, classification_report
import matplotlib.pyplot as plt# 数据读取X = []
Y = []
path = os.listdir(r'C:\\Users\\Administrator\\Desktop\\classification\\originals')
classes = {'benign':0, 'malignant':1}
for cls in classes:pth = r'C:\\Users\\Administrator\\Desktop\\classification\\originals\\'+clsfor j in os.listdir(pth):img = cv2.imread(pth+'/'+j, 0)img = cv2.resize(img, (200,200))
#         img = img.reshape(1,-1)/255X.append(img)Y.append(classes[cls])
X = np.array(X)
Y = np.array(Y)# 数据处理
#求自相关系数
def autocorrelation_2D(image):'''Computes the 2D autocorrelation'''data = image-image.mean() data = data / data.std()temp = real(fftshift(ifft2(fft2(data)*conj(fft2(data)))))temp = temp / (data.shape[0]*data.shape[1])return tempimport numpy as np
import pickle
import scipy
from scipy.interpolate import interp1dfrom numpy.fft import fftshift, fft2, ifft2, fft, ifft
from numpy import real
from numpy import conj
aoto = autocorrelation_2D(img)
plt.imshow(aoto, cmap='gray')
T = []
for i in range(X.shape[0]):T.append(autocorrelation_2D(img))T = np.array(T)#单级小波分解，返回分别为低频分量，水平高频，竖直高频，对角高频，分别相当于上图中LL,HL,LH,HH
from pywt import dwt2,wavedec2Z = []
for i in range(X.shape[0]):cA, (cH, cV, cD) = dwt2(X[i], 'db1')Z.append((cH+cV+cD) / 3)Z = np.array(Z)
Z.shape
plt.imshow(Z[2], cmap='gray')#对数据PCA降维
from sklearn.decomposition import PCA
pca = PCA(5)
P = []
# CON.reshape(-1,1)
for i in range(T.shape[0]):pca = pca.fit(T[i])P.append(pca.transform(T[i]))P = np.array(P)
P.shape
plt.imshow(P[2], cmap='gray')from sklearn.decomposition import PCA
pca = PCA(5)
S = []
# CON.reshape(-1,1)
for i in range(Z.shape[0]):pca = pca.fit(Z[i])S.append(pca.transform(Z[i]))S = np.array(S)
S.shape
plt.imshow(S[2], cmap='gray')#将两个数据concat
S = S.reshape(len(S), -1)
P = P.reshape(len(P), -1)
print(P.shape)
print(S.shape)
CON = []
CON = np.hstack((S,P))
CON = np.array(CON)
CON.shape# 交叉验证+roc曲线
from sklearn import svm
from sklearn.metrics import roc_curve, auc
from sklearn.model_selection import StratifiedKFold# 交叉验证
cv = StratifiedKFold(n_splits=5)  # 导入该模型，后面将数据划分6份
classifier = svm.SVC(kernel='rbf', gamma=0.0001, C=1000, probability=True, random_state=0)  # SVC模型 可以换作AdaBoost模型试试# 画平均ROC曲线的两个参数
mean_tpr1 = 0.0  # 用来记录画平均ROC曲线的信息
mean_fpr1 = np.linspace(0, 1, 100)
cnt = 0
for i, (train, test) in enumerate(cv.split(CON, Y)):  # 利用模型划分数据集和目标变量 为一一对应的下标cnt += 1probas_ = classifier.fit(CON[train], Y[train]).predict_proba(CON[test])  # 训练模型后预测每条样本得到两种结果的概率fpr, tpr, thresholds = roc_curve(Y[test], probas_[:, 1])  # 该函数得到伪正例、真正例、阈值，这里只使用前两个mean_tpr1 += np.interp(mean_fpr1, fpr, tpr)  # 插值函数 interp(x坐标,每次x增加距离,y坐标)  累计每次循环的总值后面求平均值mean_tpr1[0] = 0.0  # 将第一个真正例=0 以0为起点roc_auc = auc(fpr, tpr)  # 求auc面积plt.plot(fpr, tpr, lw=1, label='ROC fold {0:.2f} (area = {1:.2f})'.format(i, roc_auc))  # 画出当前分割数据的ROC曲线print(classifier.score(CON[test], Y[test])) pre_y = classifier.predict(CON[test])pre_y = np.array(pre_y)
#     confusion_matrix = confusion_matrix(Y[test],pre_y)
#     print(confusion_matrix)print(classification_report(Y[test], pre_y))plt.plot([0, 1], [0, 1], '--', color=(0.6, 0.6, 0.6), label='Luck')  # 画对角线mean_tpr1 /= cnt  # 求数组的平均值
mean_tpr1[-1] = 1.0  # 坐标最后一个点为（1,1）  以1为终点
mean_auc1 = auc(mean_fpr1, mean_tpr1)plt.plot(mean_fpr1, mean_tpr1, 'k--', label='Mean ROC (area = {0:.2f})'.format(mean_auc1), lw=2)plt.xlim([-0.05, 1.05])  # 设置x、y轴的上下限，设置宽一点，以免和边缘重合，可以更好的观察图像的整体
plt.ylim([-0.05, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')  # 可以使用中文，但需要导入一些库即字体
plt.title('Receiver operating characteristic example')
plt.legend(loc="lower right")
plt.show()