我自己的原文哦~~https://blog.csdn.net/weixin_49587977/article/details/143161192

一、一维插值

对现有数据进行拟合或插值是数学分析中常见的方式。

  • 通过分析现有数据,得到一个连续的函数(也就是曲线);或者更密集的离散方程与已知数据互相吻合,这个过程叫做拟合
  • 通过已知的、离散的数据点,在范围内推求新数据点的过程或方法则叫做插值

简单来说,插值与拟合最大的区别就是,插值所获得的曲线一定要通过数据点,而拟合需要的是总体上最接近的结果。

w~python合集1_python

编辑

所以说要实现插值算法,我们的目标是通过给定的x值和y值,创建一个函数y=f(x),可以在该函数中插入想要的任何值a并获得相应的值y=f(a)。

例如给定x=[0:5:5],y=x^2,在python中画出散点图为

import numpy as np
import matplotlib.pyplot as plt
x_data = np.linspace(0,5,5)
y_data = x_data**2
plt.scatter(x_data,y_data)
plt.show()


w~python合集1_python_02

编辑

 假如我们想获取x=2的值,可以看到图中没有对应的数据点,因此无法直接获得,所以需要用到插值,而最简单的插值方式为线性插值,就是通过直线连接数据点。

w~python合集1_python_03

编辑

这样通过scipy.interpolate.interp1d()的函数即可进行一维插值。

from scipy.interpolate import interp1d
x_data = np.linspace(0,5,5)
y_data = x_data**2
y_f = interp1d(x_data, y_data, 'linear')
print(y_f(2))


通过上述代码使用线性插值,运行代码后求出x=2时y的值为4.375,由于线性插值是通过相邻两点的连线来进行插值,无法考虑到其他数据点,我们都知道x=2时y=x^2的值为4,这里的插值结果明显有问题。

于是使用2次多项式的方式进行插值,可以得到y(2)=4

y_f = interp1d(x_data, y_data,'quadratic')
x = np.linspace(0,5,100)
y = y_f(x)
plt.plot(x,y,'--')
plt.show()
print(y_f(2))


w~python合集1_python_04

上述案例的数据较为简单,我们可以轻易地判断插值后的数据是否准确,但是对于复杂的数据,选择合适的插值方式就变得十分重要。 




二、二维插值

一维插值就是一维插值就是给出y=f(x)上的点(x1,y1),(x2,y2),…,(xn,yn),由此求出y=f(x)在点xa处的值ya的值。

在此基础上进行扩展,二维插值就是给出z=f(x,y)上的点(x1,y1,z1),…,(xn,yn,zn),由此求出在(xa,ya)处求出za的值。

使用Scipy中的interpolate.interp2d函数可以实现二维插值。

class scipy.interpolate.interp2d(x, y, z, kind='linear', copy=True, 
bounds_error=False, fill_value=None)


定义数据点坐标的数组。如果这些点位于规则网格上,x 可以指定列坐标,y 可以指定行坐标,例如:

x = [0,1,2];  y = [0,3]; z = [[1,2,3], [4,5,6]]


否则,x 和 y 必须指定每个点的完整坐标,例如:

x = [0,1,2,0,1,2];  y = [0,0,0,3,3,3]; z = [1,4,2,5,3,6]


如果 x 和 y 是多维的数组,则会在使用前将它们展平。

我们使用以下代码绘图进行示例:

import numpy as np
from scipy.interpolate import interp2d
import matplotlib.pyplot as plt
x = np.linspace(0, 4, 13)
y = np.array([0, 2, 3, 3.5, 3.75, 3.875, 3.9375, 4])
X, Y = np.meshgrid(x, y)
Z = np.sin(np.pi*X/6) * np.exp(Y/3)
x2 = np.linspace(0, 4, 65)
y2 = np.linspace(0, 4, 65)
f1 = interp2d(x, y, Z, kind='cubic')
Z2 = f1(x2, y2)
f2 = interp2d(x,y,Z,kind='linear')
Z3 = f2(x2,y2)
fig, ax = plt.subplots(nrows=1, ncols=3)
ax[0].pcolormesh(X, Y, Z,cmap='rainbow')
X2, Y2 = np.meshgrid(x2, y2)
ax[1].pcolormesh(X2, Y2, Z2,cmap='rainbow')
X3, Y3 = np.meshgrid(x2, y2)
ax[2].pcolormesh(X3, Y3, Z3,cmap='rainbow')
ax[0].set_title('origin')
ax[1].set_title('cubic')
ax[2].set_title('linear')
plt.show()


