机器学习简单小总结 解决糖尿病分类问题和sin函数回归问题

记录

数据分析

分类任务

数据信息

数据条数	标签为1	标签为0	数据维度
576	197	379	8

数据分析

由图像可分析,糖尿病患者其属性大部分都略高于正常人，1、2、8属性尤为明显
3、4,、5属性与正常人差异几乎没有。

回归任务

5000个数据点,由图像可知显然是一个正弦函数回归数据

BP分类任务

• 结构 使用三层Linear神经网络,使用$\sigma (x)$激活函数，均方差损失函数
  $$\sigma (x)=\frac{1}{1+e^{-x}}{E}_{loss}=\frac{1}{n}\sum _{i=1}^n(\hat{y}-y{)}^2$$

• 参数 每层的参数为(8,8),(8,5),(5,1)，学习率0.02，迭代次数10000次
  
  

• 测试结果
  
  
  训练loss变化情况如下：
  

大致看出在训练20个左右单位圈数时得到最佳效果，在之后发生过拟合现象。(一个单位是100epoch)


迭代10000次后,在576个样本的训练集中，482 个数据分类正确，准确率0.837，在192个样本的测试集中144个数据分类正确，准确率0.75,显然发生过拟合。


在训练过程中产生的最佳效果（未发生过拟合时）192个测试集中，156个数据分类正确，准确率0.8125

结果

混淆矩阵（尚未过拟合时的最好效果）

真实情况	预测结果
	正例	反例
正例	36(TP)	30(FN)
反例	6（FP)	120(TN)

查准率 $P=0.857$
查全率 $R=0.545$

AUC 为 0.8541366041366042

SVM分类任务

• 结构 使用核函数和软间隔加持的支持向量机
  
  

• 参数 使用线性核函数 ，$C=3000$的软间隔惩罚，10次
  
  

• 测试结构在576个训练集中，432个数据分类成功，准确率0.75，192个测试集中，150个数据分类正确，准确率0.78125

混淆矩阵

真实情况	预测结果
	正例	反例
正例	49(TP)	17(FN)
反例	25（FP)	101(TN)

查准率 $P=0.662$
查全率 $R=0.742$

AUC is : 0.8733766233766234

beyesian分类任务

• 结构 使用高斯分布的朴素贝叶斯分类器
  
  

• 参数 没有任何参数，直接统计也不需要迭代
  
  

• 测试结构在576个训练集中，434个数据分类成功，准确率0.7535，在192个测试集中，151个数据分类正确，准确率0.7865

混淆矩阵

真实情况	预测结果
	正例	反例
正例	41(TP)	25(FN)
反例	16（FP)	110(TN)

查准率 $P=0.719$
查全率 $R=0.621$

AUC is : 0.8464405964405965

BP回归任务

• 结构 使用两层Linear神经网络,使用$\sigma (x)$激活函数，均方差损失函数
  $$\sigma (x)=\frac{1}{1+e^{-x}}{E}_{loss}=\frac{1}{n}\sum _{i=1}^n(\hat{y}-y{)}^2$$

• 参数 每层的参数为(1,15),(15,1)，学习率0.0015，迭代次数14000次
  
  

• 测试结果
  
  
  训练loss变化图
  

测试集对比结果

线性回归

分析数据可知，该问题很难直接使用线性回归问题，也无法使用广义线性回归求解。
$$y=wx+b$$
先对数据中的x先升维处理(暂且不易引入核函数升维)。

数据本身为1维度，暂把数据升维4维数据,升维方式如：
$$(x{)}_i\to (x^2,x^4,x^6,x^8{)}_i$$
再安装多元线性回归求解
$$W=(X^{T}X{)}^{-1}{X}^{T}Y$$

小结

	支持向量机分类器	朴素贝叶斯分类器	BP神经网络分类器
训练集	0.750	0.754	0.781
测试集	0.782	0.787	0.813
AUC	0.873	0.846	0.854

相关代码

读入数据及其分析

import numpy as np


def getClassificationData(str):
    f = open(str, 'r')
    lines = f.readlines()
    f.close()
    dataSet = []
    Label = []
    for line in lines:
        s = line.split()
        data = [float(ch) for ch in s]
        Label.append([int(data[0])])
        data = data[1:]
        dataSet.append(data)
    return map(np.array, (dataSet, Label))


def getRegressionData(str):
    f = open(str, 'r')
    lines = f.readlines()
    f.close()
    dataSet = []
    Label = []
    for line in lines:
        s = line.split()
        dataSet.append([float(s[1])])
        Label.append([float(s[0])])
    return map(np.array, (dataSet, Label))


def viewRegressionData():
    data, Label = getRegressionData('regression/sinc_train')
    plt.plot(data.reshape(-1), Label.reshape(-1), 'o')
    plt.show()


