NNDL 实验三线性回归

最新推荐文章于 2022-10-14 22:00:18 发布

红肚兜

最新推荐文章于 2022-10-14 22:00:18 发布

阅读量498

点赞数 4

文章标签：线性回归 python 机器学习

本文链接：https://blog.csdn.net/weixin_63539289/article/details/126871666

版权

2.2 线性回归

这次作业用到了较多的数据可视化方面的东西，不了解的同学可以去看一看这一块的知识点。这里附一个链接，这位博主讲的很详细。Gingkens

2.2.1 数据库构建

构造一个小的回归数据集：

生成 150 个带噪音的样本，其中 100 个训练样本，50 个测试样本，并打印出训练数据的可视化分布。

import torch as tc
from matplotlib import pyplot as plt  # matplotlib 是 Python 的绘图库
from random import *


def linear_func(x, w=1.2, b=0.5):
    y = w * x + b
    return y


def create_toy_data(func, interval, sample_num, noise, add_outlier, outlier_ratio):
    """
    根据给定的函数，生成样本
    输入：
       - func：函数
       - interval： x的取值范围
       - sample_num： 样本数目
       - noise： 噪声均方差
       - add_outlier：是否生成异常值
       - outlier_ratio：异常值占比
    输出：
       - X: 特征数据，shape=[n_samples,1]
       - y: 标签数据，shape=[n_samples,1]
    """
    # 均匀采样
    # 使用torch.rand在生成sample_num个随机数，这些数范围在interval[0]到interval[1]之间，并组成一个张量
    X = tc.rand(sample_num) * (interval[1] - interval[0]) + interval[0]
    y = func(X)

    # 生成高斯分布的标签噪声
    # 使用torch.normal生成0均值，noise标准差的数据
    epsilon = tc.normal(0, noise, size=y.shape)
    y = y + epsilon

    if add_outlier:  # 生成额外的异常点
        outlier_num = int(len(y) * outlier_ratio)
        if outlier_num != 0:
            # 使用torch.randint生成服从均匀分布的、范围在[0, len(y))的随机Tensor
            outlier_idx = tc.randint(outlier_num, (0, len(y)))
            y[outlier_idx] = y[outlier_idx] * 5
    return X, y


func = linear_func
interval = (-10, 10)
train_num = 100  # 训练样本数目
test_num = 50  # 测试样本数目
noise = 2
add_outlier = False
outlier_ratio = 0.001
X_train, y_train = create_toy_data(func, interval, train_num, noise, add_outlier, outlier_ratio)
X_test, y_test = create_toy_data(func, interval, test_num, noise, add_outlier, outlier_ratio)
X_train_large, y_train_large = create_toy_data(func, interval, 5000, noise, add_outlier, outlier_ratio)

# paddle.linspace返回一个Tensor，Tensor的值为在区间start和stop上均匀间隔的num个值，输出Tensor的长度为num
X_underlying = tc.linspace(interval[0], interval[1], train_num)  # 从interval[0]开始，到interval[1]结束，均匀生成train_num个点组成一个数组
y_underlying = linear_func(X_underlying)

# 绘制数据
plt.scatter(X_train, y_train, marker='*', color='r', s=20, label="train data")
plt.scatter(X_test, y_test, color='y', s=20, label="test data")
plt.plot(X_underlying, y_underlying, c='k', label=r"underlying distribution")
plt.legend(fontsize='x-large')  # 给图像加图例
# plt.savefig('ml-vis.pdf')  # 保存图像到PDF文件中
plt.show()

运行结果如下：
在这里插入图片描述
红色星星代表的是训练集，黄色圆圈代表的是测试集。可以看到，训练集在直线上下浮动，测试集也在直线上下浮动，这正是我们想要的结果。

2.2.2 模型构建

在这里插入图片描述
在实践中，为了提高预测样本的效率，我们通常会将N样本归为一组进行成批地预测，这样可以更好地利用GPU设备的并行计算能力。

其中X为N个样本的特征矩阵，y为N个预测值组成的列向量。

注意：在实践中，样本的矩阵X是由N个x的行向量组成。

import torch as tc

tc.manual_seed(11)  # 设置随机种子


# 线性算子
class Linear():
    def __init__(self, input_size):
        """
        输入：
           - input_size:模型要处理的数据特征向量长度
        """

        self.input_size = input_size

        # 模型参数
        self.params = {}
        self.params['w'] = tc.randn(self.input_size, 1, dtype=tc.float32)
        self.params['b'] = tc.zeros(1, dtype=tc.float32)

    def __call__(self, X):
        return self.forward(X)

    # 前向函数
    def forward(self, X):
        """
        输入：
           - X: tensor, shape=[N,D]
           注意这里的X矩阵是由N个x向量的转置拼接成的，与原教材行向量表示方式不一致
        输出：
           - y_pred： tensor, shape=[N]
        """
        N, D = X.shape

        if self.input_size == 0:
            return tc.full((N, 1), fill_value=self.params['b'])

        assert D == self.input_size  # 输入数据维度合法性验证
        # 使用paddle.matmul计算两个tensor的乘积
        y_pred = tc.matmul(X, self.params['w']) + self.params['b']
        return y_pred


