笔触狂放9

这个屌丝很懒，什么也没留下！

热门标签

多元线性回归，岭回归，lasso回归（具体代码（包括调用库代码和手写代码实现）+一点点心得)_多元线性回归代码

作者：笔触狂放9 | 2024-05-07 22:29:25

踩

多元线性回归代码

最近数据挖掘导论老师布置了一项作业，主要就是线性回归的实现，笔者之前听过吴恩达的线性回归的网课，但一直没有进行代码的实现，这次正好相对系统的整理一下，方便各位同学的学习，也希望能够对其进行优化，优化的点最后再说。笔者写这篇博客也为了给实验报告打底稿，各位小伙伴2021年9月30号提交报告的时候别跟我实验报告一样啊，打回的话苦的是自己人，到时候我直接一波举报，哈哈哈。不过，发表这篇文章笔者是真的希望给没有实现的同学提供一点小小的帮助或者提供一点小小的想法。

一.实验指导书：

数据挖掘实验指导书

数据预处理与线性回归

实验目的

掌握数据预处理方法
掌握线性回归预测基本原理与实现。

实验问题背景

鲍鱼的年龄可以通过鲍鱼壳的“环数”来判断，但是获取这个“环数”是十分耗时的，需要锯开壳，然后在显微镜下观察得到。

可以通过鲍鱼的其他特征比如性别、长度、直径、高度、整体重量、去壳后重量、脏器重量、壳的重量等，通过机器学习的方法来预测其环数，从而得到年龄，具有很大的应用价值。

实验问题描述

现有一份鲍鱼的数据集abalone.csv，该数据集有4177个数据样本，每个样本有9个特征，具体信息如下：

属性	数据类型	单位	内容描述
性别（Sex）	标称		M, F,I(infant)
长度（Length）	连续	毫米
直径（Diameter）	连续	毫米
高度（Height）	连续	毫米
整体重量（Whole weight）	连续	克
去壳后重量（Shucked weight）	连续	克
脏器重量（Viscera weight）	连续	克
壳的重量（Shell weight）	连续	克
环数（Rings）	连续

实验要求

通过数据集的前面8个特征预测环数，具体包括：

对数据做适当必要的预处理，包括（1）标称数据的连续化（2）未定义数据NaN的处理（3）非常大或者非常小的，明显错误的数据处理。
将数据集按照3:1分为训练集和验证集，在验证集上计算预测值和真值误差，作为评价模型的标准。
编程实现线性回归和岭回归算法，调用第三方软件包（比如Sklearn）实现线性回归、岭回归、LASSO回归，并将5种方法的结果比较，分析其优缺点。
分析算法实现中关键参数的影响，以曲线等可视化的方式展示。

二.要求了解pandas的dataframe等函数，numpy的array,ravel,mat等函数，sklean.linear_model库（主要就调用这个库函数实现功能），matplotlib.pyplot(画折线图）。

三.流程：

1.数据预处理，然后将数据集按照3:1分为训练集和验证集，在验证集上计算预测值和真值误差，作为评价模型的标准。


import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
import seaborn as sns
from sklearn.model_selection import train_test_split 
from sklearn.linear_model import LinearRegression, Ridge,RidgeCV, Lasso, LassoCV  
from sklearn import metrics
from sklearn import linear_model
 
def readcsv(path):
    datas = pd.read_csv(path)
    datas_Sex = datas.Sex
    datas_length = len(datas_Sex)
    head = datas.columns
    return datas ,datas_Sex, head ,datas_length
def main():
    filePath = "abalone.csv"
    file_datas ,file_datas_Sex, file_head ,file_datas_length = readcsv(filePath)
    print("该表头为：",file_head)
    print("共有列数：",len(file_head))
    print("共有行数：",file_datas_length)
    df = pd.DataFrame(data=file_datas)
    sex = []
    for i in file_datas_Sex:
        if i == "M":
            sex.append(2)
        elif i == "F":
            sex.append(1)
        else:sex.append(0)
    df.Sex = sex
    dffr = df.fillna(method="ffill")
    for i in range(file_datas_length):
        for j in range(1,len(file_head)):
            if dffr.values[i][j] >=30 or dffr.values[i][j] <=0:
                dffr.replace(dffr.iloc[i][j],dffr.iloc[i-1][j],inplace = True)
                # days_sorted = sorted(dffr.Diameter.items(), key=lambda d: d[0], reverse=False)
                # print(i,j,dffr.iloc[i][j])
    X =dffr.loc[:, ('Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
       'Visceraweight', 'Shellweight')]
    y = dffr.loc[:, 'Rings']
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=100)
    # print ('X_train.shape={}\ny_train.shape ={}\nX_test.shape={}\ny_test.shape={}'.format(X_train.shape,y_train.shape, X_test.shape,y_test.shape))
    return dffr, X_train, X_test, y_train, y_test