def analyseClassificationData():
    data, Label = getClassificationData('classification/diabetes_train')
    data_a = []
    data_b = []
    for da, la in zip(data, Label):
        if la == 1:
            data_a.append(da)
        else:
            data_b.append(da)
    print("1:", len(data_a), "0:", len(data_b))
    data_a = np.array(data_a)
    data_b = np.array(data_b)
    for i in range(8):
        plt.subplot(1, 2, 1)
        plt.hist([data_a[:, i], data_b[:, i]], bins=20, range=[0, 1], density=True, cumulative=False)
        plt.title("histogram of attribute" + str(i+1))
        plt.subplot(1, 2, 2)
        plt.hist([data_a[:, i], data_b[:, i]], bins=20, range=[0, 1], density=True, cumulative=True)
        plt.title("cumulative of attribute" + str(i+1))
        plt.savefig("fig"+str(i+1))
        plt.show()

朴素贝叶斯分类器

import numpy as np
import ReadData

def Norm_distribute(x, sigma, mean):
    return 1 / ((2 * np.pi) ** 0.5 * sigma) * np.exp(-(x - mean) ** 2 / (2 * sigma ** 2))


class Beyesian_continue:
    def __init__(self, dataSet, label):
        dic = {}
        for la in label:
            if dic.__contains__(la) is False:
                dic[la] = 0
            dic[la] = dic[la] + 1
        lable_kind = len(dic.keys())
        dim = len(dataSet[0])
        totNum = len(dataSet)

        self.dim = dim
        self.C_num = lable_kind
        self.PC = np.zeros(lable_kind)
        self.EX = np.zeros([self.C_num, self.dim])
        self.sig = np.zeros([self.C_num, self.dim])
        for c in range(self.C_num):
            if dic.__contains__(c):
                self.PC[c] = dic[c] / totNum

        for u in range(totNum):
            for d in range(dim):
                self.EX[label[u], d] = self.EX[label[u], d] + dataSet[u, d]
                self.sig[label[u], d] = self.sig[label[u], d] + dataSet[u, d] ** 2

        for c in range(self.C_num):
            for d in range(dim):
                self.EX[c, d] = self.EX[c, d] / dic[c]
                self.sig[c, d] = (self.sig[c, d] / dic[c] - self.EX[c, d] ** 2) ** 0.5

    def predict(self, x):
        t = []
        for c in range(self.C_num):
            p = self.PC[c]
            for d in range(self.dim):
                p = p * Norm_distribute(x[d], self.sig[c, d], self.EX[c, d])
            t.append(p)
        return np.argmax(t)

    def __call__(self, *args, **kwargs):
        return self.predict(*args)


train_x, train_y = ReadData.getClassificationData('classification/diabetes_train')
test_x, test_y = ReadData.getClassificationData('classification/diabetes_test')
beyes = Beyesian_continue(train_x, train_y.reshape(-1))

ans = 0
for (tx, ty) in zip(train_x, train_y):
    ans = ans + (beyes(tx) == ty)[0]
print(ans, ans / len(train_x))

ans = 0
TP = 0
TN = 0
FN = 0
FP = 0
for (tx, ty) in zip(test_x, test_y):
    haty = beyes(tx)
    ans = ans + (haty == ty)[0]
    TP = TP + (haty == 1) * (ty == 1)[0]
    TN = TN + (haty == 0) * (ty == 0)[0]
    FN = FN + (haty == 0) * (ty == 1)[0]
    FP = FP + (haty == 1) * (ty == 0)[0]
print(ans, ans / len(test_x))

print('TP\tFN\t[', TP, ' ', FN, ']')
print('FP\tTN\t[', FP, ' ', TN, ']')
print('P =', TP / (TP + FP), ' R=', TP / (TP + FN))
Ttot = np.sum(test_y == 1)
Ftot = np.sum(test_y == 0)

point = np.array([[beyes.getP(test_x[i], 1), test_y[i][0]] for i in range(len(test_x))])

l = len(point)
for i in range(l):
    j = np.argmax(point[i:, 0]) + i
    t = point[i].copy()
    point[i] = point[j]
    point[j] = t

cnt_T = 0
TPR = []
FPR = []
for i in range(l):
    cnt_T = cnt_T + point[i][1]
    TPR.append(cnt_T / Ttot)
    FPR.append((i - cnt_T + 1) / Ftot)
TPR = np.array(TPR)
FPR = np.array(FPR)

plt.plot(FPR, TPR)
font = FontProperties(fname=r"C:\Windows\Fonts\simhei.ttf", size=14)
plt.title('贝叶斯分类器 ROC曲线', fontproperties=font)
plt.xlabel("假正例率", fontproperties=font)
plt.ylabel("真正例率", fontproperties=font)
print('AUC is :', np.sum((TPR[1:] + TPR[:-1]) * 0.5 * (FPR[1:] - FPR[:-1])))
plt.show()