input_size = 3
N = 2
X = tc.randn(N, input_size, dtype=tc.float32)  # 生成2个维度为3的数据
model = Linear(input_size)
y_pred = model(X)
print("y_pred:", y_pred)  # 输出结果的个数也是2个

输出结果如下：
在这里插入图片描述

2.2.3 损失函数

回归任务是对连续值的预测，希望模型能根据数据的特征输出一个连续值作为预测值。因此回归任务中常用的评估指标是均方误差。
在这里插入图片描述

import torch as tc


def mean_squared_error(y_true, y_pred):
    """
    输入：
       - y_true: tensor，样本真实标签
       - y_pred: tensor, 样本预测标签
    输出：
       - error: float，误差值
    """

    assert y_true.shape[0] == y_pred.shape[0]

    # tc.square计算输入的平方值
    # tc.mean沿 axis 计算 x 的平均值，默认axis是None，则对输入的全部元素计算平均值。
    error = tc.mean(tc.square(y_true - y_pred))/2

    return error


# 构造一个简单的样例进行测试:[N,1], N=2
y_true = tc.tensor([[-0.2], [4.9]], dtype=tc.float32)
y_pred = tc.tensor([[1.3], [2.5]], dtype=tc.float32)

error = mean_squared_error(y_true=y_true, y_pred=y_pred).item()
print("error:", error)

输出结果如下：
在这里插入图片描述

2.2.4 模型优化

在这里插入图片描述

def optimizer_lsm(model, X, y, reg_lambda=0):
    """
      输入：
         - model: 模型
         - X: tensor, 特征数据，shape=[N,D]
         - y: tensor,标签数据，shape=[N]
         - reg_lambda: float, 正则化系数，默认为0
      输出：
         - model: 优化好的模型
      """

    N, D = X.shape

    # 对输入特征数据所有特征向量求平均
    x_bar_tran = tc.mean(X, axis=0).T

    # 求标签的均值,shape=[1]
    y_bar = tc.mean(y)

    # tc.subtract通过广播的方式实现矩阵减向量
    x_sub = tc.subtract(X, x_bar_tran)

    # 使用tc.all判断输入tensor是否全0
    if tc.all(x_sub) == 0:
        model.params['b'] = y_bar
        model.params['w'] = tc.zeros(D)
        return model

    # tc.inverse求方阵的逆
    tmp = tc.inverse(tc.matmul(x_sub.T, x_sub) + reg_lambda * tc.eye(D))

    w = tc.matmul(tc.matmul(tmp, x_sub.T), (y - y_bar))

    b = y_bar - tc.matmul(x_bar_tran, w)

    model.params['b'] = b
    model.params['w'] = tc.squeeze(w, axis=-1)

    return model

2.2.5 模型训练

在准备了数据、模型、损失函数和参数学习的实现之后，我们开始模型的训练。在回归任务中，模型的评价指标和损失函数一致，都为均方误差。

通过上文实现的线性回归类来拟合训练数据，并输出模型在训练集上的损失。

input_size = 1
model = Linear(input_size)
model = optimizer_lsm(model,X_train.reshape([-1,1]),y_train.reshape([-1,1]))
print("w_pred:",model.params['w'].item(), "b_pred: ", model.params['b'].item())

y_train_pred = model(X_train.reshape([-1,1])).squeeze()
train_error = mean_squared_error(y_true=y_train, y_pred=y_train_pred).item()
print("train error: ",train_error)

运行结果如下：
在这里插入图片描述

model_large = Linear(input_size)
model_large = optimizer_lsm(model_large,X_train_large.reshape([-1,1]),y_train_large.reshape([-1,1]))
print("w_pred large:",model_large.params['w'].item(), "b_pred large: ", model_large.params['b'].item())

y_train_pred_large = model_large(X_train_large.reshape([-1,1])).squeeze()
train_error_large = mean_squared_error(y_true=y_train_large, y_pred=y_train_pred_large).item()
print("train error large: ",train_error_large)

运行结果如下：
在这里插入图片描述

2.2.6 模型评估

下面用训练好的模型预测一下测试集的标签，并计算在测试集上的损失。

y_test_pred = model(X_test.reshape([-1,1])).squeeze()
test_error = mean_squared_error(y_true=y_test, y_pred=y_test_pred).item()
print("test error: ",test_error)