可以通过图片对比各插值方式的区别。

w~python合集1_python_05





三、几种一维插值方式对比

使用python的scipy库可进行一维插值,并且对比了线性(linear)和二次多项式(quadratic)两种方法,这里主要来介绍一下scipy.interpolate.interp1d函数的其他插值选项。

首先,interp1d函数由以下部分所组成

class scipy.interpolate.interp1d(x, y, kind='linear', axis=- 1, 
copy=True, bounds_error=None,fill_value=nan, assume_sorted=False)


  • x为一维数组
  • y可以是N维数组,但y 沿插值轴的长度必须等于 x 的长度
  • kind表示插值的方式(默认为‘linear’),包括 ‘linear’, ‘nearest’, ‘nearest-up’, ‘zero’, ‘slinear’, ‘quadratic’, ‘cubic’, ‘previous’, or ‘next’.
    “zero”、“linear”、“quadratic”和“cubic”是指零、一、二或三阶的样条插值;
    ‘previous’ 和 ‘next’ 只是返回该点的上一个或下一个值
    'nearest-up' 和 'nearest' 在对半整数(例如 0.5、1.5)进行插值时有所不同,“nearest-up”向上取整,“nearest”向下取整。
  • axis用于指定沿其进行插值的 y 轴。默认为最后一个y轴。
  • copy(布尔类型)如果为 True,则该类函数创建 x 和 y 的内部副本。如果为 False,则使用对 x 和 y 的引用。默认是True。
  • bounds_error(布尔类型)如果为 True,则在任何时候尝试对 x 范围之外的值进行插值时都会引发 ValueError(需要外插)。如果为 False,则为越界值分配 fill_value。默认情况下会报错,除非fill_value="extrapolate"。
  • fill_value,如果是 ndarray(或浮点数),则此值将用于填充数据范围之外的请求点。如果未提供,则默认值为 NaN。类数组必须正确广播到非插值轴的尺寸。如果是双元素元组,则第一个元素用作 x_new < x[0] 的填充值,第二个元素用于 x_new > x[-1]。任何不是 2 元素元组的东西(例如,列表或 ndarray,无论形状如何)都被视为一个类似数组的参数,用于两个边界,如下所示,上面 = fill_value,fill_value。使用二元素元组或 ndarray 需要 bounds_error=False。

我们使用下述代码绘图进行对比

import numpy as np
from scipy.interpolate import interp1d
import pylab
A, nu, k = 10, 4, 2
def f(x, A, nu, k):
    return A * np.exp(-k*x) * np.cos(2*np.pi * nu * x)
xmax, nx = 0.5, 8
x = np.linspace(0, xmax, nx)
y = f(x, A, nu, k)
f_nearest = interp1d(x, y,'nearest')
f_linear  = interp1d(x, y,'linear')
f_cubic   = interp1d(x, y,'cubic')
f_next    = interp1d(x, y, 'next')
x2 = np.linspace(0, xmax, 100)
pylab.plot(x, y, 'o', label='data points')
pylab.plot(x2, f(x2, A, nu, k), label='exact')
pylab.plot(x2, f_nearest(x2), label='nearest')
pylab.plot(x2, f_linear(x2), label='linear')
pylab.plot(x2, f_cubic(x2), label='cubic')
pylab.plot(x2, f_next(x2), label='next')
pylab.legend(loc=1)
pylab.show()


w~python合集1_python_06

编辑