支持向量机分类器

import numpy as np
import matplotlib.pyplot as plt
from ML import ReadData
from matplotlib.font_manager import FontProperties

class SVM:
    def __init__(self, dim, kennel=None, Name='svm'):
        self.dim = dim
        self.w = np.zeros([1, dim])
        self.b = np.zeros([1, 1])
        self.kennel = kennel
        self.LamAndY = np.zeros([1, 1])
        self.data = np.zeros([1, 1])
        self.Name = Name

    def train(self, x, y, epochs=20, callback=None, Lam=None, C_val=np.inf, Log=True):
        N = len(x)
        Mat_dot = np.zeros([N, N])
        Mat_dot_y = np.zeros([N, N])
        for i in range(N):
            for j in range(N):
                if self.kennel is None:
                    Mat_dot[i, j] = np.dot(x[i], x[j])
                else:
                    Mat_dot[i, j] = self.kennel(x[i], x[j])
                Mat_dot_y[i, j] = Mat_dot[i, j] * y[i] * y[j]
        if Lam is None:
            Lam = np.random.rand(N)
            delta = np.dot(Lam, y)
            for i in range(N):
                if delta * y[i] < 0:
                    Lam[i] = Lam[i] - y[i] * delta
                    delta = 0
                    if Lam[i] > C_val:
                        delta = (Lam[i] - C_val) * y[i]
                        Lam[i] = C_val
        if Log:
            print('初始化完毕')
        for epoch in range(epochs):
            for a in range(N):
                for b in range(N):
                    subb = Mat_dot[a, a] - Mat_dot[b, a] * 2 + Mat_dot[b, b]
                    if subb == 0:
                        continue

                    C = - Lam[a] * y[a] - Lam[b] * y[b]
                    Sum_Lam_a = np.dot(Lam, Mat_dot_y[:, a]) - Lam[a] * Mat_dot_y[a, a] - Lam[b] * Mat_dot_y[b, a]
                    Sum_Lam_b = np.dot(Lam, Mat_dot_y[:, b]) - Lam[a] * Mat_dot_y[a, b] - Lam[b] * Mat_dot_y[b, b]

                    suba = y[a] * C * (Mat_dot[a, b] - Mat_dot[b, b]) + 1 - y[a] * y[b] - Sum_Lam_a + y[a] * \
                           y[b] * Sum_Lam_b

                    U = suba / subb
                    V = -(C + U * y[a]) * y[b]
                    rate = 1
                    if U < 0:
                        rate = min(rate, (0 - Lam[a]) / (U - Lam[a]))
                    if V < 0:
                        rate = min(rate, (0 - Lam[b]) / (V - Lam[b]))
                    if U > C_val:
                        rate = min(rate, (C_val - Lam[a]) / (U - Lam[a]))
                    if V > C_val:
                        rate = min(rate, (C_val - Lam[b]) / (V - Lam[b]))
                    Lam[a] = Lam[a] + (U - Lam[a]) * rate
                    Lam[b] = Lam[b] + (V - Lam[b]) * rate
            if callback is not None:
                callback(x, y, Lam, self.dim, self.kennel)
        if self.kennel is None:
            w = np.zeros([1, self.dim])
            for i in range(N):
                w = w + Lam[i] * y[i] * x[i]
            h = []
            for i in range(N):
                if Lam[i] > 1e-9:
                    h.append(np.dot(w, x[i]))
            h = np.array(h)
            self.w = w
            self.b[0][0] = (np.max(h) + np.min(h)) / 2
        else:
            self.LamAndY = []
            for i in range(N):
                self.LamAndY.append(Lam[i] * y[i])
            self.data = x
            h = []
            for i in range(N):
                if Lam[i] < 1e-9:
                    continue
                res = 0
                for lay, vec in zip(self.LamAndY, self.data):
                    res = res - lay * self.kennel(vec, self.data[i])
                h.append(res)
            h = np.array(h)
            self.b = (np.max(h) + np.min(h)) / 2
        if Log:
            print('训练结束')
        return Lam

    def forward(self, x):
        if self.kennel is None:
            return self.w @ x.T + self.b
        else:
            ret = np.zeros([1, len(x)])
            for lay, vec in zip(self.LamAndY, self.data):
                ret = ret + lay * self.kennel(vec, x)
            return ret + self.b

    def save(self):
        np.save('save/' + self.Name + '_LamAndY', np.array(self.LamAndY))
        np.save('save/' + self.Name + '_w', np.array(self.w))
        np.save('save/' + self.Name + '_b', np.array(self.b))

    def load(self, data=None):
        try:
            self.LamAndY = np.load('save/' + self.Name + '_LamAndY.npy')
            self.w = np.load('save/' + self.Name + '_w.npy')
            self.b = np.load('save/' + self.Name + '_b.npy')
        finally:
            self.data = data

    def __call__(self, *args, **kwargs):
        return self.forward(*args)


