ex1线性回归

单变量回归问题

​ 我们将使用一个变量实现线性回归，根据城市人口数量，预测开小吃店的利润，数据在ex1data1.txt里，第一列是城市人口数量，第二列是该城市小吃店利润。

读取数据

利用pandas读取数据，并显示前5行数据检查是否读取成功。

import pandas as pd
data = pd.read_csv('ex1data1.txt', header=None, names=['Population', 'Profit'])
print(data.head())

	Population	Profit
0	6.1101	17.5920
1	5.5277	9.1302
2	8.5186	13.6620
3	7.0032	11.8540
4	5.8598	6.8233

数据可视化

第一步我们需要对数据进行可视化，初步判断数据之间是否存在线性关系，这对后面的分析十分重要。

import matplotlib.pyplot as plt
plt.scatter(x=data['Population'], y=data['Profit'], label='data', marker='.', color='b')
plt.legend()
plt.show()

数据处理

根据吴恩达课程所讲，我们需要在x的第一列加上一列x_0=1,初始化初值为0的theta。

import numpy as np
data.insert(0, 'ones', 1)  # 在第0列插入1，表头为ones 
x = data.iloc[:, 0:-1]  # 切片工具[几行到几行, 几列到几列]
y = data.iloc[:, -1:]
x, y = np.matrix(x.values), np.matrix(y.values)
theta = np.matrix(np.zeros(x.shape[1]))  # 初始化theta

代价函数

采用回归问题中常用的平方误差代价函数。

$$J（\Theta）= \frac{1}{2m}\times\sum_{i=1}^m(h_\Theta(x^i)-y^i)^2$$

其中

$$h_\Theta(x)=\Theta_0+\Theta_1\times x_1$$

注意在前面初始化的时候添加了一列x_0,故

$$h_\Theta(x)=x\cdot\Theta^T$$

我们的目标是求出能够最小化J(θ)的θ，这样预测值才能更加接近y，即:

$$argmin_\Theta=J(\Theta)$$

import numpy as np
# 计算代价函数
def cost_function(x_, y_, theta_):
    ans = np.power(x_.dot(theta_.T)-y_, 2)
    return np.sum(ans)/(2*x_.shape[0])

print(cost_function(x, y, theta))

得出的结果应为theta=32.072733877455676。

梯度下降

$$\Theta_j:=\Theta_j-\alpha\times\frac{\partial J(\Theta)}{\partial \Theta_j} :=\Theta_j-\frac{1}{m}\times \sum_{i=1}^m(h_\Theta(x^i)-y^i)x^i$$

其中α为learning rate,用来控制梯度下降的幅度。

def gradient_descent(x_, y_, theta_, alpha_, inters_):
    temp = np.matrix(np.zeros(theta.shape))
    nums = int(theta_.shape[1])
    cost = np.zeros(inters_)
    for inter in range(inters_):
        ans = x_.dot(theta_.T)-y_
        for i in range(nums):
            temp[0, i] = theta_[0, i] - alpha_*np.sum(np.multiply(ans, x_[:, i]))/len(x_)
        theta_ = temp
        cost[inter] = cost_function(x_, y_, theta_)
    return theta_, cost

预测函数

设置α以及梯度下降的迭代次数，对θ进行求解。

# 设置参数
alpha = 0.01
inters = 1500
g, cost = gradient_descent(x, y, theta, alpha, inters)
# 画出预测函数
x_test = np.linspace(data['Population'].min(), data['Population'].max(), data.shape[0])
# linspace(start,end,num=50)从间隔取出50个数
predict = g[0, 0]+g[0, 1]*x_test
plt.scatter(x=data['Population'], y=data['Profit'], label='data', marker='.', color='b')
plt.plot(x_test, predict, c='r', label='predict')
plt.xlabel('Population')
plt.ylabel('Profit')
plt.legend()
plt.show()

得出最后的预测函数如下图所示。

多变量回归问题

​ 在这一部分中，我们将使用多个变量实现线性回归以预测房屋价格。假设你正在出售房屋，并且想知道一个好的市场价格。一种方法是首先收集最近有关出售房屋的信息，并建立房屋价格模型。
​ 文件ex1data2.txt包含某地区房屋价格的训练集。第一列是房屋的大小（以平方英尺为单位），第二列是卧室的数量，第三列是房屋的价格。

读取数据

data2 = pd.read_csv('ex1data2.txt', header=None, names=['Size', 'Bedrooms', 'Price'])
print(data2.head())

显示前5行数据如下，

	Size	Bedrooms	Price
0	2104	3	399900
1	1600	3	329900
2	2400	3	369000
3	1416	2	232000
4	3000	4	539900

这时我们需要进行均值归一化让数据统一量级，让梯度下降算法更快的收敛。均值归一化就是将每个特征的值减去平均值再除以标准差。

data2 = (data2-data2.mean())/data2.std()
print(data2.head())

	Size	Bedrooms	Price
0	0.130010	-0.223675	0.475747
1	-0.504190	-0.223675	-0.084074
2	0.502476	-0.223675	0.228626
3	-0.735723	-1.537767	-0.867025
4	1.257476	1.090417	1.595389

数据处理

data2.insert(0, 'ones', 1)
x2 = data2.iloc[:, 0:-1]  
y2 = data2.iloc[:, -1:]
x2, y2 = np.matrix(x2.values), np.matrix(y2.values)
theta2 = np.matrix(np.zeros(x2.shape[1]))  

梯度下降

g2, cost2 = gradient_descent(x2, y2, theta2, alpha, inters)
print(g2)

输出结果为[[-1.10892383e-16 8.84042349e-01 -5.24551809e-02]]

特征方程

正规方程是通过求解下面的方程来找出使得代价函数最小的参数的：

$$\frac{\partial J(\Theta_j)}{\partial\Theta_j}=0$$

通过特征方程，解得：

$$\Theta=(X^TX)^{-1}X^TY$$

python实现代码如下：

def normal_eqn(x_, y_):  # 正规方程
    return np.linalg.inv(x_.T.dot(x_)).dot(x_.T).dot(y_)

梯度下降与正规方程的比较：

梯度下降：需要选择学习率α，需要多次迭代，当特征数量n大时也能较好适用，适用于各种类型的模型

正规方程：不需要选择学习率α，一次计算得出，需要计算逆矩阵，如果特征数量n较大则运算代价大。

下面使用特征方程求解θ与梯度下降的结果进行对比。

print(normal_eqn(x2, y2))

结果为如下，可以看出与梯度下降所得的结果很相近。

[[-9.36750677e-17]
[ 8.84765988e-01]
[-5.31788197e-02]]

最后给出预测函数图像：