为了方便大家的体验，我把所有的库都先粘进来啦。

好吧，我承认，注释写的确实不咋地（其实我是怕老师认为我需要看注释才能讲出来，嘿嘿）

首先定义了一个readcsv()函数，用来读取csv文件，返回的是整表数据，数据中的Sex列，表头，行数，这样定义会在下面用的时候方便很多。

然后是一个main(）函数，先把该打印的打印出来，然后把数据用pd.DataFrame变成一个dataframe框架（用这个也是很方便）。把Sex列中“M"变为2，"F"变为1，“I”变为0，df.fillna()是把NAN向上填充，下面判断异常值，>=30是自己试了几次之后订的，<=0的数据肯定是不对的，然后得到异常值，因为异常值比较少，所以观察得知他们变量值跟上一行差不多，所以我用上一行的数据进行替代了，一般用平均值或者中位数替代，当然，这都比较好实现，老师也着重问了我这点，我这方法确实也没正常人用过，哈哈哈。

接下来用loc函数把自变量数据和因变量数据分开，再调用train_test_split()把训练集和测试集按7:分开，random_state是随机数种子，100是随便给的，值不同，训练集和测试集的数据就不同，但数据数量相同。终于，得到了我们想要的X_train, X_test, y_train, y_test，dffr是处理好的dataframe框架数据，到此，准备工作全部完成。

异常值（不算NAN)如下

2.调用系统库方法实现线性回归


def xxhg(X_train, X_test, y_train, y_test):             #多元线性回归模型的训练,测试,评价
    # mpl.rcParams['font.sans-serif'] = ['SimHei']  #配置显示中文，否则乱码
    # mpl.rcParams['axes.unicode_minus']=False
    # sns.pairplot(data, x_vars=['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
    #    'Visceraweight', 'Shellweight'], y_vars='Rings',kind="reg", size=5, aspect=0.7)
    # plt.show()  #各变量与Rings的大致关系
    linreg = LinearRegression()
    model=linreg.fit(X_train, y_train)
    # print (model)
    # 训练后模型截距
    print ('The value of intercept:{}'.format(linreg.intercept_))
    # 训练后模型权重
    # print ('训练后模型权重:',linreg.coef_)
    # feature_cols = ['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
    #    'Visceraweight', 'Shellweight']
    # B = list(zip(feature_cols, linreg.coef_))
    # print(B)
    # 预测
    fd = pd.DataFrame({'variable': ['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
                                    'Visceraweight', 'Shellweight'],
                       'weights': linreg.coef_})
    print(fd)
    y_pred = linreg.predict(X_test)
    print('The prediction of test:\n{}'.format(y_pred)) # 1254个样本的预测结果
    #评价
    # 均方根误差(Root Mean Squared Error, RMSE)
    sum_mean=0
    for i in range(len(y_pred)):
        sum_mean+=(y_pred[i]-y_test.values[i])**2
    sum_erro=np.sqrt(sum_mean/1254)
    # calculate RMSE by hand
    print ("the value of RMSE:",sum_erro)
    #做预测值与真实值的折线图
    plt.figure()
    plt.title('test&prediction')
    plt.plot(range(len(y_pred)),y_pred,'b',label="prediction")
    plt.plot(range(len(y_pred)),y_test,'r',label="test")
    plt.legend(loc="upper right") #图中标签的位置
    plt.xlabel("the number of test")
    plt.ylabel('the value of Rings')
    plt.show()

