[机器学习西瓜书]3.线性回归

介于看完之后老是不知道看了些啥,所以开这么一部分blog对看到的内容进行记录整理
当然知识介绍肯定没有书上写的详细,仅限于自己对内容的整理

0. 线性模型概述

  1. 试图学习一个通过属性的线性组合来进行预测的函数
  2. 这里学到的主要内容
    线性回归——最小二乘法
    局部加权线性回归
    对数几率回归(逻辑斯蒂回归)
    感知机(书上是神经网络部分的内容,但是讲的时候把他放到这部分讲啦)

1. 线性回归

总体内容

线性回归:学习一个线性模型尽可能准确地预测实值输出标记
提到线性回归就需要涉及另外两个概念:均方误差和最小二乘法

  • 均方误差: ( w ∗ , b ∗ ) = a r g w , b m i n ∑ ( y i − w x i − b ) 2 (w^*,b^*)=arg_{w,b}min\sum(y_i-wx_i-b)^2 (w,b)=argw,bmin(yiwxib)2
    均方误差是回归任务中最常用的性能度量
  • 最小二乘法:
    基于均方误差来进行模型求解的方法
    最小二乘法试图找到一条直线,使所有样本到直线上的欧氏距离之和最小

线性回归模型的最小二乘“参数估计”是求解w和b使下面公式最小:
E w , b = ∑ i = 1 m ( y i − w x i − b ) 2 E_{w,b}=\sum _{i=1}^m(y_i-wx_i-b)^2 Ew,b=i=1m(yiwxib)2
分别对w和b求偏导并令偏导=0,可以得到w和b的最优闭式解:
w = ∑ i = 1 m y i ( x i − x ‾ ) ∑ i = 1 m x i 2 − 1 m ( ∑ i = 1 m x i ) 2 = ∑ i = 1 m ( y i x i − y i x ‾ ) ∑ i = 1 m ( x i 2 − x i x ‾ ) b = 1 m ∑ i = 1 m ( y i − w x i ) w=\frac{\sum_{i=1}^m y_i(x_i-\overline{x})}{\sum_{i=1}^m x_i^2-\frac{1}{m}(\sum_{i=1}^mx_i)^2} \\=\frac{\sum_{i=1}^m (y_ix_i-y_i\overline{x})}{\sum_{i=1}^m(x_i^2-x_i\overline{x})} \\ b=\frac{1}{m}\sum_{i=1}^m(y_i-wx_i) w=i=1mxi2m1(i=1mxi)2i=1myi(xix)=i=1m(xi2xix)i=1m(yixiyix)b=m1i=1m(yiwxi)

使用最小二乘法进行参数估计

f ( x i ) = w T x i + b 中,可以为 x 扩展一列常数 1 变为增广矩阵, b 可以视为乘上这列常数 1 f(x_i)=w^Tx_i+b中,可以为x扩展一列常数1变为增广矩阵,b可以视为乘上这列常数1 f(xi)=wTxi+b中,可以为x扩展一列常数1变为增广矩阵,b可以视为乘上这列常数1
常将w和b吸收入向量形式 w ‾ = ( w ; b ) \overline{w}=(w;b) w=(w;b)
将均方误差写作向量的形式,即有:
w ^ ∗ = a r g   m i n w ^ ( y − X w ^ ) T ( y − X w ^ ) \hat{w}^* = arg\ min_{\hat{w}}(y-X\hat{w})^T(y-X\hat{w}) w^=arg minw^(yXw^)T(yXw^)
E ^ x = ( y − X w ^ ) T ( y − X w ^ ) E_{\hat{}x}=(y-X\hat{w})^T(y-X\hat{w}) E^x=(yXw^)T(yXw^)
对w求导且XTX为满秩矩阵时令导数=0可得:
w ^ ∗ = ( X T X ) − 1 X T y \hat{w}^*=(X^TX)^{-1}X^Ty w^=(XTX)1XTy

令xi=(xi;1)最终学得的多元线性回归模型为:
f ( x i ^ ) = x i ^ T ( X T X ) − 1 X T y f(\hat{x_i})=\hat{x_i}^T(X^TX)^{-1}X^Ty f(xi^)=xi^T(XTX)1XTy