运行结果如下：
在这里插入图片描述

y_test_pred_large = model_large(X_test.reshape([-1,1])).squeeze()
test_error_large = mean_squared_error(y_true=y_test, y_pred=y_test_pred_large).item()
print("test error large: ",test_error_large)

运行结果如下：
在这里插入图片描述

2.2.7 样本数量 & 正则化系数

（1）调整训练数据的样本数量，由 100 调整到 5000，观察对模型性能的影响。

（2）调整正则化系数，观察对模型性能的影响。

2.3 多项式回归

在这里插入图片描述
（借老师个图）

2.3.1 数据集构建

假设我们要拟合的非线性函数为一个缩放后的sin函数。

构建训练和测试数据，其中：
训练数样本 15 个，测试样本 10 个，高斯噪声标准差为 0.1，自变量范围为 (0,1)。

import torch as tc
from matplotlib import pyplot as plt  # matplotlib 是 Python 的绘图库
import math
from random import *

# sin函数: sin(2 * pi * x)
def sin(x):
    y = tc.sin(2 * math.pi * x, )
    return y

def create_toy_data(func, interval, sample_num, noise, add_outlier, outlier_ratio):
    """
    根据给定的函数，生成样本
    输入：
       - func：函数
       - interval： x的取值范围
       - sample_num： 样本数目
       - noise： 噪声均方差
       - add_outlier：是否生成异常值
       - outlier_ratio：异常值占比
    输出：
       - X: 特征数据，shape=[n_samples,1]
       - y: 标签数据，shape=[n_samples,1]
    """
    # 均匀采样
    # 使用torch.rand在生成sample_num个随机数，这些数范围在interval[0]到interval[1]之间，并组成一个张量
    X = tc.rand(sample_num) * (interval[1] - interval[0]) + interval[0]
    y = func(X)

    # 生成高斯分布的标签噪声
    # 使用torch.normal生成0均值，noise标准差的数据
    epsilon = tc.normal(0, noise, size=y.shape)
    y = y + epsilon

    if add_outlier:  # 生成额外的异常点
        outlier_num = int(len(y) * outlier_ratio)
        if outlier_num != 0:
            # 使用torch.randint生成服从均匀分布的、范围在[0, len(y))的随机Tensor
            outlier_idx = tc.randint(outlier_num, (0, len(y)))
            y[outlier_idx] = y[outlier_idx] * 5
    return X, y


# 生成数据
func = sin
interval = (0, 1)
train_num = 15
test_num = 10
noise = 0.5  # 0.1
add_outlier = False
outlier_ratio = 0.001
X_train, y_train = create_toy_data(func, interval, train_num, noise, add_outlier, outlier_ratio)
X_test, y_test = create_toy_data(func, interval, test_num, noise, add_outlier, outlier_ratio)

X_underlying = tc.linspace(interval[0], interval[1], train_num)
y_underlying = sin(X_underlying)

# 绘制图像
plt.rcParams['figure.figsize'] = (8.0, 6.0)
plt.scatter(X_train, y_train, color='y', s=20, label="train data")
plt.scatter(X_test, y_test, color="r", s=20, label="test data")
plt.plot(X_underlying, y_underlying, c='#000000', label=r"$\sin(2\pi x)$")
plt.legend(fontsize='x-large')
plt.savefig('ml-vis2.pdf')
plt.show()

运行结果如下：
在这里插入图片描述
我们可以看到训练集数据点在sin图像上下波动，这正是我们想要的结果。

这里这条sin（2Πx）函数并不是很光滑，因为步长太大的缘故，我这里取得是训练集中样点的个数来分段，可以将X_underlying的取值细化，这样得出来的图像就平滑许多了。

2.3.2 模型构建

通过多项式的定义可以看出，多项式回归和线性回归一样，同样学习参数w，只不过需要对输入特征ϕ(x)根据多项式阶数进行变换。因此，我们可以套用求解线性回归参数的方法来求解多项式回归参数。
首先，我们实现 多项式基函数 polynomial_basis_function对原始特征x进行转换。

# 多项式转换
def polynomial_basis_function(x, degree):
    """
    输入：
       - x: tensor, 输入的数据，shape=[N,1]
       - degree: int, 多项式的阶数
       example Input: [[2], [3], [4]], degree=2
       example Output: [[2^1, 2^2], [3^1, 3^2], [4^1, 4^2]]
       注意：本案例中,在degree>=1时不生成全为1的一列数据；degree为0时生成形状与输入相同，全1的Tensor
    输出：
       - x_result： tensor
    """

    if degree == 0:
        return tc.ones(x.shape, dtype=tc.float32)

    x_tmp = x
    x_result = x_tmp

    for i in range(2, degree + 1):
        x_tmp = tc.multiply(x_tmp, x)  # 逐元素相乘
        x_result = tc.concat((x_result, x_tmp), dim=-1)

    return x_result