第一个被注释掉的部分是看各个自变量和因变量(Rings)之间的大致关系，如下

第二个被注释掉的是权重的dict(zip())输出，返回值是字典类型。

linreg = LinearRegression() #线性回归模型的实例化
model=linreg.fit(X_train, y_train) #拟合，求权重的计算公式在这个里面

所以上面这两行代码已经解决截距和自变量系数。下面是一些基本函数方法。

linreg.intercept_：求截距

linreg.coef_：求自变量系数

linreg.predict(X_test)：用求得的权重得到的相应的方程对X_test预测,返回其对应的预测值，结果为列表类型。

下面是手写的均方根误差代码部分

#评价
# 均方根误差(Root Mean Squared Error, RMSE)
sum_mean=0
for i in range(len(y_pred)):
    sum_mean+=(y_pred[i]-y_test.values[i])**2
sum_erro=np.sqrt(sum_mean/1254)
# calculate RMSE by hand
print ("the value of RMSE:",sum_erro)

画图部分是画的线性图:

最后结果显示：

3.调用系统库方法实现岭回归


def linghg(X_train, X_test, y_train, y_test):       #岭回归模型的训练,测试,评价
    n_alphas = 200
    alphas = np.logspace(-3,4, n_alphas) #对数等比数列
    clf = linear_model.Ridge(fit_intercept=True)
    coefs = []
    score = []
    intercept = []
    c = 0
    for a in alphas:
        clf.set_params(alpha=a)
        clf.fit(X_train, y_train)
        coefs.append(clf.coef_)
        m = clf.score(X_test,y_test)
        n = clf.intercept_
        score.append(m)
        intercept.append(n)
    for i in range(len(score)):
        if score[i] == max(score):
            c = alphas[i]
            fd = pd.DataFrame({'variable': ['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
                                            'Visceraweight', 'Shellweight'],
                               'weights': coefs[i]})
            print(fd)
            print('intercept:{}\nalphas:{}\nmax(score):{}'.format(intercept[i],alphas[i],score[i]))
    clf1 = linear_model.Ridge(fit_intercept=True)
    clf1.set_params(alpha=c)
    clf1.fit(X_train, y_train)
    y_pre = clf1.predict(X_test)
    # 评价
    # 均方根误差(Root Mean Squared Error, RMSE)
    sum_mean = 0
    for i in range(len(y_pre)):
        sum_mean += (y_pre[i] - y_test.values[i]) ** 2
    sum_erro = np.sqrt(sum_mean / len(y_pre))
    print('The value of RMSE in max(score):', sum_erro)
    print('***************************************************')
    # 交叉验证
    clf2 = linear_model.RidgeCV(alphas=alphas)
    clf2.fit(X_train, y_train)
    fr = pd.DataFrame({'variable': ['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
                                    'Visceraweight', 'Shellweight'],
                       'weights': clf2.coef_})
    print(fr)
    print('intercept:{}\nalphas:{}'.format(clf2.intercept_,clf2.alpha_))
    y_pred = clf2.predict(X_test)
    # print(y_pred)  # 1254个样本的预测结果
    # 评价
    # 均方根误差(Root Mean Squared Error, RMSE)
    sum_mean = 0
    for i in range(len(y_pred)):
        sum_mean += (y_pred[i] - y_test.values[i]) ** 2
    sum_erro = np.sqrt(sum_mean / len(y_pred))
    print('the value of RMSE in cross validation:', sum_erro)
    #图形展示
    #设置刻度
    ax = plt.gca()
    #设置刻度的映射
    ax.plot(alphas, score)
    # ax.plot(alphas,coefs)
    #设置x轴的刻度显示方式
    ax.set_xscale('log')
    # #翻转x轴
    # ax.set_xlim(ax.get_xlim()[::-1])
    #设置x、y标签以及标题
    plt.xlabel('alpha')
    plt.ylabel('score')
    # plt.ylabel(('weights'))
    plt.title('Ridge score&alphas')
    # plt.title('Ridge coefficients as a function of the regularization')
    #使得坐标轴最大值和最小值与数据保持一致
    plt.axis('tight')
    plt.show()

