机器学习基于加州房价的线性回归实验

柔情只为你懂 2023-09-30 11:40 6阅读 0赞

### 1.线性回归闭合形式参数求解的原理 ###

如果定义X为m\*(n+1)的矩阵，Y为m*1的矩阵，θ为(n+1)1维的矩阵，那么在之前的定义中就可以表示为h(x)=Xθ。则代价函数可以表示为J(θ)=1/2*(X*θ-y)Т*(X*θ-y),J(θ)为凹函数，我们要让其值最小化，只需对该函数求导，然后令导数为0即可求得θ。对其求导后得到XT*X*θ-XT*y，令其等于0，得到θ=(XT*X)^-1*XT\*y。

### 2.线性回归梯度下降参数求解的原理 ###

我们构造了拟合函数h(θ),并且得到了损失函数J（θ），我们要求得使J（θ）取得最小值的θ，其原理还是求偏导然后使导数为0，我们对J（θ）求导得到(hθ(x) − y) xj，然后可以得到对θj的更新公式

![在这里插入图片描述][c34a0d97b9ad4af6a2a8854284a92412.png]

由于数据量较大，所以采用了随机梯度下降，但是准确度相较于批量梯度下降来说会下降。在我的程序里，由于数据采用矩阵形式存储，所以更新过程可以替换为  
![在这里插入图片描述][80148e47b67a46e7b359297d5d5bc06f.png]  
其中θ为(n+1)*1维，X为m*(n+1)维，Y为m\*1维。梯度为0用损失函数差值小于1e-18来表示，说明这个点是损失函数的极小值点，但并不一定是最小值点。

### 3.[相关文件][Link 1] ###

### 4.程序清单 ###

相关包：

from sklearn.preprocessing import StandardScaler
    from sklearn.model_selection import train_test_split
    import pandas as pd
    import numpy as np
    import os
    import time
    from numpy import median
    from sklearn.preprocessing import OneHotEncoder

（一）读取数据：

# 读取数据
    HOUSE_PATH = './'
    
    def load_housing_data(housing_path=HOUSE_PATH):
        csv_path = os.path.join(housing_path, 'housing.csv')
        return pd.read_csv(csv_path)
    housing = load_housing_data()

（二）数据处理：

# 将中位数补全空位
    median = housing["total_bedrooms"].median()
    housing["total_bedrooms"].fillna(median, inplace=True)
    
    # 独热编码
    housing_category = housing[["ocean_proximity"]]
    cat_encoder = OneHotEncoder()
    housing_category_onehot = cat_encoder.fit_transform(housing_category)
    housing = housing.drop("ocean_proximity", axis=1)
    housing_values = np.c_[housing.values, housing_category_onehot.toarray()]
    housing_fixed = pd.DataFrame(
        housing_values,
        columns=list(housing.columns) +
        ['<1H OCEAN', 'INLAND', 'ISLAND', 'NEAR BAY', 'NEAR OCEAN'],
        index=housing.index
    )

（三）分析数据相关性

# 分析数据相关性
    corr_matrix = housing_fixed.corr()   # 用corr计算两两特征之间的相关性系数
    correlation = corr_matrix["median_house_value"].sort_values(ascending=False)   # 跟街区价格中位数特征的其他特征的相关系数
    print(correlation)

（四）将数据集分类

# 将数据集分类
    train_set, test_set = train_test_split(
        housing_fixed, test_size=0.3, random_state=42)
    X_train = train_set[['longitude', 'latitude', 'housing_median_age',
                         'total_rooms', 'total_bedrooms', 'population', 'households', 'median_income', '<1H OCEAN', 'INLAND', 'ISLAND',
                         'NEAR BAY', 'NEAR OCEAN']]
    y_train = train_set[["median_house_value"]]
    X_test = test_set[['longitude', 'latitude', 'housing_median_age',
                       'total_rooms', 'total_bedrooms', 'population', 'households', 'median_income', '<1H OCEAN', 'INLAND', 'ISLAND',
                       'NEAR BAY', 'NEAR OCEAN']]
    y_test = test_set[["median_house_value"]]
    X = np.hstack([np.ones((len(X_train_std), 1)), X_train_std])  # 训练集X
    Y = np.array(y_train_std)  # 训练集Y
    x = np.hstack([np.ones((len(X_test_std), 1)), X_test_std])  # 测试集x
    y = np.array(y_test_std)  # 测试集y
    y_var = np.var(y)  # 标准差

（五）特征标准化