# 简单测试
data = [[2], [3], [4]]
X = tc.tensor(data=data, dtype=tc.float32)
degree = 3
transformed_X = polynomial_basis_function(X, degree)
print("转换前：", X)
print("阶数为", degree, "转换后：", transformed_X)

输出结果如下：
在这里插入图片描述

2.3.3 模型训练

对于多项式回归，我们可以同样使用前面线性回归中定义的LinearRegression算子、训练函数train、均方误差函数mean_squared_error。拟合训练数据的目标是最小化损失函数，同线性回归一样，也可以通过矩阵运算直接求出 w 的值。

我们设定不同的多项式阶，M的取值分别为0、1、3、8，之前构造的训练集上进行训练，观察样本数据对sin⁡曲线的拟合结果。

plt.rcParams['figure.figsize'] = (12.0, 8.0)

for i, degree in enumerate([0, 1, 3, 8]):  # []中为多项式的阶数
    model = Linear(degree)
    X_train_transformed = polynomial_basis_function(X_train.reshape([-1, 1]), degree)
    X_underlying_transformed = polynomial_basis_function(X_underlying.reshape([-1, 1]), degree)

    model = optimizer_lsm(model, X_train_transformed, y_train.reshape([-1, 1]))  # 拟合得到参数

    y_underlying_pred = model(X_underlying_transformed).squeeze()

    print(model.params)

    # 绘制图像
    plt.subplot(2, 2, i + 1)
    plt.scatter(X_train, y_train, facecolor="none", edgecolor='#e4007f', s=50, label="train data")
    plt.plot(X_underlying, y_underlying, c='#000000', label=r"$\sin(2\pi x)$")
    plt.plot(X_underlying, y_underlying_pred, c='#f19ec2', label="predicted function")
    plt.ylim(-2, 1.5)
    plt.annotate("M={}".format(degree), xy=(0.95, -1.4))

# plt.legend(bbox_to_anchor=(1.05, 0.64), loc=2, borderaxespad=0.)
plt.legend(loc='lower left', fontsize='x-large')
plt.savefig('ml-vis3.pdf')
plt.show()

输出结果如下：
在这里插入图片描述

2.3.4 模型评估

通过均方误差来衡量训练误差、测试误差以及在没有噪音的加入下sin函数值与多项式回归值之间的误差，更加真实地反映拟合结果。多项式分布阶数从0到8进行遍历。

对于模型过拟合的情况，可以引入正则化方法，通过向误差函数中添加一个惩罚项来避免系数倾向于较大的取值。

下面通过均方误差来衡量训练误差、测试误差以及在没有噪音的加入下sin函数值与多项式回归值之间的误差，更加真实地反映拟合结果。多项式分布阶数从0到8进行遍历。

# 训练误差和测试误差
training_errors = []
test_errors = []
distribution_errors = []

# 遍历多项式阶数
for i in range(9):
    model = Linear(i)

    X_train_transformed = polynomial_basis_function(X_train.reshape([-1, 1]), i)
    X_test_transformed = polynomial_basis_function(X_test.reshape([-1, 1]), i)
    X_underlying_transformed = polynomial_basis_function(X_underlying.reshape([-1, 1]), i)

    optimizer_lsm(model, X_train_transformed, y_train.reshape([-1, 1]))

    y_train_pred = model(X_train_transformed).squeeze()
    y_test_pred = model(X_test_transformed).squeeze()
    y_underlying_pred = model(X_underlying_transformed).squeeze()

    train_mse = mean_squared_error(y_train, y_train_pred).item()
    training_errors.append(train_mse)

    test_mse = mean_squared_error(y_test, y_test_pred).item()
    test_errors.append(test_mse)

    # distribution_mse = mean_squared_error(y_true=y_underlying, y_pred=y_underlying_pred).item()
    # distribution_errors.append(distribution_mse)

print("train errors: \n", training_errors)
print("test errors: \n", test_errors)
# print ("distribution errors: \n", distribution_errors)

# 绘制图片
plt.rcParams['figure.figsize'] = (8.0, 6.0)
plt.plot(training_errors, '-.', mfc="none", mec='#e4007f', ms=10, c='#e4007f', label="Training")
plt.plot(test_errors, '--', mfc="none", mec='#f19ec2', ms=10, c='#f19ec2', label="Test")
# plt.plot(distribution_errors, '-', mfc="none", mec="#3D3D3F", ms=10, c="#3D3D3F", label="Distribution")
plt.legend(fontsize='x-large')
plt.xlabel("degree")
plt.ylabel("MSE")
plt.savefig('ml-mse-error.pdf')
plt.show()