这部分跟线性回归相比多加了alpha*第二范式,其实很差不多，函数方法功能一样，名称也差不多。

clf = linear_model.Ridge(fit_intercept=True)
clf.set_params(alpha=a)
clf.fit(X_train, y_train)

三行解决权重。

m = clf.score(X_test,y_test) #其实就是下面这个方程啦

上面的循环目的是找到当m最大时，得到相应的权重，alpha值，max（score)值（即最大的m值)，然后再实例化clf1,把此时的alphas值带入，得到X_text的预测值，再求一个均方误差。按理说到这就结束了，但我想再来点，于是乎：

clf2 = linear_model.RidgeCV(alphas=alphas)
clf2.fit(X_train, y_train)

linear_model.RidgeCV()实例化之后，拟合之后能够直接在alpas（注意是一组数）中找到模型最优时相应的权重，alpha值，（当然它是用了他自己个的方法来判断模型是不是最优，但不是score方法，因为两者得到的权重有差异，也不是均方误差最小,后面的手写代码会有验证。）然后得到对应的预测值，均方误差值。

结果如下：

4.调用系统库方法实现lasso回归


def lassohg(X_train, X_test, y_train, y_test):          #lasso回归模型的训练,测试,评价
    n_alphas = 200
    alphas = np.logspace(-5, 0, n_alphas)  # 对数等比数列
    coefs = []
    score = []
    intercept = []
    c = 0
    las = Lasso(normalize=True)
    for a in alphas:
        las.set_params(alpha=a)
        las.fit(X_train, y_train)
        coefs.append(las.coef_)
        m = las.score(X_test, y_test)
        n = las.intercept_
        score.append(m)
        intercept.append(n)
    for i in range(len(score)):
        if score[i] == max(score):
            c = alphas[i]
            fd = pd.DataFrame({'variable': ['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
                                            'Visceraweight', 'Shellweight'],
                               'weights': coefs[i]})
            print(fd)
            print('intercept:{}\nalphas:{}\nmax(score):{}'.format(intercept[i], alphas[i], score[i]))
    las1 = Lasso(normalize=True)
    las1.set_params(alpha = c)
    las1.fit(X_train, y_train)
    y_pre = las1.predict(X_test)
    # 评价
    # 均方根误差(Root Mean Squared Error, RMSE)
    sum_mean = 0
    for i in range(len(y_pre)):
        sum_mean += (y_pre[i] - y_test.values[i]) ** 2
    sum_erro = np.sqrt(sum_mean / len(y_pre))
    print('The value of RMSE in max(score):', sum_erro)
    print('***************************************************')
    # 交叉验证，找到模型最优的alphas值
    lasso_cv = LassoCV(alphas=alphas, normalize=True, max_iter=1000, cv=None)
    lasso_cv.fit(X_train, y_train)
    fr = pd.DataFrame({'variable': ['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
                                    'Visceraweight', 'Shellweight'],
                       'weights': lasso_cv.coef_})
    print(fr)
    print('intercept:{}\nalphas:{}'.format(lasso_cv.intercept_, lasso_cv.alpha_))
    y_pred = lasso_cv.predict(X_test)
    # print(y_pred)  # 1254个样本的预测结果
    # 评价
    # 均方根误差(Root Mean Squared Error, RMSE)
    sum_mean = 0
    for i in range(len(y_pred)):
        sum_mean += (y_pred[i] - y_test.values[i]) ** 2
    sum_erro = np.sqrt(sum_mean / 1254)
    print('The value of RMSE in cross validation:', sum_erro)
 
    # 图形展示
    # 设置刻度
    ax = plt.gca()
    # 设置刻度的映射
    ax.plot(alphas, score)
    # ax.plot(alphas,coefs)
    # 设置x轴的刻度显示方式
    ax.set_xscale('log')
    # # 翻转x轴
    # ax.set_xlim(ax.get_xlim()[::-1])
    # 设置x、y标签以及标题
    plt.xlabel('alpha')
    plt.ylabel('score')
    # plt.ylabel(('weights'))
    plt.title('Lasso score&alphas')
    # plt.title('Lasso coefficients as a function of the regularization')
    # 使得坐标轴最大值和最小值与数据保持一致
    plt.axis('tight')
    plt.show()