可以通过上图对比各插值方式与原函数(exact)的曲线区别。




四、曲线拟合

在处理数据时经常需要进行曲线拟合,在拟合过程中我们采用了插值的思想。

对于给定的数据x=[...]和y=[...],插值的目的是找到参数β的最优集合,使函数

w~python合集1_python_07

编辑

能够与原数据最为相似。

  • 其中一种方法是通过调整β使:

  • 能够最小化,这种方法称之为最小二乘法。
  • 如果yi的值有相应的误差,那么使下式最小化
  • 称之为β的最大似然估计。在给定xi与yi的情况下,这样得到的β值最为准确。

我们给定以下随机数据点并绘制散点图:

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from scipy.optimize import curve_fit
x_data = np.array([0.        , 0.15789474, 0.31578947, 0.47368421, 0.63157895,
       0.78947368, 0.94736842, 1.10526316, 1.26315789, 1.42105263,
       1.57894737, 1.73684211, 1.89473684, 2.05263158, 2.21052632,
       2.36842105, 2.52631579, 2.68421053, 2.84210526, 3.        ])
y_data = np.array([  2.95258285,   3.49719803,  -2.1984975 ,  -4.88744346,
        -7.41326345,  -8.44574157, -10.01878504, -13.83743553,
       -12.91548145, -16.41149046, -14.93516299, -13.42514157,
       -14.12110495, -17.6412464 , -16.1275509 , -17.11533771,
       -15.66076021, -12.48938865, -11.33918701, -11.70467566])
plt.scatter(x_data,y_data)
plt.show()


w~python合集1_python_08

编辑

 我们使用下式的模型进行拟合:

w~python合集1_python_09

编辑

我们需要求解出符合上述数据点的最优的a、b、c的值,需要在python中:

(1)需要定义模型方程;

(2)使用scipy库中的curve_fit函数,使用该函数时需要定义一个β的初始值,对于一些复杂的模型,初始值决定了函数能否顺利运算。

下面定义函数并进行拟合

def fitfun(x, a, b, c):
    return a*(x-b)**2 + c
popt, pcov = curve_fit(fitfun, x_data, y_data, p0=[3,2,-16])


其中:

popt为我们定义的函数fitfun的最佳参数。

pcov为协方差矩阵,用于给定误差。

接下来我们就可以绘制出拟合好的曲线

a, b, c = popt
x_model = np.linspace(min(x_data), max(x_data), 100)
y_model = fitfun(x_model, a, b, c)
plt.scatter(x_data,y_data)
plt.plot(x_model,y_model, color='r')
plt.show()


w~python合集1_python_10





五、python支持向量机

本文将首先推导SVM的主要公式,接着基于Platt-SMO算法,从零开始实现支持多种核函数的SVM,然后基于One-Versus-One策略实现多分类,最后在MNIST和CIFAR-10数据集上进行性能测试。从零开始实现支持向量机

这篇文章诞生于机器学习课程无聊的大作业,既然已经为此浪费了不少时间,不妨再多花点时间写一篇文章,借此记录一下实现过程。支持向量机的数学形式简约而直观,但一旦涉及具体实现,各种问题就会接踵而来。在本文中,笔者将首先推导SVM的主要公式,接着基于Platt-SMO算法,从零开始实现支持多种核函数的SVM,然后基于One-Versus-One策略实现多分类,最后在MNIST和CIFAR-10数据集上进行性能测试

数学推导基本形式

w~python合集1_python_11

编辑

w~python合集1_python_12

编辑

w~python合集1_python_13

编辑

这是一个二次规划问题,我们可以采用梯度下降或者坐标下降等方法求解。

核函数与核技巧

w~python合集1_python_14

编辑

Platt-SMO算法

w~python合集1_python_15

编辑

w~python合集1_python_16

编辑

w~python合集1_python_17

编辑

算法实现核函数

我们对上述四种核函数进行实现,这里将核函数封装成类,通过实现__call__方法,使其实例可以像函数一样被调用