class TestSVM:
    def __init__(self, N=50):
        x = np.random.rand(N, 2)
        x[:, 1] = x[:, 1] * 2
        y = ((x[:, 0] * 3 - 1) > x[:, 1]) * 2 - 1
        for i in range(len(x)):
            if y[i] < 0:
                x[i][1] = x[i][1] - 0.5
        # y = (np.random.rand(N) > 0.5) * 2 - 1
        self.x = x
        self.y = y
        self.point_x = []
        self.point_y = []
        for a, b in zip(x, y):
            if b == 1:
                self.point_x.append(a.tolist())
            else:
                self.point_y.append(a.tolist())
        self.point_x, self.point_y = map(np.array, (self.point_x, self.point_y))
        self.svm = SVM(2, kennel=self.kennel, Name='small_svm')

    def kennel(self, x, y):
        if len(x.shape) == 1 and len(y.shape) == 1:
            return np.dot(x, y)
        return np.dot(x.reshape([1, -1]), y.T)

    def test(self, train=True):
        plt.ion()
        if train:
            self.svm.train(self.x, self.y, 30, self.callback, C_val=1000)
        conlose = (self.svm(self.x) > 0) * 2 - 1
        print(conlose, self.y)
        print(np.sum(conlose == self.y), '/', len(self.x))
        plt.ioff()

    def callback(self, x, y, Lam, dim, kennel):
        w = np.zeros([1, dim])
        for i in range(len(x)):
            w = w + Lam[i] * y[i] * x[i]

        h = []
        for i in range(len(x)):
            if Lam[i] < 1e-9:
                continue
            res = 0

            for la_, y_, vec in zip(Lam, y, x):
                res = res - la_ * y_ * kennel(vec, x[i])
            h.append(res)
        h = np.array(h)
        b = (np.max(h) + np.min(h)) / 2
        lx = np.arange(0, 1, 0.01)
        ly = -(w[0][0] * lx + b) / w[0][1]
        plt.cla()
        plt.plot(self.point_x[:, 0], self.point_x[:, 1], 'o')
        plt.plot(self.point_y[:, 0], self.point_y[:, 1], 'o')
        plt.plot(lx, ly)
        plt.pause(0.05)
        cal = np.sum(Lam)

        for i in range(len(x)):
            for j in range(len(x)):
                cal = cal - 0.5 * Lam[i] * Lam[j] * y[i] * y[j] * np.dot(x[i], x[j])
        print(cal, w, b)