大家一看代码，是不是觉得代码似曾相识，哈哈哈哈哈，答对了，因为lasso回归和岭回归除了一些函数方法名称不一样，其他都是一样的，我在写代码的时候也是直接复制粘贴，然后改一些方法名。lasso回归方法跟线性回归相比多加了alpha*第一范式。代码解释可参考岭回归解释。

结果如下：

5.对三种方法的实现


if __name__ == '__main__':
    data,X_train, X_test, y_train, y_test = main()
    xxhg(X_train, X_test, y_train, y_test)
    print('############################################################')
    linghg(X_train, X_test, y_train, y_test)
    print('############################################################')
    lassohg(X_train, X_test, y_train, y_test)

放到一个文件里面直接运行，打印出来的效果还是可以的，出个结果出个图，非常的方便。

岭回归和lasso回归画图部分还隐藏了score&alphas的图像哦。

6.手写代码实现线性回归


import numpy as np
import pandas as pd
from sklearn.model_selection import train_test_split
 
class LinearRegression_SelfDefined():
    def __init__(self):
        self.w = None
        self.i = None
        self.y_pred = None
 
    def fit(self, x, y):
        m = x.shape[0]
        X = np.concatenate((np.ones((m, 1)), x), axis=1)
        xMat = np.mat(X)
        yMat = np.mat(y.values.reshape(-1, 1))
        xTx = xMat.T * xMat
        # xTx.I为xTx的逆矩阵
        allparameters = xTx.I * xMat.T * yMat
        self.i = allparameters[0]
        self.w = allparameters[1:]
 
    def get(self):
        W = np.array(self.w)
        WW = np.ravel(W)
        I = np.array(self.i)
        II = np.ravel(I)
        fd = pd.DataFrame({'variable': ['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight','Visceraweight', 'Shellweight'],
                           'weights': WW})
        print(fd)
        print('intercept:', II[0])
 
    def predict(self,x):
        # m = x.shape[0]
        # X = np.concatenate((np.ones((m, 1)), x), axis=1)
        self.y_pred = x.dot(self.w)  # X与self.w所表示的矩阵相乘
        y = np.array(self.y_pred)
        i = np.array(self.i)
        ii = np.ravel(i)
        j = []
        yy = np.ravel(y)
        for i in yy:
            j.append(i + ii[0])
        self.y_pred = np.array(j)
        return self.y_pred
 
    def linearScore(self,y):
        u = 0
        v = 0
        sum = 0
        y = np.array(y)
        y_pre = np.array(self.y_pred)
        for i in y:
            sum += i
        mean = sum/len(y)
        for i in range(len(y)):
            u += (y[i]-y_pre[i])**2
            v += (y[i]-mean)**2
        R = 1-u/v
        return R
 
    def linearRMSE(self,y):
        # 评价
        # 均方根误差(Root Mean Squared Error, RMSE)
        y = np.array(y)
        y_pre = np.array(self.y_pred)
        sum_mean = 0
        for i in range(len(y_pre)):
            sum_mean += (y_pre[i] - y[i]) ** 2
        sum_erro = np.sqrt(sum_mean / len(y))
        return sum_erro

对于这部分呢，咱们已经知道系统库中的fit(),predict()方法，所以以此为基础，我们也来写相应的方法。首先定义LinearRegression_SelfDefined()类，心中明确fit()函数功能，predict()函数功能，整体框架写完之后，就开始补充了，开干。

self.w = None   #各变量系数
self.i = None   #截距
self.y_pred = None   #X_test的预测值