运行结果如下：
在这里插入图片描述
观察可视化结果：

当阶数较低的时候，模型的表示能力有限，训练误差和测试误差都很高，代表模型欠拟合；
当阶数较高的时候，模型表示能力强，但将训练数据中的噪声也作为特征进行学习，一般情况下训练误差继续降低而测试误差显著升高，代表模型过拟合。

对于模型过拟合的情况，可以引入正则化方法，通过向误差函数中添加一个惩罚项来避免系数倾向于较大的取值。下面加入l2正则化项，查看拟合结果。

degree = 8 # 多项式阶数
reg_lambda = 0.0001 # 正则化系数

X_train_transformed = polynomial_basis_function(X_train.reshape([-1,1]), degree)
X_test_transformed = polynomial_basis_function(X_test.reshape([-1,1]), degree)
X_underlying_transformed = polynomial_basis_function(X_underlying.reshape([-1,1]), degree)

model = Linear(degree)

optimizer_lsm(model,X_train_transformed,y_train.reshape([-1,1]))

y_test_pred=model(X_test_transformed).squeeze()
y_underlying_pred=model(X_underlying_transformed).squeeze()

model_reg = Linear(degree)

optimizer_lsm(model_reg,X_train_transformed,y_train.reshape([-1,1]),reg_lambda=reg_lambda)

y_test_pred_reg=model_reg(X_test_transformed).squeeze()
y_underlying_pred_reg=model_reg(X_underlying_transformed).squeeze()

mse = mean_squared_error(y_test, y_test_pred).item()
print("mse:",mse)
mes_reg = mean_squared_error(y_test, y_test_pred_reg).item()
print("mse_with_l2_reg:",mes_reg)

# 绘制图像
plt.scatter(X_train, y_train, facecolor="none", edgecolor="#e4007f", s=50, label="train data")
plt.plot(X_underlying, y_underlying, c='#000000', label=r"$\sin(2\pi x)$")
plt.plot(X_underlying, y_underlying_pred, c='#e4007f', linestyle="--", label="$deg. = 8$")
plt.plot(X_underlying, y_underlying_pred_reg, c='#f19ec2', linestyle="-.", label="$deg. = 8, \ell_2 reg$")
plt.ylim(-1.5, 1.5)
plt.annotate("lambda={}".format(reg_lambda), xy=(0.82, -1.4))
plt.legend(fontsize='large')
plt.savefig('ml-vis4.pdf')
plt.show()

运行结果如下：
加粗样式 在这里插入图片描述
观察可视化结果，其中浅红色曲线为加入l2正则后多项式分布拟合结果，红色曲线为未加入l2正则的拟合结果，浅红色曲线的拟合效果明显好于红色曲线。

2.4 Runner类介绍

机器学习方法流程包括数据集构建、模型构建、损失函数定义、优化器、模型训练、模型评价、模型预测等环节。

为了更方便地将上述环节规范化，我们将机器学习模型的基本要素封装成一个Runner类。

除上述提到的要素外，再加上模型保存、模型加载等功能。

Runner类的成员函数定义如下：

__init__函数：实例化Runner类，需要传入模型、损失函数、优化器和评价指标等；
train函数：模型训练，指定模型训练需要的训练集和验证集；
evaluate函数：通过对训练好的模型进行评价，在验证集或测试集上查看模型训练效果；
predict函数：选取一条数据对训练好的模型进行预测；
save_model函数：模型在训练过程和训练结束后需要进行保存；
load_model函数：调用加载之前保存的模型。

（借老师一个图）

class Runner():
    def __init__(self,model,model_loss,model_optimizer,model_evaluate):
        self.model = model
        self.model_loss = model_loss
        self.model_optimizer = model_optimizer
        self.model_evaluate = model_evaluate
    def train(self,Train_feature,Test_feature,**kwargs):
        self.Train_feature = Train_feature
        pass
        return #各个参数
    def evaluate(self,label,**kwargs):
        loss = self.model_loss(label,self.predict(self.Train_feature))
        # accuracy =
        pass
        return #返回精度和误差
    def predict(self,X,**kwargs):
        #用train 的各个参数来对输入序列X进行预测，返回预测值
        pass
        return #返回预测值
    def save_model(self,save_path):
        pass
    def load_model(self,model_path):
        pass

2.5 基于线性回归的波士顿房价预测

使用线性回归来对马萨诸塞州波士顿郊区的房屋进行预测。

实验流程主要包含如下5个步骤：