现实中XTX往往不满秩,此时可以解出多个w,他们都能使均方误差最小化,选择哪一个解作为输出将由学习算法的归纳偏好决定

实验代码


import matplotlib.pyplot as plt
import numpy as np
from sklearn.linear_model import LinearRegression

"""
函数说明:加载数据

Parameters:
    filename - 文件名
    
Returns:
    xArr - x数据集
    yArr - y数据集

"""
def loadDataSet(filename):
    
    
    xArr = [] # 输入数据
    yArr = []
    fr = open(filename)
    # 计算特征个数,由于最后一列为y值所以减1
    numFeat = len(fr.readline().split('\t')) - 1

    for line in fr.readlines():
        lineArr = [] # 每一行的输入数据
        curLine = line.strip().split('\t')
        for i in range(numFeat):
            lineArr.append(float(curLine[i]))
        xArr.append(lineArr)
        yArr.append(float(curLine[-1]))
    return xArr, yArr


"""
函数说明:计算回归系数w

Parameters:
    xArr - x数据集
    yArr - y数据集
    
Returns:
    ws - 回归系数

"""
def standRegres(xArr, yArr):
    xMat = np.mat(xArr)

    #学习xMat生成方式,生成yMat。  yMat为yArr转成矩阵后再转置的结果。
    yMat = np.mat(yArr).T
    #计算 xTx
    xTx = xMat.T * xMat
    
    # 求矩阵的行列式
    if np.linalg.det(xTx) == 0.0: # linalg: 线性代数函数库
        print("矩阵为奇异矩阵,不能求逆")
        return
    #求系数矩阵ws   提示:     .I  :求逆矩阵
    ws = xTx.I * (xMat.T * yMat)
    return ws


"""
函数说明:绘制数据集以及回归直线

"""
def plotDataSet():
    xArr, yArr = loadDataSet('./实验3线性回归/ex0.txt')
    ws = standRegres(xArr, yArr)
    xMat = np.mat(xArr)
    yMat = np.mat(yArr)
    xCopy = xMat.copy()
    # 排序
    xCopy.sort(0) # 0表示按列排序,1表示按行排序
    yHat = xCopy * ws     
    
    plt.plot(xCopy[:,1],yHat,color='red')

    # 矩阵.A(等效于矩阵.getA())变成了数组
    plt.scatter(xMat[:, 1].flatten().A[0], yMat.flatten().A[0], s=20, c='blue', alpha=.5)
    plt.title('DataSet')
    plt.xlabel('X')
    plt.show()    
    
if __name__ == '__main__':
    plotDataSet()
    

局部加权线性回归

线性回归的一个问题是有可能出现欠拟合现象,因为它求的是具有最小均方误差的无偏估计。然鹅如果模型欠拟合将不能取得最好的预测效果。所以有些方法允许在估计中引人一些偏差,从而降低预测的均方误差。其中的一个方法是局部加权线性回归LWLR
在LWLR中,给待测点附近的每个店赋予一定的权重,于是公式就变成了:
f ( x i ^ ) = x i ^ T ( X T W X ) − 1 X T W y f(\hat{x_i})=\hat{x_i}^T(X^TWX)^{-1}X^TWy f(xi^)=xi^T(XTWX)1XTWy
LWLR使用“核”来对附近的点赋予更高的权重,核的类型可以自由选择,最常用的就是高斯核
高斯核对应权重:
w ( i , i ) = e x p ( ∣ x ( i ) − x ∣ − 2 k 2 ) w(i,i)=exp(\frac{|x^{(i)}-x|}{-2k^2}) w(i,i)=exp(2k2x(i)x)

实验代码

# -*- coding: utf-8 -*-

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

"""
函数说明:加载数据

Parameters:
    filename - 文件名
    
Returns:
    xArr - x数据集
    yArr - y数据集


"""
def loadDataSet(filename):
    # 计算特征个数,由于最后一列为y值所以减一
    numFeat = len(open(filename).readline().split('\t')) - 1
    xArr = []
    yArr = []
    fr = open(filename)
    for line in fr.readlines():
        lineArr = []
        curLine = line.strip().split('\t')
        for i in range(numFeat):
            lineArr.append(float(curLine[i]))
        xArr.append(lineArr)
        yArr.append(float(curLine[-1]))
    return xArr, yArr