fit()函数中实现求权值的方程，这个需要了解一些np库中相应的矩阵方法，返回值是权重（截距+各变量系数），此时返回值类型为矩阵形式。然后是predict()函数，首先将各变量系数与测试集数据矩阵求和，得到每行数据的一个预测结果（未加截距），此时数据还是矩阵形式，然后我用np.array(),np.reval()把矩阵形式变成数组，然后降维，再通过数组计算得到最终预测结果。按理说已经可以结束了，但我，还是想再加点，毕竟写到这了，再把评价函数加进去吧，反正也不难，所以，我加了linearScore(self,y)函数，linearRMSE(self,y)函数，方法中的公式已经是老生常谈了。然后，可以结束了，但我觉得输出还能再美化点，所以，我加了一个get()函数。不得不说，dataframe格式的数据输出出来就是强一点。结果如下：

手写代码岭回归（*****上面），手写代码lasso回归（*******下面）的最终结果

这张图片****上面是手写代码的线性回归结果，*****下面是手写代码的岭回归结果。

对比调用系统函数的线性回归模型结果，发现，竟然！竟然！一模一样！好吧，也没什么可惊讶的，说明库函数中的方法就是用的咱们刚写的公式法。

8.手写代码实现岭回归


class RidgeRegression_SelfDefined():
    def __init__(self):
        self.w = None
        self.i = None
        self.y_pred = None
        self.rmse = None
        self.alphas = None
 
    def fit(self, x, y,alphas):
        m = x.shape[0]
        X = np.concatenate((np.ones((m, 1)), x), axis=1)
        xMat = np.mat(X)
        yMat = np.mat(y.values.reshape(-1, 1))
        xTx = xMat.T * xMat
        q = 2.5
        for alpha in alphas:
            denom = xTx + np.eye(np.shape(xMat)[1]) * alpha
            # if np.linalg.det(denom) == 0.0:
            #     print("This matrix is singular, cannot do inverse")
            #     return
            allparameters = denom.I * (xMat.T * yMat)
            o = allparameters[0]
            p = allparameters[1:]
            a = RidgeRegression_SelfDefined()
            y_pred = a.predictn(x,o,p)
            if a.RidgeRMSEn(y,y_pred)<q:
                q = a.RidgeRMSEn(y,y_pred)
                self.i = o
                self.w = p
                self.y_pred = y_pred
                self.rmse = q
                self.alphas = alpha
    def get(self):
                W = np.array(self.w)
                WW = np.ravel(W)
                I = np.array(self.i)
                II = np.ravel(I)
                alp = self.alphas
                # feature_cols = ['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight',
                #    'Visceraweight', 'Shellweight']
                # B = dict(zip(feature_cols, WW))
                # print(B)
 
                fd = pd.DataFrame({'variable':['Sex', 'Length', 'Diameter', 'Height', 'Wholeweight', 'Shuckedweight','Visceraweight', 'Shellweight'],
                'weights':WW})
                print(fd)
                print('intercept:',II[0])
                print('alphas:',alp)
 
    def predictn(self,x,o,p):
        y_pred = x.dot(p)  # X与self.w所表示的矩阵相乘
        y = np.array(y_pred)
        i = np.array(o)
        ii = np.ravel(i)
        j = []
        yy = np.ravel(y)
        for m in yy:
            j.append(m + ii[0])
        y_pred = np.array(j)
        return y_pred
 
    def predict(self,x):
        # m = x.shape[0]
        # X = np.concatenate((np.ones((m, 1)), x), axis=1)
        self.y_pred = x.dot(self.w) # X与self.w所表示的矩阵相乘
        y = np.array(self.y_pred)
        i = np.array(self.i)
        ii = np.ravel(i)
        j = []
        yy = np.ravel(y)
        for i in yy:
            j.append(i + ii[0])
        self.y_pred = np.array(j)
        return self.y_pred
 
    def RidgeRMSEn(self,y,y_pred):
        # 评价
        # 均方根误差(Root Mean Squared Error, RMSE)
        y = np.array(y)
        sum_mean = 0
        for i in range(len(y_pred)):
            sum_mean += (y_pred[i] - y[i]) ** 2
        sum_erro = np.sqrt(sum_mean / len(y))
        return sum_erro
 
    def RidgeScore(self,y):
        u = 0
        v = 0
        sum = 0
        y = np.array(y)
        y_pre = np.array(self.y_pred)
        for i in y:
            sum += i
        mean = sum/len(y)
        for i in range(len(y)):
            u += (y[i]-y_pre[i])**2
            v += (y[i]-mean)**2
        R = 1-u/v
        return R
 
    def RidgeRMSE(self,y):
        # 评价
        # 均方根误差(Root Mean Squared Error, RMSE)
        y = np.array(y)
        y_pre = np.array(self.y_pred)
        sum_mean = 0
        for i in range(len(y_pre)):
            sum_mean += (y_pre[i] - y[i]) ** 2
        sum_erro = np.sqrt(sum_mean / len(y))
        return sum_erro