数据处理：包括数据清洗（缺失值和异常值处理）、数据集划分，以便数据可以被模型正常读取，并具有良好的泛化性;
模型构建：定义线性回归模型类；
训练配置：训练相关的一些配置，如：优化算法、评价指标等；
组装训练框架Runner：Runner用于管理模型训练和测试过程；
模型训练和测试：利用Runner进行模型训练和测试。

2.5.1 数据处理

2.5.1.1 数据清洗

仿佛是第一次听说这个词。。搜索后才有些眉目。这里附上该博主的链接数据清洗

2.5.1.1.1 重复观测处理

import pandas as pd # 开源数据分析和操作工具


# 利用pandas加载波士顿房价的数据集
data = pd.read_csv("boston_house_prices.csv")
print("是否存在重复值：", any(data.duplicated()))
data.drop_duplicates(inplace=True)
data.to_csv('boston_house_prices.csv', index=False)   # 保存移除重复值后的数据集

运行结果如下：
在这里插入图片描述

2.5.1.1.2 缺失值处理

print("查看各个属性的缺失值：\n", data.isna().sum())

输出结果如下：
在这里插入图片描述
由此可以看到，该数据库无缺失值。

2.5.1.1.3 异常值处理

异常值是指那些远离正常值的观测值，异常值的出现会给模型的常见和预测产生严重的后果，但有时也会利用异常值进行异常数据查找。
对于异常值的检测，一般采用两种方法，一种是标准差法，另一种是箱线图判别法。标准差法的判别公式是outlier > x+nδ 或者outlier < x-nδ,其中x为样本均值，δ为样本标准差。当n=2时，满足条件的观测就是异常值；当n=3时，满足条件的观测就是极端异常值。箱线图的判别公式是outlier > Q3+nIQR 或者outlier < Q1-nIQR,其中Q1为下四分位数，Q3为上四分位数，IQR为上四分位数和下四分位数的差，当n=1.5时，满足条件的观测为异常值，当n=3时，满足条件的观测为极端异常值。

箱线图介绍：
在这里插入图片描述

# 箱线图查看异常值分布
def boxplot(data, fig_name):
    # 绘制每个属性的箱线图
    data_col = list(data.columns)

    # 连续画几个图片
    plt.figure(figsize=(5, 5), dpi=300)
    # 子图调整
    plt.subplots_adjust(wspace=0.6)
    # 每个特征画一个箱线图
    for i, col_name in enumerate(data_col):
        plt.subplot(3, 5, i + 1)
        # 画箱线图
        plt.boxplot(data[col_name],
                    showmeans=True,
                    meanprops={"markersize": 1, "marker": "D", "markeredgecolor": '#f19ec2'},  # 均值的属性
                    medianprops={"color": '#e4007f'},  # 中位数线的属性
                    whiskerprops={"color": '#e4007f', "linewidth": 0.4, 'linestyle': "--"},
                    flierprops={"markersize": 0.4},
                    )
        # 图名
        plt.title(col_name, fontdict={"size": 5}, pad=2)
        # y方向刻度
        plt.yticks(fontsize=4, rotation=90)
        plt.tick_params(pad=0.5)
        # x方向刻度
        plt.xticks([])
    plt.savefig(fig_name)
    plt.show()


boxplot(data, 'ml-vis5.pdf')

在这里插入图片描述
由图可以直观的看出，数据中存在很多的异常值，即超出上下两条杠的部分中的黑圆圈，我们将这些异常值认为是数据集中的“噪声”，并将临界值取代噪声点，从而完成对数据集异常值的处理。

用于替换的代码：

# 四分位处理异常值
num_features = data.select_dtypes(exclude=['object', 'bool']).columns.tolist()

for feature in num_features:
    if feature == 'CHAS':
        continue

    Q1 = data[feature].quantile(q=0.25)  # 下四分位
    Q3 = data[feature].quantile(q=0.75)  # 上四分位

    IQR = Q3 - Q1
    top = Q3 + 1.5 * IQR  # 最大估计值
    bot = Q1 - 1.5 * IQR  # 最小估计值
    values = data[feature].values
    values[values > top] = top  # 临界值取代噪声
    values[values < bot] = bot  # 临界值取代噪声
    data[feature] = values.astype(data[feature].dtypes)

# 再次查看箱线图，异常值已被临界值替换（数据量较多或本身异常值较少时，箱线图展示会不容易体现出来）
boxplot(data, 'ml-vis6.pdf')

运行结果如下：
在这里插入图片描述

2.5.1.2 数据集划分

数据集划分为测试集和训练集：

import torch

torch.manual_seed(10)