class TestSVM_Classification:
    def __init__(self, TestNum=1000, ):
        self.train_x, self.train_y = ReadData.getClassificationData('classification/diabetes_train')
        self.test_x, self.test_y = ReadData.getClassificationData('classification/diabetes_test')
        self.train_y = self.train_y * 2 - 1
        self.test_y = self.test_y * 2 - 1
        if TestNum < len(self.test_x):
            self.train_x = self.train_x[:TestNum, ]
            self.train_y = self.train_y[:TestNum, ]
        self.svm = SVM(len(self.train_x[0]),kennel=self.linearKennel, Name='svm_Classification')
        self.gamma = 1
        self.Lam = np.zeros([1, len(self.train_x)])

    def train(self, epochs=20, C=1e9):  # 训练模型
        self.Lam = self.svm.train(self.train_x, self.train_y, epochs=epochs, callback=self.callback, C_val=C)
        print(np.sum((self.svm(self.train_x) > 0) * 2 - 1 == self.train_y.T), len(self.train_x))

    def trainAgain(self, epochs=20, C=1e9, log=False):  # 在原理基础上继续训练
        self.svm.train(self.train_x, self.train_y, epochs=epochs, callback=self.callback, Lam=self.Lam, C_val=C,
                       Log=log)
        right = np.sum((self.svm(self.train_x) > 0) * 2 - 1 == self.train_y.T)
        print(right, right / len(self.train_x))

    def test(self):  # 用测试集测试模型
        self.svm.data = self.train_x
        hhy = self.svm(self.test_x)
        haty = (hhy > 0) * 2 - 1
        print(haty)
        print(self.test_y.T)
        TP = np.sum((haty == 1) * (self.test_y.T == 1))
        TN = np.sum((haty == -1) * (self.test_y.T == -1))
        FN = np.sum((haty == -1) * (self.test_y.T == 1))
        FP = np.sum((haty == 1) * (self.test_y.T == -1))
        print((TP + TN), (TP + TN) / len(self.test_x))
        print('TP\tFN\t[', TP, ' ', FN, ']')
        print('FP\tTN\t[', FP, ' ', TN, ']')
        print('P =', TP / (TP + FP), ' R=', TP / (TP + FN))
        print("F1=",2*TP/(TP-TN+len(self.svm.data)))
        Ttot = np.sum(self.test_y == 1)
        Ftot = np.sum(self.test_y == -1)
        point = np.array([[p[0], y[0]] for p, y in zip(hhy.T, self.test_y)])
        l = len(point)
        for i in range(l):
            j = np.argmax(point[i:, 0]) + i
            t = point[i].copy()
            point[i] = point[j]
            point[j] = t
        cnt_T = 0
        TPR = []
        FPR = []
        for i in range(l):
            cnt_T = cnt_T + (point[i][1] == 1) * 1
            TPR.append(cnt_T / Ttot)
            FPR.append((i - cnt_T + 1) / Ftot)
        TPR = np.array(TPR)
        FPR = np.array(FPR)
        plt.plot(FPR, TPR)
        font = FontProperties(fname=r"C:\Windows\Fonts\simhei.ttf", size=14)
        plt.title('SVM的ROC曲线', fontproperties=font)
        plt.xlabel("假正例率", fontproperties=font)
        plt.ylabel("真正例率", fontproperties=font)
        print('AUC is :', np.sum((TPR[1:] + TPR[:-1]) * 0.5 * (FPR[1:] - FPR[:-1])))
        plt.show()
        return FPR, TPR

    def save(self):  # 保存模型
        np.save('save/' + self.svm.Name + '_Lam.npy', self.Lam)
        self.svm.save()

    def load(self):  # 重载模型
        self.svm.load()
        self.Lam = np.load('save/' + self.svm.Name + '_Lam.npy')

    def kennel(self, x, y):  # 高斯核函数
        if len(x.shape) == 1 and len(y.shape) == 1:
            return np.exp(-np.dot(x - y, x - y) / (2 * self.gamma ** 2))
        return np.exp(- np.sum((x - y) ** 2, axis=1) / (2 * self.gamma ** 2))

    def linearKennel(self, x, y):
        if len(x.shape) == 1 and len(y.shape) == 1:
            return np.dot(x, y)
        return np.dot(x.reshape([1, -1]), y.T)

    def callback(self, x, y, Lam, dim, kennel):  # 回调函数
        cal = np.sum(Lam)
        for i in range(len(x)):
            for j in range(len(x)):
                cal = cal - 0.5 * Lam[i] * Lam[j] * y[i] * y[j] * np.dot(x[i], x[j])
        print(cal)

#测试的小例子
#test = TestSVM(30)
#test.test()

csvm = TestSVM_Classification()
csvm.train(5, 2100)
csvm.test()

BP神经网络分类器

import matplotlib.pyplot as plt
import numpy as np
import ReadData
from matplotlib.font_manager import FontProperties


class Sigmoid:
    def __call__(self, x):
        return 1 / (1 + np.exp(-x))

    def grad(self, d):
        return d * (1 - d)


class Linear:
    def __init__(self, in_size, out_size, Func=Sigmoid(), lr=0.02, Name='BpNet'):
        self.Name = Name
        self.weight = np.random.rand(in_size, out_size)  # shape=(in_size,out_size)
        self.bias = np.random.rand(1, out_size)  # shape = (1,out_size)
        self.data_in = []
        self.data_out = []
        self.Func = Func
        self.lr = lr

    def save(self):
        np.save('save/' + self.Name + '_weight', np.array(self.weight))
        np.save('save/' + self.Name + '_bias', np.array(self.bias))

    def load(self):
        self.weight = np.load('save/' + self.Name + '_weight.npy')
        self.bias = np.load('save/' + self.Name + '_bias.npy')

    def forward(self, x):
        self.data_in = np.reshape(x, [-1, 1, self.weight.shape[0]])  # shape = (-1,1,in_size)
        self.data_out = self.Func(self.data_in @ self.weight + self.bias)  # shape =(-1,1,out_size)
        return self.data_out

    def backward(self, grad):
        grad = grad.reshape(self.data_out.shape) * self.Func.grad(self.data_out)  # shape =(-1,1,out_size)
        xgrad = grad @ self.weight.T  # shape =(-1,1,in_size)
        self.bias = self.bias - np.sum(grad, axis=0) * self.lr
        self.weight = self.weight - np.sum(self.data_in.swapaxes(1, 2) @ grad, axis=0) * self.lr
        return xgrad

    def __call__(self, *args, **kwargs):
        return self.forward(*args)