对于这部分，我想多说点，因为这个是我花时间最多的，虽然很大一部分是因为我弄错了alpha的值，我还以为输入alpha的值是一个数组，气呼呼。不过，显然更大一部分原因是我代码水平太差（小白努力中.......）。刚开始，还是老套路，先把fit()函数,predict()函数要实现的功能整明白，写上大致框架，然后，开干！

在fit()函数中公式部分实现没啥问题，等到实例化时我意识到一个问题：我现在还不知道alpha的值，对于这个alphas数组，我不知道alpha要取其中的哪个值，这咋办？我只能开始循环，调用自身predict()函数，RidgeRMSE()函数，找到让均方误差最小的那个alpha,然后得到权重，进而得到想要的模型。但结果发现，不行啊，参数和方法都是带self的，一个改变全都改变，有种按着葫芦起了瓢的感觉。所以，索性，我又定义了两个函数：predictn(),RidgeRMSEn()函数，把参数和返回值的self统统去掉，调用之后得到让均方误差最小的那个alpha，从而得到模型。

结果看上图。lasso手写代码跟岭回归除了名称外，还有一点不一样，就是fit()函数的矩阵公式，其他都是一样的大家可以尝试一下，这里就不展示了。

9.两种方法的实现


if __name__ == '__main__':
    data,X_train, X_test, y_train, y_test = main()
 
    B = LinearRegression_SelfDefined()
    B.fit(X_train,y_train)
    B.get()
    print('The prediction of test:\n{}'.format(B.predict(X_test)))
    print('The score of this model:{}'.format(B.linearScore(y_test)))
    print('The value of RMSE:{}'.format(B.linearRMSE(y_test)))
    print('******************************************************')
    n_alphas = 200
    alphas = np.logspace(-3,4, n_alphas) #对数等比数列
    C = RidgeRegression_SelfDefined()
    C.fit(X_train,y_train,alphas)
    C.get()
    print('The prediction of test:\n{}'.format(C.predict(X_test)))
    print('The score of this model:{}'.format(C.RidgeScore(y_test)))
    print('The value of RMSE:{}'.format(C.RidgeRMSE(y_test)))

6-9的代码放一块运行，结果还是不错的，实验到此结束。

最后，还记得开头说的优化点吗？

1.两个手写代码的pridict()方法中,可以直接用X_test数据矩阵（得在X_test矩阵中第一列加入一列1）*权重矩阵，就不用化成数组倒腾了。

2.手写岭回归代码的fit()函数的q我定义的是2.5，是因为我知道最后的比较小的RMSE结果是2.22，这里可以先从alphas数组中取出一个alpha，得到一个RMSE赋给q，这样模型就比较通用了。

3.还是想想怎么调用自身函数找出alphas数组中最好的alpha值，为了一个fit()再去定义两个功能函数不好，以后要调用更多的功能函数怎么办？还是每个功能函数写两份？哈哈哈哈，个人感觉不是很可取。

OK，分析最终结果，我们已经得到了最优的那个模型，也能证实多元线性回归模型用的是公式法，但RidgeCV,LassoCV如何找到alphas的公式还不知道。所以，笔者就写到这吧，去看库函数源码了，肯定比自己手写的好得多。手写代码之前没看源码哈，是自己个根据库函数方法，参数一点点写的，不然感觉像糊弄老师，哈哈哈，其实也是想自己先实现一下，毕竟，这是机器学习里面最简单的模型了。

声明：本文内容由网友自发贡献，转载请注明出处：【wpsshop博客】