# 划分训练集和测试集
def train_test_split(X, y, train_percent=0.8):
    n = len(X)
    shuffled_indices = torch.randperm(n)  # 返回一个数值在0到n-1、随机排列的1-D Tensor
    train_set_size = int(n * train_percent)
    train_indices = shuffled_indices[:train_set_size]
    test_indices = shuffled_indices[train_set_size:]

    X = X.values
    y = y.values

    X_train = X[train_indices]
    y_train = y[train_indices]

    X_test = X[test_indices]
    y_test = y[test_indices]

    return X_train, X_test, y_train, y_test


X = data.drop(['MEDV'], axis=1)
y = data['MEDV']

X_train, X_test, y_train, y_test = train_test_split(X, y)  # X_train每一行是个样本，shape[N,D]

2.5.1.3 特征工程

为了消除纲量对数据特征之间影响，在模型训练前，需要对特征数据进行归一化处理，将数据缩放到[0, 1]区间内，使得不同特征之间具有可比性。

import torch
X_train = torch.as_tensor(X_train,dtype=torch.float32)
X_test = torch.as_tensor(X_test,dtype=torch.float32)
y_train = torch.as_tensor(y_train,dtype=torch.float32)
y_test = torch.as_tensor(y_test,dtype=torch.float32)
X_min = torch.min(X_train,axis=0)[0]
X_max = torch.max(X_train,axis=0)[0]
X_train = (X_train-X_min)/(X_max-X_min)
X_test  = (X_test-X_min)/(X_max-X_min)
# 训练集构造
train_dataset=(X_train,y_train)
# 测试集构造
test_dataset=(X_test,y_test)

2.5.2 模型构建

from nndl.op import Linear

# 模型实例化
input_size = 12
model=Linear(input_size)

其中nndl.op代码如下：

import torch
torch.seed()  # 设置随机种子

class Op(object):
    def __init__(self):
        pass

    def __call__(self, inputs):
        return self.forward(inputs)

    def forward(self, inputs):
        raise NotImplementedError

    def backward(self, inputs):
        raise NotImplementedError

# 线性算子
class Linear(Op):
    def __init__(self, dimension):
        """
        输入：
           - dimension:模型要处理的数据特征向量长度
        """

        self.dim = dimension

        # 模型参数
        self.params = {}
        self.params['w'] = torch.randn(self.dim, 1, dtype=torch.float32)
        self.params['b'] = torch.zeros([1], dtype=torch.float32)

    def __call__(self, X):
        return self.forward(X)

    # 前向函数
    def forward(self, X):
        """
        输入：
           - X: tensor, shape=[N,D]
           注意这里的X矩阵是由N个x向量的转置拼接成的，与原教材行向量表示方式不一致
        输出：
           - y_pred： tensor, shape=[N]
        """

        N, D = X.shape

        if self.dim == 0:
            return torch.full([N, 1], fill_value=self.params['b'])
        
        assert D == self.dim  # 输入数据维度合法性验证

        # 使用torch.matmul计算两个tensor的乘积
        y_pred = torch.matmul(X,self.params['w'])+self.params['b']
        
        return y_pred

2.5.3 完善Runner类

在测试集上使用MSE对模型性能进行评估。

import torch.nn as nn
mse_loss = nn.MSELoss()

完整实现如下：

import torch.nn as nn
import torch
import os
from nndl.op import Linear
from nndl.opitimizer import optimizer_lsm
# 模型实例化
input_size = 12
model=Linear(input_size)
mse_loss = nn.MSELoss()



class Runner(object):
    def __init__(self, model, optimizer, loss_fn, metric):
        # 优化器和损失函数为None,不再关注

        # 模型
        self.model = model
        # 评估指标
        self.metric = metric
        # 优化器
        self.optimizer = optimizer

    def train(self, dataset, reg_lambda, model_dir):
        X, y = dataset
        self.optimizer(self.model, X, y, reg_lambda)

        # 保存模型
        self.save_model(model_dir)

    def evaluate(self, dataset, **kwargs):
        X, y = dataset

        y_pred = self.model(X)
        result = self.metric(y_pred, y)

        return result

    def predict(self, X, **kwargs):
        return self.model(X)

    def save_model(self, model_dir):
        if not os.path.exists(model_dir):
            os.makedirs(model_dir)

        params_saved_path = os.path.join(model_dir, 'params.pdtensor')
        torch.save(model.params, params_saved_path)

    def load_model(self, model_dir):
        params_saved_path = os.path.join(model_dir, 'params.pdtensor')
        self.model.params = torch.load(params_saved_path)


optimizer = optimizer_lsm
runner = Runner(model, optimizer=optimizer,loss_fn=None, metric=mse_loss)

2.5.4 模型训练

optimizer = optimizer_lsm
runner = Runner(model, optimizer=optimizer,loss_fn=None, metric=mse_loss)