"""
函数说明:使用局部加权线性回归计算回归系数w

Parameters:
    testPoint - 测试样本点
    xArr - x数据集
    yArr - y数据集
    k - 高斯核的k,自定义参数
    
Returns:
    testPoint上的回归结果


"""
def lwlr(testPoint, xArr, yArr, k=1.0):
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T
    m = np.shape(xMat)[0]
    # 创建加权对角阵
    weights = np.mat(np.eye((m)))
    for j in range(m):
        # 高斯核
        diffMat = testPoint - xMat[j, :]
        #计算 w[j,j] 
        weights[j, j] = np.exp(diffMat * diffMat.T / (-2 * pow(k, 2)))
    
    #计算 xTx(局部加权后)
    xTx = xMat.T * weights * xMat
    
    # 求矩阵的行列式
    if np.linalg.det(xTx) == 0.0:
        print("矩阵为奇异矩阵,不能求逆")
        return
    #求系数矩阵ws   提示:     .I  :求逆矩阵
    #
    ws = xTx.I * xMat.T * weights * yMat
    
    #返回测试样本点testPoint的预测结果
    return testPoint * ws


"""
函数说明:局部加权线性回归测试

Parameters:
    testArr - 测试数据集
    xArr - x数据集
    yArr - y数据集
    k - 高斯核的k,自定义参数
    
Returns:
    yHat - testArr上的回归结果


"""
def lwlrTest(testArr, xArr, yArr, k=1.0):
    m = np.shape(testArr)[0]
    yHat = np.zeros(m)
    for i in range(m):
        yHat[i] = lwlr(testArr[i], xArr, yArr, k)
    return yHat


"""
函数说明:绘制多条局部加权回归曲线

Parameters:
    None
    
Returns:
    None


"""
def plotlwlrRegression():
    # 设置中文字体
    font = FontProperties(fname=r"C:\Windows\Fonts\simsun.ttc", size=14)
    xArr, yArr = loadDataSet("./机器学习基础实验/实验4局部加权线性回归/ex1.txt")
    yHat_1 = lwlrTest(xArr, xArr, yArr, 1.0)
    yHat_2 = lwlrTest(xArr, xArr, yArr, 0.01)
    yHat_3 = lwlrTest(xArr, xArr, yArr, 0.003)
    xMat = np.mat(xArr)
    yMat = np.mat(yArr)
    
    yMat_1 = np.mat(yHat_1)
    yMat_2 = np.mat(yHat_2)
    yMat_3 = np.mat(yHat_3)
    
    
    
    
    ##########################   计算相关性 ############################
  
    print(np.corrcoef(yHat_1.T, yArr))
    print(np.corrcoef(yHat_2.T, yArr))
    print(np.corrcoef(yHat_3.T, yArr))
    
    # 排序,返回索引值
    srtInd = xMat[:, 1].argsort(0)
    xSort = xMat[srtInd][:, 0, :]
   
    fig, axs = plt.subplots(nrows=6, ncols=1, sharex=False, sharey=False, figsize=(10, 30))
    axs[0].plot(xSort[:, 1], yHat_1[srtInd], c='red')
    axs[1].plot(xSort[:, 1], yHat_2[srtInd], c='red')
    axs[2].plot(xSort[:, 1], yHat_3[srtInd], c='red')
    
    axs[0].scatter(xMat[:, 1].flatten().A[0], yMat.flatten().A[0], s=20, c='blue', alpha=.5)
    axs[1].scatter(xMat[:, 1].flatten().A[0], yMat.flatten().A[0], s=20, c='blue', alpha=.5)
    axs[2].scatter(xMat[:, 1].flatten().A[0], yMat.flatten().A[0], s=20, c='blue', alpha=.5)
    
    axs[3].scatter(xMat[:, 1].flatten().A[0], (yMat-yMat_1).flatten().A[0], s=20, c='blue', alpha=.5)
    axs[4].scatter(xMat[:, 1].flatten().A[0], (yMat-yMat_2).flatten().A[0], s=20, c='blue', alpha=.5)
    axs[5].scatter(xMat[:, 1].flatten().A[0], (yMat-yMat_3).flatten().A[0], s=20, c='blue', alpha=.5)
    
    axs0_title_text = axs[0].set_title(u'局部加权回归曲线,k=1.0', fontproperties=font)
    axs1_title_text = axs[1].set_title(u'局部加权回归曲线,k=0.01', fontproperties=font)
    axs2_title_text = axs[2].set_title(u'局部加权回归曲线,k=0.003', fontproperties=font)
    
    plt.setp(axs0_title_text, size=8, weight='bold', color='red')
    plt.setp(axs1_title_text, size=8, weight='bold', color='red')
    plt.setp(axs2_title_text, size=8, weight='bold', color='red')
    plt.xlabel('X')
    plt.show()
    
    