# 特征标准化
    stdsc = StandardScaler()
    X_train_std = stdsc.fit_transform(X_train)
    X_test_std = stdsc.transform(X_test)
    y_train_std = stdsc.fit_transform(y_train)
    y_test_std = stdsc.fit_transform(y_test)

（六）初始化θ

theta = np.zeros((14, 1))  # 初始化theta

（七）正规方程

# 正规方程
    def nomal(X, Y):
        theta = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(Y)
        return theta

# 损失函数
    def cost_function(x, theta, y):
        cost = np.sum((np.dot(x, theta)-y)**2)
        return cost/(2*len(y))
    
    
    # 正规方程
    def nomal(X, Y):
        theta = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(Y)
        return theta
    
    
    # 梯度下降方向
    def gradient(X, theta, Y):
        return X.T.dot((X.dot(theta)-Y))/len(Y)
    
    
    # 梯度下降
    def gradient_descent(X, theta, Y, eta):
        while True:
            last_theta = theta
            grad = gradient(X, theta, Y)
            theta = theta - eta*grad
            if abs(cost_function(X, last_theta, Y) - cost_function(X, theta, Y)) < 1e-18:
                break
        return theta
    
    
    # theta = nomal(X, Y)  # 闭合形式求解
    theta = gradient_descent(X, theta, Y, 0.001)  # 梯度下降
    
    
    # 评估项(R2)
    def evaluation(x, theta, y, y_var):
        return 1 - ((np.sum((np.dot(x, theta)-y)**2))/(y_var*len(y)))
    
    
    MSE = np.sum(np.power((np.dot(x, theta)-y), 2))/len(y)
    cost = cost_function(x, theta, y)  # 损失函数值
    R2 = evaluation(x, theta, y, y_var)  # 评估值
    
    end = time.time()
    
    # print("The normal equations:")
    print("Gradient descent:")
    print("theta=")
    print(theta)
    print("MSE=", MSE)
    print("cost=", cost)
    print("R2=", R2)
    print('Running time: %s Seconds' % (end-start))

实验结果：  
（一）正规方程求解  
![在这里插入图片描述][watermark_type_d3F5LXplbmhlaQ_shadow_50_text_Q1NETiBA6IKH5LqL5Y-45py66LW15Zub_size_18_color_FFFFFF_t_70_g_se_x_16]  
（二）梯度下降求解  
![在这里插入图片描述][watermark_type_d3F5LXplbmhlaQ_shadow_50_text_Q1NETiBA6IKH5LqL5Y-45py66LW15Zub_size_18_color_FFFFFF_t_70_g_se_x_16 1]  
可以看到两个方法得出的结果差别不大，用测试集进行测试时候，损失函数值均为0.18，评估项R²均为0.6多，梯度下降的拟合效果会比正规方程的好一点。在运算过程中，能很明显看到正规方程的计算速度要比梯度下降快很多，原因在于梯度下降在更新θ时候需要迭代很多次才能得到较优解。但是梯度下降在特征数量n较大时也能很好使用，而正规方程需要计算(X\*X)-1，如果特征数量太多则运算代价较大因为矩阵的运算时间复杂度为O(n³)，而且只适用于线性模型，不适用于逻辑回归模型等其他模型。在这个模型里面，由于特征数量不是很多，因此用正规式求解比较合理。

[c34a0d97b9ad4af6a2a8854284a92412.png]: https://img-blog.csdnimg.cn/c34a0d97b9ad4af6a2a8854284a92412.png
[80148e47b67a46e7b359297d5d5bc06f.png]: https://img-blog.csdnimg.cn/80148e47b67a46e7b359297d5d5bc06f.png
[Link 1]: https://download.csdn.net/download/weixin_45750572/85049514
[watermark_type_d3F5LXplbmhlaQ_shadow_50_text_Q1NETiBA6IKH5LqL5Y-45py66LW15Zub_size_18_color_FFFFFF_t_70_g_se_x_16]: https://img-blog.csdnimg.cn/f39436c7cf5e433fa46d82e22017464d.png?x-oss-process=image/watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBA6IKH5LqL5Y-45py66LW15Zub,size_18,color_FFFFFF,t_70,g_se,x_16
[watermark_type_d3F5LXplbmhlaQ_shadow_50_text_Q1NETiBA6IKH5LqL5Y-45py66LW15Zub_size_18_color_FFFFFF_t_70_g_se_x_16 1]: https://img-blog.csdnimg.cn/09dc7f112e6542058b426d2ea8ff44e9.png?x-oss-process=image/watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBA6IKH5LqL5Y-45py66LW15Zub,size_18,color_FFFFFF,t_70,g_se,x_16