# 模型保存文件夹
saved_dir = 'D:/models_'

# 启动训练
runner.train(train_dataset,reg_lambda=0,model_dir=saved_dir)
columns_list = data.columns.to_list()
weights = runner.model.params['w'].tolist()
b = runner.model.params['b'].item()

for i in range(len(weights)):
    print(columns_list[i],"weight:",weights[i])

print("b:",b)

运行结果如下：
在这里插入图片描述

2.5.5 模型测试

# 加载模型权重
runner.load_model(saved_dir)
mse = runner.evaluate(test_dataset)
print('MSE:', mse.item())

运行结果如下：
在这里插入图片描述

2.5.6 模型预测

runner.load_model(saved_dir)
pred = runner.predict(X_test[:1])
print("真实房价：",y_test[:1].item())
print("预测的房价：",pred.item())

运行结果如下：
在这里插入图片描述
这里附上optimizer文件

import torch

def optimizer_lsm(model, X, y, reg_lambda=0):
  """
    输入：
       - model: 模型
       - X: tensor, 特征数据，shape=[N,D]
       - y: tensor,标签数据，shape=[N]
       - reg_lambda: float, 正则化系数，默认为0
    输出：
       - model: 优化好的模型
    """

  N, D = X.shape

  # 对输入特征数据所有特征向量求平均
  x_bar_tran = torch.mean(X,axis=0).T 
  
  # 求标签的均值,shape=[1]
  y_bar = torch.mean(y)
  
  # torch.subtract通过广播的方式实现矩阵减向量
  x_sub = torch.subtract(X,x_bar_tran)

  # 使用torch.all判断输入tensor是否全0
  if torch.all(x_sub == 0):
    model.params['b'] = y_bar
    model.params['w'] = torch.zeros(shape=[D])
    return model
  
  # torch.inverse求方阵的逆
  tmp = torch.inverse(torch.matmul(x_sub.T, x_sub) + reg_lambda*torch.eye(D))

  w = torch.matmul(torch.matmul(tmp,x_sub.T),(y-y_bar))
  
  b = y_bar-torch.matmul(x_bar_tran,w)
  
  model.params['b'] = b
  model.params['w'] = torch.squeeze(w,axis=-1)

  return model

from abc import abstractmethod

#新增优化器基类
class Optimizer(object):
    def __init__(self, init_lr, model):
        """
        优化器类初始化
        """
        #初始化学习率，用于参数更新的计算
        self.init_lr = init_lr
        #指定优化器需要优化的模型
        self.model = model

    @abstractmethod
    def step(self):
        """
        定义每次迭代如何更新参数
        """
        pass

#新增梯度下降法优化器
class SimpleBatchGD(Optimizer):
    def __init__(self, init_lr, model):
        super(SimpleBatchGD, self).__init__(init_lr=init_lr, model=model)

    def step(self):
        #参数更新
        #遍历所有参数，按照公式(3.8)和(3.9)更新参数
        if isinstance(self.model.params, dict):
            for key in self.model.params.keys():
                self.model.params[key] = self.model.params[key] - self.init_lr * self.model.grads[key]

问题1：使用类实现机器学习模型的基本要素有什么优点？
答：使用类实现机器学习模型基本要素方便整洁，更容易被人理解。让整体思路更加清晰，优化了许多后续的代码。同时也加强了代码的稳定性和正确性。
问题2：算子op、优化器opitimizer放在单独的文件中，主程序在使用时调用该文件。这样做有什么优点？
答：简化操作，方便备份，以后再用到算子op,优化器opitimizer时不用再写一遍，直接调用该文件就可以了。同时使原代码更简洁，减少了很多空间和时间。
问题3：线性回归通常使用平方损失函数，能否使用交叉熵损失函数？为什么？
答：交叉熵损失函数主要用于离散情况，假设误差是二值分部，平方损失函数是实际值与目标值之间的误差，假设误差服从正态分布。因此线性回归通常使用平方损失函数。

NNDL 实验三
 机器学习实现简单的线性回归方程

红肚兜

关注

4
点赞
踩
1

收藏

觉得还不错? 一键收藏
打赏
0
评论
NNDL 实验三线性回归

机器学习方法流程包括数据集构建、模型构建、损失函数定义、优化器、模型训练、模型评价、模型预测等环节。为了更方便地将上述环节规范化，我们将机器学习模型的基本要素封装成一个Runner类。除上述提到的要素外，再加上模型保存、模型加载等功能。__init__函数：实例化Runner类，需要传入模型、损失函数、优化器和评价指标等；train函数：模型训练，指定模型训练需要的训练集和验证集；evaluate函数：通过对训练好的模型进行评价，在验证集或测试集上查看模型训练效果；predict。
复制链接

扫一扫