class ClasssificationBPNet:
    def __init__(self, Name='BP_Classification', lr=0.02):
        self.l1 = Linear(8, 8, Name=Name + '_layer1', lr=lr)
        self.l2 = Linear(8, 5, Name=Name + '_layer2', lr=lr)
        self.l3 = Linear(5, 1, Name=Name + '_layer3', lr=lr)

    def forward(self, x):
        x = self.l3(self.l2(self.l1(x)))
        return x.reshape((-1, x.shape[2]))

    def backward(self, grad):
        x = self.l1.backward(self.l2.backward(self.l3.backward(grad)))

    def save(self):
        self.l1.save()
        self.l2.save()
        self.l3.save()

    def load(self):
        self.l1.load()
        self.l2.load()
        self.l3.load()

    def __call__(self, *args, **kwargs):
        return self.forward(*args)


def Classification():
    train_x, train_y = ReadData.getClassificationData('classification/diabetes_train')
    test_x, test_y = ReadData.getClassificationData('classification/diabetes_test')
    Bp = ClasssificationBPNet()
    epochs = 10000
    a = []
    b = []
    maxval_test = 0
    maxval_train = 0
    for epoch in range(epochs):
        dy = Bp(train_x) - train_y
        Bp.backward(2 * dy)
        if (epoch + 1) % 100 == 0:
            loss = np.mean(dy * dy)
            test_dy = Bp(test_x) - test_y
            testloss = np.mean(test_dy * test_dy)
            a.append(loss)
            b.append(testloss)
            print('epoch=%5d loss %.8f testloss %.8f' % (epoch + 1, loss, testloss), end='\r')
        te = np.sum((Bp(test_x) > 0.5) == test_y)
        tr = np.sum((Bp(train_x) > 0.5) == train_y)
        if maxval_test < te:
            maxval_test = te
            maxval_train = tr
            Bp.save()
        print(tr, '\t', tr / len(train_x), te, '\t', te / len(test_x), '\t', end='\r')
    print('best recode')
    print(maxval_train, '\t', maxval_train / len(train_x))
    print(maxval_test, '\t', maxval_test / len(test_x))
    Bp.load()
    print('*' * 20)
    Bp.load()
    haty = (Bp(test_x) > 0.5) * 1
    TP = np.sum((haty == 1) * (test_y == 1))
    TN = np.sum((haty == 0) * (test_y == 0))
    FN = np.sum((haty == 0) * (test_y == 1))
    FP = np.sum((haty == 1) * (test_y == 0))
    print('TP\tFN\t[', TP, ' ', FN, ']')
    print('FP\tTN\t[', FP, ' ', TN, ']')
    print('P =', TP / (TP + TP), ' R=', TP / (TP + FN))
    print('*' * 20)

    plt.figure(1)
    plt.plot(a)
    plt.plot(b)
    plt.show()


def ROC_AUC_Classification():
    test_x, test_y = ReadData.getClassificationData('classification/diabetes_test')
    Bp = ClasssificationBPNet()
    Bp.load()
    predict = Bp(test_x)
    Ttot = np.sum(test_y == 1)
    Ftot = np.sum(test_y == 0)
    point = np.array([[p[0], y[0]] for p, y in zip(predict, test_y)])
    l = len(point)
    for i in range(l):
        j = np.argmax(point[i:, 0]) + i
        t = point[i].copy()
        point[i] = point[j]
        point[j] = t
    cnt_T = 0
    TPR = []
    FPR = []
    for i in range(l):
        cnt_T = cnt_T + point[i][1]
        TPR.append(cnt_T / Ttot)
        FPR.append((i - cnt_T + 1) / Ftot)
    TPR = np.array(TPR)
    FPR = np.array(FPR)
    plt.figure(2)
    plt.plot(FPR, TPR)
    font = FontProperties(fname=r"C:\Windows\Fonts\simhei.ttf", size=14)
    plt.title('Bp神经网络ROC曲线', fontproperties=font)
    plt.xlabel("假正例率", fontproperties=font)
    plt.ylabel("真正例率", fontproperties=font)
    print('AUC is :', np.sum((TPR[1:] + TPR[:-1]) * 0.5 * (FPR[1:] - FPR[:-1])))
    plt.show()


Classification()
ROC_AUC_Classification()

支持向量机cpp版

python跑的太慢了，cpp跑的还快一些，方便调参