线性核代码如下

class LinearKernel(object):
    def __init__(self):
        self.name = 'linear'
 
    def __call__(self, X, y):
        return X @ y.T


多项式核代码如下

class PolynomialKernel(object):
    def __init__(self, gamma=1.0, degree=3):
        self.name = 'polynomial'
        self.gamma = gamma
        self.degree = degree
    
    def __call__(self, X, y):
        return np.power(self.gamma * (X @ y.T) + 1, self.degree)


高斯核代码如下

class GaussianKernel(object):
    def __init__(self, gamma=1.0):
        self.name = 'gaussian'
        self.gamma = gamma
    
    def __call__(self, X, y):
        return np.exp(-self.gamma * np.sum(np.square(X - y), axis=1))


Sigmoid核代码如下

class SigmoidKernel(object):
    def __init__(self, gamma=1.0, bias=0.0):
        self.name = 'sigmoid'
        self.gamma = gamma
        self.bias = bias
    
    def __call__(self, X, y):
        return np.tanh(self.gamma * (X @ y.T) + self.bias)


另外,我们定义一个工具函数,方便核函数的创建

def CreateKernel(entry):
    if entry['name'] == 'linear':
        return LinearKernel()
    elif entry['name'] == 'polynomial':
        return PolynomialKernel(entry['gamma'], entry['degree'])
    elif entry['name'] == 'gaussian':
        return GaussianKernel(entry['gamma'])
    elif entry['name'] == 'sigmoid':
        return SigmoidKernel(entry['gamma'], entry['bias'])
    raise AttributeError('invalid kernel')


支持向量机

参考scikit-learn的封装,我们定义一个类,提供fit和predict两种方法,参数包括最大迭代次数、惩罚系数、误差精度和核函数类型,利用私有函数实现  和  的选择和单步更新,对于线性核,我们提供weight属性,用于获取线性核的分类超平面参数,除了一些简化以外,代码基本按照Platt-SMO算法进行实现

class SupportVectorMachine(object):
    def __init__(self, iteration=100, penalty=1.0, epsilon=1e-6, kernel=None):
        self.iteration = iteration
        self.penalty = penalty
        self.epsilon = epsilon
        if kernel is None:
            kernel = {'name': 'linear'}
        self.kernel = CreateKernel(kernel)
    
    def __compute_w(self):
        return (self.a * self.y) @ self.X

    def __compute_e(self, i):
        return (self.a * self.y) @ self.K[:, i] + self.b - self.y[i]
    
    def __select_j(self, i):
        j = np.random.randint(1, self.m)
        return j if j > i else j - 1
    
    def __step_forward(self, i):
        e_i = self.__compute_e(i)
        if ((self.a[i] > 0) and (e_i * self.y[i] > self.epsilon)) or ((self.a[i] < self.penalty) and (e_i * self.y[i] < -self.epsilon)):
            j = self.__select_j(i)
            e_j = self.__compute_e(j)
            a_i, a_j = np.copy(self.a[i]), np.copy(self.a[j])
            if self.y[i] == self.y[j]:
                L = max(0, a_i + a_j - self.penalty)
                H = min(self.penalty, a_i + a_j)
            else:
                L = max(0, a_j - a_i)
                H = min(self.penalty, self.penalty + a_j - a_i)
            if L == H:
                return False
            d = 2 * self.K[i, j] - self.K[i, i] - self.K[j, j]
            if d >= 0:
                return False
            self.a[j] = np.clip(a_j - self.y[j] * (e_i - e_j) / d, L, H)
            if np.abs(self.a[j] - a_j) < self.epsilon:
                return False
            self.a[i] = a_i + self.y[i] * self.y[j] * (a_j - self.a[j])
            b_i = self.b - e_i - self.y[i] * self.K[i, i] * (self.a[i] - a_i) - self.y[j] * self.K[j, i] * (self.a[j] - a_j)
            b_j = self.b - e_j - self.y[i] * self.K[i, j] * (self.a[i] - a_i) - self.y[j] * self.K[j, j] * (self.a[j] - a_j)
            if 0 < self.a[i] < self.penalty:
                self.b = b_i
            elif 0 < self.a[j] < self.penalty:
                self.b = b_j
            else:
                self.b = (b_i + b_j) / 2
            return True
        return False
    
    def setup(self, X, y):
        self.X, self.y = X, y
        self.m, self.n = X.shape
        self.b = 0.0
        self.a = np.zeros(self.m)
        self.K = np.zeros((self.m, self.m))
        for i in range(self.m):
            self.K[:, i] = self.kernel(X, X[i, :])
    
    def fit(self, X, y):
        self.setup(X, y)
        entire = True
        for _ in range(self.iteration):
            change = 0
            if entire:
                for i in range(self.m):
                    change += self.__step_forward(i)
            else:
                index = np.nonzero((0 < self.a) * (self.a < self.penalty))[0]
                for i in index:
                    change += self.__step_forward(i)
            if entire:
                entire = False
            elif change == 0:
                entire = True

    def predict(self, X):
        m = X.shape[0]
        y = np.zeros(m)
        for i in range(m):
            y[i] = np.sign((self.a * self.y) @ self.kernel(self.X, X[i, :]) + self.b)
        return y
    
    @property
    def weight(self):
        if self.kernel.name != 'linear':
            raise AttributeError('non-linear kernel')
        return self.__compute_w(), self.b


