用scikit-learn训练线性回归模型

---

1. 读取数据和预处理函数

import pandas as pd

# 从excel文件中读取数据
def loaddata(datafile):
    # 读取自变量和因变量数据  index:温度数据0-12   DYB数据13-38  挠度数据39-40
    df = pd.read_excel(datafile, header=None)
    data = df.values
    data = np.delete(data, 0, 0)
    data = np.delete(data, [0, 1], 1)
    data = data.astype(np.float)
    # 温度数据
    X1 = data[:, 0:13]
    # 动应变数据
    t1 = []
    t2 = []
    for i in range(26):
        if i % 2 == 0:
            t1.append(13 + i)
        else:
            t2.append(13 + i)
    X21 = data[:, t1]
    X22 = data[:, t2]
    # 挠度数据
    Y11 = data[:, 39:40]
    Y12 = data[:, 40:41]
    return X1, X21, X22, Y11, Y12

# 根据数据集Data.xlsx
def main():
    datafile = "./data/Data.xlsx"
    X1, X21, X22, Y11, Y12 = loaddata(datafile)
    # 完整正值数据
    X11 = np.concatenate([X1, X21], axis=1)
    # 完整负值数据
    X12 = np.concatenate([X1, X22], axis=1)
    return X11, X12, Y11, Y12

---

2. 准备算法的输入数据，划分训练集和测试集

if __name__ == "__main__":
    X11, X12, Y11, Y12 = main()
    X_data_1 = X11  # (677,26)
    Y_data_1 = Y11  # (677,1)

    X_train, X_test, y_train, y_test = train_test_split(X_data_1, Y_data_1, test_size=0.3, random_state=2018)
    print(np.array(X_train).shape)  # (473,26)
    print(np.array(X_test).shape)  # (204,26)

train_test_split()的参数解释：

• train_data：被划分的样本特征集
• train_target：被划分的样本标签
• test_size：如果是浮点数，在0-1之间，表示样本占比；如果是整数的话就是样本的数量
• random_state：是随机数的种子。其实就是该组随机数的编号，在需要重复试验的时候，保证得到一组一样的随机数。比如你每次都填1，其他参数一样的情况下你得到的随机数组是一样的。但填0或不填，每次都会不一样。
  随机数的产生取决于种子，随机数和种子之间的关系遵从以下两个规则：
  种子不同，产生不同的随机数；种子相同，即使实例不同也产生相同的随机数。

---

3. 训练scikit-learn的线性模型

from sklearn.linear_model import LinearRegression

# 下面为sklearn线性回归过程
    linreg = LinearRegression()
    linreg.fit(X_train, y_train)
# 查看模型的系数结果
    print(linreg.intercept_)
    print(linreg.coef_)

如果输出如下：
[ 447.06297099]
[[-1.97376045 -0.23229086 0.0693515 -0.15806957]]
那么我们就得到了线性回归模型：
Y=447.06297099−1.97376045∗X1−0.23229086∗X2+0.0693515∗X3−0.15806957∗X4

---

4. 模型预测和评价

    # 预测
    y_pred = linreg.predict(X_test)

    # 用scikit-learn计算MSE
    print("一次交叉验证之后结果---------------------")
    print("MSE:",metrics.mean_squared_error(y_test, y_pred))
    # 用scikit-learn计算RMSE
    print("RMSE:", np.sqrt(metrics.mean_squared_error(y_test, y_pred)))
    print("R2_score:", metrics.r2_score(y_test, y_pred))

    # 我们可以通过交叉验证来持续优化模型，代码如下，我们采用10折交叉验证，即cross_val_predict中的cv参数为10：
    predicted = cross_val_predict(linreg, X_data_1, Y_data_1, cv=10)
    # 用scikit-learn计算MSE
    print("十折交叉验证之后结果---------------------")
    print("MSE:",metrics.mean_squared_error(Y_data_1, predicted))
    # 用scikit-learn计算RMSE
    print("RMSE:", np.sqrt(metrics.mean_squared_error(Y_data_1, predicted)))
    print("explained_variance_score:", metrics.explained_variance_score(Y_data_1, predicted))
    print("mean_absolute_error:", metrics.mean_absolute_error(Y_data_1, predicted))
    print("median_absolute_error:", metrics.median_absolute_error(Y_data_1, predicted))
    print("R2_score:", metrics.r2_score(Y_data_1, predicted))

自己实现的验证：

''' 均方误差根 '''
def rmse(y_test, y_true):
    return sp.sqrt(sp.mean((y_test - y_true) ** 2))


''' 与均值相比的优秀程度，介于[0~1]。0表示不如均值。1表示完美预测.这个版本的实现是参考scikit-learn官网文档 '''
def R2(y_test, y_true):
    return 1 - ((y_test - y_true) ** 2).sum() / ((y_true - np.mean(y_true)) ** 2).sum()


''' 这是Conway&White《机器学习使用案例解析》里的版本 '''
def R22(y_test, y_true):
    y_mean = np.array(y_true)
    y_mean[:] = y_mean.mean()
    return 1 - rmse(y_test, y_true) / rmse(y_mean, y_true)

    # print("RMSE_1:",rmse(y_pred , y_test))
    # print("R2_score_1:",R2(y_pred , y_test))
    # print("R22_score_1:",R22(y_pred , y_test))

sklearn中的回归器性能评估方法：

1. 可释方差得分（explained_variance_score）

Var是方差，y是真实值，y扒是预测值。
可释方差得分越接近1，代表预测的目标响应越好，越小则代表结果越差。

2. 平均绝对值误差（mean_absolute_error）

3. 均方误差（mean-squared-error）

4. 中值绝对误差（Median absolute error）

5. R2 决定系数（拟合优度）

它是表征回归方程在多大程度上解释了因变量的变化，或者说方程对观测值的拟合程度如何。 因为如果单纯用残差平方和会受到你因变量和自变量绝对值大小的影响，不利于在不同模型之间进行相对比较.而用拟合优度就可以解决这个问题。例如一个模型中的因变量:10000、20000…..，而另一个模型中因变量为1、2……，这两个模型中第一个模型的残差平方和可能会很大,而另一个会很小，但是这不能说明第一个模型就别第二个模型差。

5. 画图观察效果

    fig, ax = plt.subplots()
    ax.scatter(Y_data_1, predicted)
    ax.plot([np.min(Y_data_1), np.max(Y_data_1)], [np.min(Y_data_1), np.max(Y_data_1)], 'k--', lw=4)
    ax.set_xlabel('Measured')
    ax.set_ylabel('Predicted')
    plt.show()

    plt.plot(Y_data_1, label="True")
    plt.plot(predicted, label="Predicted")
    plt.legend()
    plt.show()