if __name__ == '__main__':
    plotlwlrRegression()

感知机perceptron

感知机perceptron是二分类的线性分类模型,输入是实例的特征向量,输出为实例的类别,取值一般为+1和-1二值,属于判别模型
西瓜书把这一部分放到了后面神经网络哪里,但是由于她也是线性模型,这里就一起学习了一下

实验代码

# -*- coding: utf-8 -*-


# 二维数据集决策边界可视化
from matplotlib.colors import ListedColormap
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np

"""
类说明:构建感知机

Parameters:
    eta - 学习率(0,1]
    n_iter - 迭代次数
    w_ - 训练后的权重数组
    errors_ - 每轮训练后的误差
    
Returns:
    None

"""
class Perceptron(object):
    def __init__(self, eta=0.01, n_iter=10):
        self.eta = eta
        self.n_iter = n_iter
    
    """
    Fit training data
    Parameters
    X: {array-like}, shape=[n_samples, n_features]
       Trainig vectors, where n_samples is the number of samples and n_features is the number of features.
    y: array-like, shape = [n_samples]
       Target values.
    Returns
    self: object
    """
    def fit(self, X, y):        # 拟合训练数据
        # self.w_中的权值初始化为一个零向量R(m+1),其中m是数据集中维度(特征)的数量
        # 我们在此基础上增加一个0权重列(也就是阈值)
        self.w_ = np.zeros(1 + X.shape[1])  # X.shape[1]:特征数,self.w_[0]是b
        self.errors_ = []
        for _ in range(self.n_iter):
            errors = 0
            for xi, target in zip(X, y):

                update = self.eta * (target - self.predict(xi))/2

                # 更新 w b
                self.w_[1:] += update * xi
                self.w_[0] += update

                # 每轮中错分类样本的数量
                errors += int(update != 0.0)    
            self.errors_.append(errors)
        return self
    
    def net_input(self, X):         # 激活函数
        """Calculate net input"""
        # 计算X和w_的点积 
        return np.dot(X, self.w_[1:]) + self.w_[0]
    
    def predict(self, X):
        """Return class label after unit step"""
        return np.where(self.net_input(X) >= 0.0, 1, -1)