多分类

基于One-Versus-One策略, 我们构造  个SVM, 其中  为类别数, 训练每个分类器时, 选取相应类别的样本作为训练集, 并将标签映射到 -1 和 1 , 在预测时, 用每个分类器的预测结果进行投票, 从而得到最终结果

我们采用与支持向量机完全相同的封装,提供fit和predict两种方法,使该类成为通用的分类模型

class SupportVectorClassifier(object):
    def __init__(self, iteration=100, penalty=1.0, epsilon=1e-6, kernel=None):
        self.iteration = iteration
        self.penalty = penalty
        self.epsilon = epsilon
        self.kernel = kernel
        self.classifier = []

    def __build_model(self, y):
        self.label = np.unique(y)
        for i in range(len(self.label)):
            for j in range(i+1, len(self.label)):
                model = SupportVectorMachine(self.iteration, self.penalty, self.epsilon, self.kernel)
                self.classifier.append((i, j, model))

    def fit(self, X, y):
        self.__build_model(y)
        for i, j, model in tqdm(self.classifier):
            index = np.where((y == self.label[i]) | (y == self.label[j]))[0]
            X_ij, y_ij = X[index], np.where(y[index] == self.label[i], -1, 1)
            model.fit(X_ij, y_ij)
    
    def predict(self, X):
        vote = np.zeros((X.shape[0], len(self.label)))
        for i, j, model in tqdm(self.classifier):
            y = model.predict(X)
            vote[np.where(y == -1)[0], i] += 1
            vote[np.where(y == 1)[0], j] += 1
        return self.label[np.argmax(vote, axis=1)]


性能测试

首先,我们在二维平面上构造两组简单的正态分布数据,用于可视化支持向量机的分类效果,首先构造数据并训练模型

X = np.concatenate((np.random.randn(500, 2) - 2, np.random.randn(500, 2) + 2))
y = np.concatenate((np.ones(500), -np.ones(500)))
C = SupportVectorMachine(iteration=100)
C.fit(X, y)
w, b = C.weight
u = np.linspace(-3, 3, 100)
v = (-b - w[0] * u) / w[1]


然后根据模型参数绘制分类效果

plt.scatter(X[:500, 0], X[:500, 1], label='Positive')
plt.scatter(X[500:, 0], X[500:, 1], label='Negative')
plt.plot(u, v, label='Separation', c='g')
plt.xlabel('$x$')
plt.ylabel('$y$')
plt.title('Separation Sample')
plt.grid()
plt.legend()
plt.tight_layout()
plt.savefig('./figure/separation.png')
plt.show()