#include<bits/stdc++.h>
using namespace std;
class SVM
{
public:
    int dim,datalen;
    double b,**data,*lable,*Lamda,**Mat_dot,**Mat_doty,C_val;
    double (*kennel)(double*,double*,int);
    SVM(int _dim,double (*_kennel)(double*,double*,int),double **x,double *y,int _datalen,double C_val)
    {
        dim=_dim;
        kennel=_kennel;
        datalen=_datalen,data=x,lable=y;
        Mat_dot = new double*[datalen];
        Mat_doty = new double*[datalen];
        for(int i=0; i<datalen; i++)Mat_dot[i]=new double[datalen];
        for(int i=0; i<datalen; i++)Mat_doty[i]=new double[datalen];
        for(int i=0; i<datalen; i++)
            for(int j=i; j<datalen; j++)
            {
                Mat_dot[i][j]=Mat_dot[j][i]=kennel(data[i],data[j],dim);
                Mat_doty[i][j]=Mat_doty[j][i]=Mat_dot[i][j]*lable[i]*lable[j];
            }
        Lamda = new double[datalen];
        double delat = 0;
        for(int i=0; i<datalen; i++)Lamda[i]=rand()%100/100.0,delat+=Lamda[i]*lable[i];
        for(int i=0; i<datalen; i++)
            if(delat*y[i]<0)
            {
                Lamda[i]=Lamda[i]-y[i]*delat;
                delat = 0;
                if(Lamda[i]>C_val)delat=(Lamda[i]-C_val)*y[i],Lamda[i]=C_val;
            }
        double res = 0;
        for(int i=0; i<datalen; i++)res+=Lamda[i]*y[i];
        this->C_val = C_val;
    }
    ~SVM()
    {
        for(int i=0; i<datalen; i++)
            delete Mat_dot[i],delete Mat_doty[i];
        delete Lamda;
        delete Mat_dot;
        double Mat_doty;
    }
    void train(int epochs)
    {
        while(epochs--)
        {
            for(int a=0; a<datalen; a++)
            {
                for(int b=0; b<datalen; b++)
                {
                    double down = Mat_dot[a][a]+Mat_dot[b][b]-2*Mat_dot[a][b];
                    if(down<1e-9)continue;
                    double sum_Lam_a = 0,sum_Lam_b=0;
                    for(int i=0; i<datalen; i++)
                        sum_Lam_a+=Lamda[i]*Mat_doty[i][a],sum_Lam_b+=Lamda[i]*Mat_doty[i][b];
                    sum_Lam_a-=Lamda[a]*Mat_doty[a][a]+Lamda[b]*Mat_doty[b][a];
                    sum_Lam_b-=Lamda[a]*Mat_doty[a][b]+Lamda[b]*Mat_doty[b][b];
                    double C = -Lamda[a]*lable[a]-Lamda[b]*lable[b];
                    double up =lable[a]*C*(Mat_dot[a][b]-Mat_dot[b][b])+1 - lable[a]*lable[b]-sum_Lam_a+lable[a]*lable[b]*sum_Lam_b;

                    double na =up/down,nb=-(C+na*lable[a])*lable[b],rate=1;
                    if(na<0)rate =min(rate,(0-Lamda[a])/(na-Lamda[a]));
                    if(nb<0)rate =min(rate,(0-Lamda[b])/(nb-Lamda[b]));
                    if(na>C_val)rate =min(rate,(C_val-Lamda[a])/(na-Lamda[a]));
                    if(nb>C_val)rate =min(rate,(C_val-Lamda[b])/(nb-Lamda[b]));
                    Lamda[a]+=(na-Lamda[a])*rate,Lamda[b]+=(nb-Lamda[b])*rate;
                }
                printf("<c_val>=%4.0f<a,b>=%d,**\r",C_val,a);
            }
        }


        double maxb=0,minb=0;
        for(int i=0; i<datalen; i++)
        {
            if(Lamda[i]<1e-9)continue;
            double t=0;
            for(int j=0; j<datalen; j++)t-=Lamda[j]*lable[j]*Mat_dot[j][i];
            if(i==0)maxb=minb=t;
            maxb=max(maxb,t);
            minb=min(minb,t);
        }
        b=(maxb+minb)/2;
    }
    double predict(double* x)
    {
        double res = b;
        for(int i =0; i<datalen; i++)
            res+=Lamda[i]*lable[i]*kennel(data[i],x,dim);
        return res>0?1:-1;
    }
    int predict(double **x,double *y,int _datalen)
    {
        int ans = 0;
        for(int i=0; i<_datalen; i++)
            ans+= abs(predict(x[i])-y[i])<1e-5?1:0;
        return ans;
    }
};
const double pi = acos(-1);
double guassKennel(double*x,double*y,int len)
{
    static double gamma = 1.0,down = 1/(sqrt(2)*gamma);
    double res = 0;
    for(int i=0; i<len; i++)res +=(x[i]-y[i])*(x[i]-y[i]);
    return down*exp(-res*down*down);
}
double LinearKennel(double*x,double*y,int len)
{
    double res = 0;
    for(int i=0; i<len; i++)res += x[i]*y[i];
    return res;
}
char ch[60000];
int p;
int readdouble(double&x)
{
    while(p<60000&!(ch[p]>='0'&&ch[p]<='9'))p++;
    if(p==60000)return 0;
    x=0;
    while(p<60000&&ch[p]>='0'&&ch[p]<='9')x=x*10+ch[p]-48,p++;
    if(ch[p]!='.')return 1;
    p++;
    double base=0.1;
    while(p<60000&&ch[p]>='0'&&ch[p]<='9')x+=(ch[p]-48)*base,base*=0.1,p++;
    return 1;
}