"""
函数说明:导入数据集

Parameters:
    None
    
Returns:
    X - 特征矩阵
    y - label列向量

Modify:
    2018-08-28
"""
def DataSet():
    # 使用pandas库直接从UCI机器学习库中将鸢尾花数据集转换为DataFrame对象并加载到内存中
    df = pd.read_csv('https://archive.ics.uci.edu/ml/machine-learning-databases/iris/iris.data', header=None)
    # 使用tail方法显示数据最后5行以保证数据正确加载
    df.tail()
    print(df)
    # 提取前100个类标,50个山鸢尾类标,50个变色鸢尾类标
    # iloc works on the positions in the index (so it only takes integers).
    y = df.iloc[0:100, 4].values
    # -1代表山鸢尾 1代表变色鸢尾,将label存到y中
    # np.where用法相当于C语言的 ? : 
    y = np.where(y == 'Iris-setosa', -1, 1) # 三元运算符,第一个参数是if的条件,第二个参数是if为真时的值,第三个参数是if为假时的值
    '''当只有一个参数时,返回的是每个符合condition条件元素的坐标,返回的是以元组的形式'''
    # 提取特征0和特征1 
    X = df.iloc[0:100, [0, 2]].values   # 切片[0,100)行,每行取第0和2列->iloc切面是前闭后开

    print("************")
    print(X)
    print("=============")
    print(y)

    # 绘制散点图
    plt.scatter(X[:50, 0], X[:50, 1], color='red', marker='o', label='setosa')
    plt.scatter(X[50:100, 0], X[50:100, 1], color='blue', marker='x', label='versicolor')
    # 花瓣长度
    plt.xlabel('petal length')
    # 萼片长度
    plt.ylabel('sepal length')
    plt.legend(loc='upper left')
    plt.show()
    # X为特征矩阵,y为label列向量(-1,1)
    return X, y


"""
函数说明:绘制迭代次数与误分点个数之间的关系

Parameters:
    None
    
Returns:
    None

Modify:
    2018-08-28
"""
def NumOfErrors():
    # 导入数据
    X, y = DataSet()
    # 实例化感知机
    ppn = Perceptron(eta=0.1, n_iter=10)
    # 训练模型
    ppn.fit(X, y)
    # 绘制迭代次数与误分点个数之间的关系
    plt.plot(range(1, len(ppn.errors_) + 1), ppn.errors_, marker='o')
    plt.xlabel('Epochs')
    plt.ylabel('Number of misclassifications')
    plt.show()
    

"""
函数说明:绘制决策区域图像

Parameters:
    X - 特征矩阵
    y - label列向量
    classifier - 分类器
    resolution - 采样间隔为0.02
    
Returns:
    None

Modify:
    2018-08-28
"""
def plot_decision_regions(X, y, classifier, resolution=0.02):
    # 散点样式
    markers = ('s', 'x', 'o', '^', 'v')
    # 颜色元组
    colors = ('red', 'blue', 'lightgreen', 'gray', 'cyan')
    # np.unique该函数是去除数组中的重复数字,并进行排序之后输出。
    # ListedColormap主要用于生成非渐变的颜色映射
    cmap = ListedColormap(colors[:len(np.unique(y))])
    # 横轴范围
    x1_min, x1_max = X[:, 0].min() - 1, X[:, 0].max() + 1
    # 纵轴范围
    x2_min, x2_max = X[:, 1].min() - 1, X[:, 1].max() + 1
    # meshgrid函数将最大值、最小值向量生成二维数组xx1和xx2
    xx1, xx2 = np.meshgrid(np.arange(x1_min, x1_max, resolution), 
                           np.arange(x2_min, x2_max, resolution))
    z = classifier.predict(np.array([xx1.ravel(), xx2.ravel()]).T)
    z = z.reshape(xx1.shape)
    plt.contourf(xx1, xx2, z, alpha=0.4, cmap=cmap)
    plt.xlim(xx1.min(), xx1.max())
    plt.ylim(xx2.min(), xx2.max())
    # 对于一个可迭代的(iterable)/可遍历的对象(如列表、字符串),enumerate将其组成一个索引序列,利用它可以同时获得索引和值
    for idx, cl in enumerate(np.unique(y)):
        plt.scatter(x=X[y==cl, 0], y=X[y==cl, 1], alpha=0.8, c=cmap(idx), marker=markers[idx], label=cl)
        

if __name__ == '__main__':
    # 导入数据
    X, y = DataSet()

    '''
    X:两个特征,图中以x[0]作横坐标,x[1]作纵坐标进行分类
    y:标签,1和-1
    '''

    # 实例化感知机
    ppn = Perceptron(eta=0.1, n_iter=10)
    # 训练模型
    ppn.fit(X, y)
    plot_decision_regions(X, y, classifier=ppn)
    # 萼片长度
    plt.xlabel('sepal length [cm]')
    # 花瓣长度
    plt.ylabel('petal length [cm]')
    plt.legend(loc='upper left')
    plt.show()
    NumOfErrors()
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值