可以看到,我们实现的SVM可以很好地将两组数据分开

w~python合集1_python_18

编辑

w~python合集1_python_19

编辑

def MNIST(path, group='train'):
    if group == 'train':
        with gzip.open(os.path.join(path, 'train-images-idx3-ubyte.gz'), 'rb') as file:
            image = np.frombuffer(file.read(), np.uint8, offset=16).reshape(-1, 1, 28, 28) / 255.0
        with gzip.open(os.path.join(path, 'train-labels-idx1-ubyte.gz'), 'rb') as file:
            label = np.frombuffer(file.read(), np.uint8, offset=8)
    elif group == 'test':
        with gzip.open(os.path.join(path, 't10k-images-idx3-ubyte.gz'), 'rb') as file:
            image = np.frombuffer(file.read(), np.uint8, offset=16).reshape(-1, 1, 28, 28) / 255.0
        with gzip.open(os.path.join(path, 't10k-labels-idx1-ubyte.gz'), 'rb') as file:
            label = np.frombuffer(file.read(), np.uint8, offset=8)
    remain = 500 if group == 'train' else 100
    image_list, label_list = [], []
    for value in range(10):
        index = np.where(label == value)[0][:remain]
        image_list.append(image[index])
        label_list.append(label[index])
    image, label = np.concatenate(image_list), np.concatenate(label_list)
    index = np.random.permutation(len(label))
    return image[index], label[index]


对于CIFAR10数据集,我们做同样的处理

def CIFAR10(path, group='train'):
    if group == 'train':
        image_list, label_list = [], []
        for i in range(1, 6):
            filename = os.path.join(path, 'data_batch_{}'.format(i))
            with open(filename, 'rb') as file:
                data = pickle.load(file, encoding='bytes')
            image_list.append(np.array(data[b'data'], dtype=np.float32).reshape(-1, 3, 32, 32) / 255.0)
            label_list.append(np.array(data[b'labels'], dtype=np.int32))
        image, label = np.concatenate(image_list), np.concatenate(label_list)
    elif group == 'test':
        filename = os.path.join(path, 'test_batch')
        with open(filename, 'rb') as file:
            data = pickle.load(file, encoding='bytes')
        image = np.array(data[b'data'], dtype=np.float32).reshape(-1, 3, 32, 32) / 255.0
        label = np.array(data[b'labels'], dtype=np.int32)
    remain = 500 if group == 'train' else 100
    image_list, label_list = [], []
    for value in range(10):
        index = np.where(label == value)[0][:remain]
        image_list.append(image[index])
        label_list.append(label[index])
    image, label = np.concatenate(image_list), np.concatenate(label_list)
    index = np.random.permutation(len(label))
    return image[index], label[index]


由于CIFAR10数据集较为困难,我们考虑利用CV方法进行特征提取,这里我们使用HOG特征提高分类效果,首先将彩色图像转换为灰度图像

def RGB2Gray(image):
    image = 0.299 * image[0] + 0.587 * image[1] + 0.114 * image[2]
    return image.reshape(1, *image.shape)


然后实现一个简单的HOG特征提取函数,这里我们没有实现区块重叠,对该函数进行改进应该可以进一步提高分类效果   

