首先，需要安装cvxopt模块。可以通过以下命令进行安装：
```
!pip install cvxopt
```
接下来，需要加载数据集。可以使用sklearn中的load_digits函数加载mnist数据集，代码如下：
```python
from sklearn.datasets import load_digits
digits = load_digits()
X = digits.data
y = digits.target
```
对于北理工的手写数据集，可以使用PIL库中的Image模块来读取图片，并使用numpy库将图片转换为数组。代码如下：
```python
from PIL import Image
import numpy as np
def read_image(file_path):
img = Image.open(file_path).convert('L')
img = img.resize((8, 8))
img_array = np.array(img)
img_array[img_array > 0] = 1
return img_array.reshape(1, -1)
X_test = read_image('test.png')
```
其中，'test.png'是北理工的手写数据集中的一张图片，需要将其放置在当前工作目录下。
接下来，可以使用cvxopt中的qp函数求解SVM的对偶问题。对于硬间隔SVM，其对偶问题可以表示为：
$$\max\limits_{\alpha}\sum_{i=1}^{n}\alpha_i-\frac{1}{2}\sum_{i=1}^{n}\sum_{j=1}^{n}y_iy_j\alpha_i\alpha_j\langle x_i,x_j\rangle$$
$$s.t.\ \sum_{i=1}^{n}y_i\alpha_i=0,\ \alpha_i\geq 0$$
其中，$\alpha$为拉格朗日乘子向量，$n$为样本数量，$x_i$为第$i$个样本的特征向量，$y_i\in\{-1,1\}$为第$i$个样本的标签。
代码如下：
```python
from cvxopt import matrix, solvers
# 构造P矩阵
n_samples, n_features = X.shape
K = np.dot(X, X.T)
P = matrix(np.outer(y, y) * K)
# 构造q、G、h、A、b向量
q = matrix(-np.ones((n_samples, 1)))
G = matrix(np.vstack((-np.eye(n_samples), np.eye(n_samples))))
h = matrix(np.hstack((np.zeros(n_samples), np.ones(n_samples) * C)))
A = matrix(y.reshape(1, -1))
b = matrix(np.zeros(1))
# 求解QP问题
sol = solvers.qp(P, q, G, h, A, b)
alphas = np.array(sol['x'])
# 计算权重向量w和截距b
w = np.dot(X.T, y * alphas)
b = y[alphas > 1e-5][0] - np.dot(X[alphas > 1e-5], w)
```
在求解QP问题之前，需要先计算样本之间的内积矩阵$K$，并构造P矩阵。同时，还需要构造G、h、A、b向量，用于表示拉格朗日乘子的约束条件。求解QP问题可以使用cvxopt中的qp函数。
接下来，可以使用求解得到的权重向量和截距，对测试数据集进行预测。代码如下：
```python
y_pred = np.sign(np.dot(X_test, w) + b)
```
对于mnist数据集，可以使用交叉验证来评估模型的性能。代码如下：
```python
from sklearn.model_selection import KFold
kf = KFold(n_splits=5, shuffle=True, random_state=42)
accs = []
for train_index, test_index in kf.split(X):
X_train, y_train = X[train_index], y[train_index]
X_test, y_test = X[test_index], y[test_index]
# 训练模型
...
# 预测测试集
y_pred = np.sign(np.dot(X_test, w) + b)
# 计算准确率
acc = np.mean(y_pred == y_test)
accs.append(acc)
print('Accuracy: %.2f%% (+/- %.2f%%)' % (np.mean(accs) * 100, np.std(accs) * 100))
```
对于北理工的手写数据集，可以分别对预处理和未预处理的数据进行预测，并比较两种情况的性能。代码如下：
```python
# 未预处理
X_test = read_image('test.png')
y_pred = np.sign(np.dot(X_test, w) + b)
print('Prediction: %d' % y_pred[0])
# 预处理
X_test = read_image('test.png')
X_test_processed = np.dot(X_test - np.mean(X_train, axis=0), pca.components_.T)
y_pred = np.sign(np.dot(X_test_processed, w_pca) + b_pca)
print('Prediction: %d' % y_pred[0])
```
其中，未预处理的代码与之前相同，预处理的代码需要先对测试数据进行降维，再进行预测。在这里，使用PCA降维，代码如下：
```python
from sklearn.decomposition import PCA
# 对训练数据进行降维
pca = PCA(n_components=32)
X_train_processed = np.dot(X_train - np.mean(X_train, axis=0), pca.components_.T)
# 训练模型
...
# 对测试数据进行降维
X_test_processed = np.dot(X_test - np.mean(X_train, axis=0), pca.components_.T)
# 预测测试集
y_pred = np.sign(np.dot(X_test_processed, w) + b)
# 计算准确率
acc = np.mean(y_pred == y_test)
print('Accuracy: %.2f%%' % (acc * 100))
```
需要注意的是，在对测试数据进行降维时，需要使用与训练数据相同的PCA模型。因此，在训练模型时需要保存PCA模型，代码如下：
```python
# 对训练数据进行降维
pca = PCA(n_components=32)
X_train_processed = np.dot(X_train - np.mean(X_train, axis=0), pca.components_.T)
# 训练模型
...
# 保存PCA模型
w_pca = np.dot(pca.components_.T, w)
b_pca = b + np.dot(np.mean(X_train, axis=0), w_pca)
# 对测试数据进行降维
X_test_processed = np.dot(X_test - np.mean(X_train, axis=0), pca.components_.T)
# 预测测试集
y_pred = np.sign(np.dot(X_test_processed, w_pca) + b_pca)
# 计算准确率
acc = np.mean(y_pred == y_test)
print('Accuracy: %.2f%%' % (acc * 100))
```
726
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