class ReadData
{
public:
    int lineNum,Maxlen;
    double **data,*label;
    ReadData(int dim,int maxlen,FILE* file)
    {
        Maxlen=maxlen;
        lineNum=0;
        label=new double[maxlen];
        data = new double*[maxlen];
        for(int i=0; i<maxlen; i++)data[i]=new double[dim];
        memset(ch,0,sizeof(ch));
        fread(ch,1,60000,file);
        p=0;
        while(readdouble(label[lineNum]))
        {
            label[lineNum]=label[lineNum]*2-1;
            for(int i=0; i<8; i++)readdouble(data[lineNum][i]);
            swap(data[lineNum][2],data[lineNum][7]);
            lineNum++;
        }
    }
    ~ReadData()
    {
        for(int i=0; i<Maxlen; i++)delete data[i];
        delete label;
        delete data;
    }
};
const double inf = 1e18;
int main()
{
    srand(time(NULL));


    FILE* train_file = fopen("classification/diabetes_train","r");
    FILE* test_file = fopen("classification/diabetes_test","r");
    ReadData train = ReadData(8,600,train_file);
    ReadData test = ReadData(8,200,test_file);
    fclose(test_file);
    fclose(train_file);

    for(int i=1; i<100; i++)
    {
        SVM* svm =new SVM(8,guassKennel,train.data,train.label,train.lineNum,i*500);
        for(int epoch=1;epoch<=5;epoch++){
            svm->train(1);
            printf("C=%d   epoch=%3d\t",500*i,epoch);
            printf("train %d/%d ",svm->predict(train.data,train.label,train.lineNum),train.lineNum);
            printf("test %d/%d\n",svm->predict(test.data,test.label,test.lineNum),test.lineNum);
        }
        delete svm;
    }
    return 0;
}

BP神经网络回归

class RegressionBPNet:
    def __init__(self, lr=0.5):
        self.l1 = Linear(1, 15, lr=lr)
        self.l2 = Linear(15, 1, lr=lr)

    def forward(self, x):
        x = self.l2(self.l1(x))
        return x.reshape((-1, x.shape[2]))

    def backward(self, grad):
        x = self.l1.backward(self.l2.backward(grad))

    def __call__(self, *args, **kwargs):
        return self.forward(*args)

def Regression():
    train_x, train_y = ReadData.getRegressionData('regression/sinc_train')
    test_x, test_y = ReadData.getRegressionData('regression/sinc_test')
    Bp = RegressionBPNet(lr=0.0015)

    maxy = np.max(train_y) + 0.2
    miny = np.min(train_y) - 0.2
    train_y = train_y / (maxy - miny) - miny / (maxy - miny)
    test_y = test_y / (maxy - miny) - miny / (maxy - miny)

    epochs = 14000
    a = []
    b = []
    for epoch in range(epochs):
        dy = Bp(train_x) - train_y
        Bp.backward(2 * dy)
        if (epoch + 1) % 100 == 0:
            loss = np.mean(dy * dy)
            test_dy = Bp(test_x) - test_y
            testloss = np.mean(test_dy * test_dy)
            a.append(loss)
            b.append(testloss)
            print('epoch=%5d loss %.12f testloss %.12f' % (epoch + 1, loss, testloss), end='\r')
    # plt.plot(a)
    # plt.plot(b)
    x = np.arange(-10, 10, 0.1)
    y = Bp(x)
    plt.plot(test_x, test_y, 'o')
    plt.plot(x, y)
    plt.show()


Regression()

多元线性回归

import numpy as np
import matplotlib.pyplot as plt
import ReadData


def linear_regression(x, y):
    w = np.linalg.inv(x.T @ x) @ x.T @ y
    return w


train_x, train_y = ReadData.getRegressionData('regression/sinc_train')
test_x, test_y = ReadData.getRegressionData('regression/sinc_test')

x = np.concatenate((train_x ** 8, train_x ** 6, train_x ** 4, train_x ** 2, np.ones([len(train_x), 1])), axis=1)
w = linear_regression(x, train_y)  # 求解模型
print(w)
# 绘制train数据
v = x @ w
plt.figure(1)
plt.plot(train_x[:, 0], train_y[:, 0], 'o', label='train_data')
plt.plot(train_x[:, 0], [_[0] for _ in v], 'o', label='predict')
plt.legend()
plt.show()
# 绘制test数据
x = np.concatenate((test_x ** 8, test_x ** 6, test_x ** 4, test_x ** 2, np.ones([len(test_x), 1])), axis=1)
v = x @ w
plt.figure(2)
plt.plot(test_x[:, 0], [_[0] for _ in test_y], 'o', label='test_data')
plt.plot(test_x[:, 0], [_[0] for _ in v], 'o', label='predict')
plt.legend()
plt.show()