使用SVM对分泌效应蛋白数据集进行分类预测

1、SVM简介 支持向量机(Support Vector Machine, SVM)是一类按监督学习(supervised learning)方式对数据进行二元

1、SVM简介
支持向量机(Support Vector Machine, SVM)是一类按监督学习(supervised learning)方式对数据进行二元分类的广义线性分类器(generalized linear classifier),其决策边界是对学习样本求解的最大边距超平面(maximum-margin hyperplane 。
在机器学习中,支持向量机是在分类与回归分析中分析数据的监督式学习模型与相关的学习算法。给定一组训练实例,每个训练实例被标记为属于两个类别中的一个或另一个,SVM训练算法创建一个将新的实例分配给两个类别之一的模型,使其成为非概率二元线性分类器。

2、Sklearn文档中代码实现
和其他分类器一样, SVC, NuSVC 和 LinearSVC 将两个数组作为输入: [n_samples, n_features] 大小的数组 X 作为训练样本, [n_samples] 大小的数组 y 作为类别标签(字符串或者整数):

from sklearn import svm
>>> X = [[0, 0], [1, 1]]
>>> y = [0, 1]
>>> clf = svm.SVC(gamma='scale')
>>> clf.fit(X, y)  
SVC(C=1.0, cache_size=200, class_weight=None, coef0=0.0,decision_function_shape='ovr', degree=3, gamma='scale', kernel='rbf',max_iter=-1, probability=False, random_state=None, shrinking=True,tol=0.001, verbose=False)

在拟合后, 这个模型可以用来预测新的值:

 clf.predict([[2., 2.]])
array([1])

SVMs 决策函数取决于训练集的一些子集, 称作支持向量. 这些支持向量的部分特性可以在 support_vectors_, support_ 和 n_support 找到:

>>> # 获得支持向量
>>> clf.support_vectors_
array([[ 0.,  0.],[ 1.,  1.]])
>>> # 获得支持向量的索引
>>> clf.support_
array([0, 1]...)
>>> # 为每一个类别获得支持向量的数量
>>> clf.n_support_
array([1, 1]...)

3、自己加工了一下的代码(对丁香型假单胞菌III型分泌效应蛋白数据集进行分类预测)

from __future__ import print_function
from sklearn.metrics  import roc_curve,auc
import matplotlib.pyplot  as plt
from sklearn.metrics import precision_score, recall_score, f1_score,matthews_corrcoef
import numpy as np
from keras.preprocessing import sequence
import pickle as cPickle
import datetime
import random
from sklearn.svm import SVC
from sklearn import metrics
# 时间差计算函数
def subtime(date1, date2):date1 = datetime.datetime.strptime(date1, "%Y-%m-%d %H:%M:%S")date2 = datetime.datetime.strptime(date2, "%Y-%m-%d %H:%M:%S")return date2 - date1def trans(str1):a = []dic = {'A':1,'B':22,'U':23,'J':24,'Z':25,'O':26,'C':2,'D':3,'E':4,'F':5,'G':6,'H':7,'I':8,'K':9,'L':10,'M':11,'N':12,'P':13,'Q':14,'R':15,'S':16,'T':17,'V':18,'W':19,'Y':20,'X':21}for i in range(len(str1)):a.append(dic.get(str1[i]))return adef createTrainData(str1):sequence_num = []label_num = []for line in open(str1,encoding='utf-8'):#print(line)proteinId, sequence,label= line.split(",")proteinId = proteinId.strip(' \t\r\n').split();sequence = sequence.strip(' \t\r\n');# print(sequence.isspace())sequence_num.append(trans(sequence))label = label.strip(' \t\r\n').split();label_num.append(label)return sequence_num,label_numa,b=createTrainData("D:\lunwen\deep1\data//T3SE//syringae.csv")
t = (a, b)
cPickle.dump(t,open("D:\lunwen\deep1\data//T3SE//syringae.pkl","wb"))def createTrainTestData(str_path, nb_words=None, skip_top=0,maxlen=None, test_split=0.20, seed=800,start_char=1, oov_char=2, index_from=3):X,labels = cPickle.load(open(str_path, "rb"))np.random.seed(seed)np.random.shuffle(X)np.random.seed(seed)np.random.shuffle(labels)if start_char is not None:X = [[start_char] + [w + index_from for w in x] for x in X]elif index_from:X = [[w + index_from for w in x] for x in X]if maxlen:new_X = []new_labels = []for x, y in zip(X, labels):if len(x) < maxlen:new_X.append(x)new_labels.append(y)X = new_Xlabels = new_labelsif not X:raise Exception('After filtering for sequences shorter than maxlen=' +str(maxlen) + ', no sequence was kept. ''Increase maxlen.')if not nb_words:nb_words = max([max(x) for x in X])if oov_char is not None:X = [[oov_char if (w >= nb_words or w < skip_top) else w for w in x] for x in X]else:nX = []for x in X:nx = []for w in x:if (w >= nb_words or w < skip_top):nx.append(w)nX.append(nx)X = nXX_train = np.array(X[:int(len(X) * (1 - test_split))])y_train = np.array(labels[:int(len(X) * (1 - test_split))])X_test = np.array(X[int(len(X) * (1 - test_split)):])y_test = np.array(labels[int(len(X) * (1 - test_split)):])return (X_train, y_train), (X_test, y_test)
print('Loading data...')
startdate = datetime.datetime.now() # 获取当前时间
startdate = startdate.strftime("%Y-%m-%d %H:%M:%S") # 当前时间转换为指定字符串格式(X_train, y_train), (X_test, y_test) = createTrainTestData("D:\lunwen\deep1\data//T3SE//syringae.pkl",nb_words=23, test_split=0.2)
print(len(X_train), 'train sequences')
print(len(X_test), 'test sequences')print('Pad sequences (samples x time)')
X_train = sequence.pad_sequences(X_train, maxlen=1000)
X_test = sequence.pad_sequences(X_test, maxlen=1000)
print('X_train shape:', X_train.shape)
print('X_test shape:', X_test.shape)
print('Build model...')svc = SVC().fit(X_train,y_train)
y_prob = svc.decision_function(X_test)                              # 决策边界距离
y_pred = svc.predict(X_test)           # 模型对测试集的预测结果fpr_svc,tpr_svc,threshold_svc = metrics.roc_curve(y_test.astype(int),y_prob)     # 获取真阳率、伪阳率、阈值
auc_svc = metrics.auc(fpr_svc,tpr_svc)                              # 模型准确率
score_svc = metrics.accuracy_score(y_test,y_pred)
tn, fp, fn, tp = metrics.confusion_matrix(y_test,y_pred).ravel()
print([score_svc,auc_svc])print('***********************************************************************')
print('SP: %.4f' % (tn/(tn+fp)))
print('SN: %.4f' % (tp/(tp+fn)))
# 正确率  (提取出的正确信息条数 / 提取出的信息条数)
print('Precision: %.4f' % precision_score(y_true=y_test, y_pred=y_pred, pos_label='1'))
# 召回率 (提出出的正确信息条数 / 样本中的信息条数)
print('Recall: %.4f' % recall_score(y_true=y_test, y_pred=y_pred, pos_label='1'))
# 马修相关系数
print('MCC: %.4f' % matthews_corrcoef(y_true=y_test, y_pred=y_pred))
# F1-score  (正确率*召回率*2 /(正确率+召回率))
print('F1: %.4f' % f1_score(y_true=y_test, y_pred=y_pred,pos_label='1'))
print('***********************************************************************')enddate = datetime.datetime.now() # 获取当前时间
enddate = enddate.strftime("%Y-%m-%d %H:%M:%S") # 当前时间转换为指定字符串格式
# 计算训练时长
print('start date ',startdate)
print('end date ',enddate)
print('Time ',subtime(startdate,enddate)) # enddate > startdate
roc_auc = auc(fpr_svc, tpr_svc)
lw = 2
plt.figure(figsize=(8, 5))
plt.plot(fpr_svc, tpr_svc, color='darkorange',lw=lw, label='ROC curve (area = %0.4f)' % roc_auc)  ###假正率为横坐标,真正率为纵坐标做曲线
plt.plot([0, 1], [0, 1], color='navy', lw=lw, linestyle='--')
plt.xlim([0.0, 1.0])
plt.ylim([0.0, 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()

最终效果图片:
在这里插入图片描述
在这里插入图片描述

自己做的一点点记录。不知道是否正确,仅作参考。