def HOG(image, block=4, partition=8):
    image = RGB2Gray(image).squeeze(axis=0)
    height, width = image.shape
    gradient = np.zeros((2, height, width), dtype=np.float32)
    for i in range(1, height-1):
        for j in range(1, width-1):
            delta_x = image[i, j-1] - image[i, j+1]
            delta_y = image[i+1, j] - image[i-1, j]
            gradient[0, i, j] = np.sqrt(delta_x ** 2 + delta_y ** 2)
            gradient[1, i, j] = np.degrees(np.arctan2(delta_y, delta_x))
            if gradient[1, i, j] < 0:
                gradient[1, i, j] += 180
    unit = 360 / partition
    vertical, horizontal = height // block, width // block
    feature = np.zeros((vertical, horizontal, partition), dtype=np.float32)
    for i in range(vertical):
        for j in range(horizontal):
            for k in range(block):
                for l in range(block):
                    rho = gradient[0, i*block+k, j*block+l]
                    theta = gradient[1, i*block+k, j*block+l]
                    index = int(theta // unit)
                    feature[i, j, index] += rho
            feature[i, j] /= np.linalg.norm(feature[i, j]) + 1e-6
    return feature.reshape(-1)


基于这些工具函数,我们可以优雅地完成图像分类任务,对于MNIST数据集,一个基于线性核的分类示例如下

X_train, y_train = MNIST('./dataset/mnist_data/', group='train')
X_test, y_test = MNIST('./dataset/mnist_data/', group='test')
X_train, X_test = X_train.reshape(-1, 28*28), X_test.reshape(-1, 28*28)

model = SupportVectorClassifier(iteration=100, penalty=0.05)
model.fit(X_train, y_train)
p_train, p_test = model.predict(X_train), model.predict(X_test)

r_train, r_test = ComputeAccuracy(p_train, y_train), ComputeAccuracy(p_test, y_test)
print('Kernel: Linear, Train: {:.2%}, Test: {:.2%}'.format(r_train, r_test))


对于CIFAR10数据集,一个基于HOG特征和高斯核的分类示例如下

X_train, y_train = CIFAR10('./dataset/cifar-10-batches-py/', group='train')
X_test, y_test = CIFAR10('./dataset/cifar-10-batches-py/', group='test')
X_train, X_test = BatchHOG(X_train, partition=16), BatchHOG(X_test, partition=16)

kernel = {'name': 'gaussian', 'gamma': 0.03}
model = SupportVectorClassifier(iteration=100, kernel=kernel)
model.fit(X_train, y_train)
p_train, p_test = model.predict(X_train), model.predict(X_test)

r_train, r_test = ComputeAccuracy(p_train, y_train), ComputeAccuracy(p_test, y_test)
print('Kernel: Gaussian, Train: {:.2%}, Test: {:.2%}'.format(r_train, r_test))


经过测试,我们实现的SVM分类器在MNIST和CIFAR10数据集上的分类精度如下表所示

w~python合集1_python_20

编辑

此外,我们对模型的收敛性和各个核函数的参数选择进行了测试,模型精度与迭代次数的关系如下图所示 

w~python合集1_python_21

编辑

w~python合集1_python_22

编辑

w~python合集1_python_23

编辑

w~python合集1_python_24

编辑

w~python合集1_python_25

编辑

上述结果揭示了各个参数对模型性能的影响,可以为调参提供一定的指导作用

写在最后

SVM从过去的炙手可热到如今的日薄西山,仅仅过去了十年的时间,无论是精度还是效率,SVM都完败于当下随处可见的神经网络,关于从零开始实现SVM的意义,我也感到迷茫,但这一过程或多或少改变了我对机器学习的认知,一个简洁优雅的多项式时间精确算法,也许只能满足理论研究者的洁癖,而优化复杂模型的近似算法,在工程上赢得了未来。





六、Python真正多线程要来了

「Python 中的 GIL 将不复存在,这是人工智能生态系统领域中的巨大胜利。」PyTorch 核心维护者 Dmytro Dzhulgakov 感慨道。

GIL 是什么?GIL 的全称是 Global Interpreter Lock(全局解释器锁),它不是 Python 独有的,而是在实现 CPython(Python 解释器)时引入的一个概念。我们可以将 GIL 理解为一个互斥锁,用来保护 Python 里的对象,防止同一时刻多个线程执行 Python 的字节码,从而确保线程安全。

w~python合集1_python_26

然而,GIL 存在一个弊端,即在同一时刻只能有一个线程在一个 CPU 上执行,无法将多个线程映射到多个 CPU 上,使得 Python 并不能实现真正的多线程并发,从而降低了执行效率。

现在,Python 团队已经正式接受了删除 GIL 的这个提议,并将其设置为可选模式,可谓是利好广大开发者。

做出这一贡献的是一位来自 Meta 的名叫 Sam Gross 的软件工程师,他花费了四年多的时间才完成这一工程。

在得知这一消息后,大家纷纷叫好,深度学习三巨头之一的 Yann LeCun 发文祝贺:没有了 GIL,现在,Python 代码可以自由的执行多线程了。

「Python 中终于没有 GIL 了!」

「这是一个里程碑式的决定,是编码社区所热切期待的。」

具体细节如何,我们接着看下文。

CPython 核心开发者 Thomas Wouters 撰文描述了 Python 中的无 GIL 细节,并对未来发展做了展望。

非常感谢所有人对无 GIL 提议的反馈,整体上都持积极的支持态度。指导委员会打算接受无 GIL 提议,并就以下具体细节与大家分享。

我们的基本设想是:

  • 长期来看(大约 5 年以上),no-GIL 构建应是唯一的构建;
  • 我们希望非常谨慎地向后兼容。我们不希望出现另一个 Python 3 的情况,所有适应 no-GIL 构建所需的任何第三方代码更改应只适用于 with-GIL 构建(尽管仍要解决更老 Python 版本的向后兼容性问题)。这不适用于 Python 4。我们仍在考虑对这两个构建的 ABI 兼容性和其他细节的要求,以及对向后兼容性的影响;
  • 在我们承诺完全转向 no-GIL 之前,需要看到社区的支持。我们不能只是更改默认设置,更希望社区弄清自己需要做什么工作来给予支持。我们核心开发团队需要获得新构建模式及相关所有内容的经验。我们要整理现有代码中的线程安全性,因而需要弄明白新的 C API 和 Python API。我们在获得这些洞见时还需要传达给 Python 社区的其他人,并确保自身想要做出的更改以及希望他们做出的更改是可取的;
  • 在我们默认 no-GIL 设置之前的任何时候,如果事实证明了,它的破坏性太大导致收益太少,我们希望能够改变主意。这也就意味着我们会回滚所有工作,因此在我们确定要将 no-GIL 设为默认方式之前,特定于 no-GIL 的代码在某种程度上应是可识别的。

目前,我们认为未来的道路分为以下三个阶段:

  • 短期内,我们会将 no-GIL 构建作为一种实验性构建模式,大概是在 3.13 版本(也有可能推迟到 3.14 版本)。之所以是实验性的,是因为我们核心开发团队虽然支持这一构建模式,但不期望整个社区都会支持它。我们需要时间弄清自己要做什么,至少在 API 设计以及打包和分发方面,从而得到社区的支持。我们也不鼓励 distributor 将实验性 no-GIL 构建作为默认解释器发布。
  • 中期来看,在我们确信得到足够的社区支持并使 no-GIL 的生产使用可行后,我们将支持 no-GIL 构建,但不是默认方式,而是在某个目标日期或某个 Python 版本中使它成为默认方式。具体的时间将取决于很多因素,比如 API 更改最终兼容性如何、社区认为他们仍然需要做多少工作等。我们预计这至少需要一至两年的时间。一旦我们宣布支持,预计将有一些 distributor 会开始默认发布 no-GIL。
  • 长期来看,我们希望 no-GIL 成为默认方式,并删除 GIL 的所有痕迹(但不会不必要地破坏向后兼容性)。我们不希望等待太长时间,毕竟两种常用的构建模式同时存在会给社区造成很大的负担(比如需要双倍测试资源和 debug 场景)。但是我们也不能急于求成。我们认为这一过程将需要花费五年的时间。

当然在整个过程中,我们整个开发团队将需要实时评估进程并对时间线进行调整。

评论区的小伙伴们,你们对 GIL 成为可选是什么看法呢?

参考链接:

https://discuss.python.org/t/a-steering-council-notice-about-pep-703-making-the-global-interpreter-lock-optional-in-cpython/30474

https://twitter.com/dzhulgakov/status/1685667015800066048