回归

回归分析用来建立方程模拟自变量和应变量之间的关系。
一元线性回归包含一个自变量和一个应变量
如果包含两个以上的自变量，则称为多元回归分析。

一元线性回归

目标就是从1，2，3……中找到最适合的一条直线，能让所有点到直线的垂直距离最短即可。这里可以采用最小二乘法：

然后对w,b两个参数进行求导，然后进行梯度更新即可。

下面采用两种方式给出代码，一种是自己写的python代码，一种是调用sklearn库

调用sklearn库

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

# 自己设定数据集
# 斜率为2，截距为4，y值服从0-20的均匀分布
x = np.array(list(range(1,100)))
y = np.array([2 * item + np.random.uniform(0,20) + 4 for item in x])
# 展示一下图形
plt.scatter(x, y)
plt.show()

model = LinearRegression()
# sklearn格式要求二维数据
x = x.reshape(-1,1)
u = y.reshape(-1,1)
model.fit(x, y)
plt.plot(x, y, 'g.') # 这里g表示绿色，.表示散点的形态
plt.plot(x, model.predict(x), 'r') # r表示红颜色

纯python写法

容易出现梯度爆炸的情况

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

# 自己设定数据集
# 斜率为2，截距为7，y值服从均匀分布
x=np.random.normal(0,10,100)
noise=np.random.normal(0,0.05,100)
y=2*x+7+noise
# 展示一下图形
plt.scatter(x, y)
plt.show()

# 计算损失
def count_loss(x, y, w, b):
    # m表示样本点的数量
    m = len(x)
    return np.sum(np.square((w * x + b)-y)) / (2.0*m)

w = 0
b = 0
epochs = 1000
lr = 0.01
for epoch in range(epochs):
    w_grad = b_grad = 0
    loss = count_loss(x, y, w, b)
    print(loss)
    m = float(len(x))
    w_grad = np.sum((w * x + b -y) * x) / m
    b_grad = np.sum((w * x + b -y) ) / m
    w = w  - lr*w_grad
    b = b  - lr*b_grad

print(w, b)

多元线性回归

在很多的情况下，变量不是唯一的。例如当考虑一栋房子的价值时，户型、学区、交通、单价都能成为评判的指标：

做法和一元线性回归是一样的，只是参数变多了，同时梯度下降的参数变多了。

代码如下：

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

# 准备好了自定义的数据集，共有2个特征，分别是x1, x2，共10组数据
x1 = [100, 50, 100, 100, 50, 80, 75, 65, 90, 90]
x2 = [4, 3, 4, 2, 2, 2, 3, 4, 3, 2]
y = [9.3, 4.8, 8.9, 6.5, 4.2, 6.2, 7.4, 6., 7.6, 6.1]

model = LinearRegression()
x = np.zeros((10, 2))
x[:, 0] = x1
x[:, 1] = x2
model.fit(x, y)

# 分别打印系数，截距
print(model.coef_, model.intercept_)

ax = plt.figure().add_subplot(111, projection='3d') # 创建了一个画布
ax.scatter(x1, x2, y, c='r', marker='o', s=100)
# 生成网格矩阵
x1, x2 = np.meshgrid(x1,x2)
'''
x1, x2 = np.meshgrid([1,2,3], [4,5,6])
print(x1)
print(x2)
可以生成网格矩阵：
(1,4)(2,4)(3,4)
(1,5)(2,5)(3,5)
(1,6)(2,6)(3,6)
可以用以下代码进行展示
x1, x2 = np.meshgrid([1,2,3], [4,5,6])
plt.scatter(x1,x2)
'''
z = model.intercept_ + x1 * model.coef_[0] + x2 * model.coef_[1]
ax.plot_surface(x1, x2, z)
plt.show()

多项式回归（非线性）

假设我们不是要找直（或者超平面），而是需要找到一个多项式所表示的曲线（或者超曲面），例如二次曲线$y=a{t}^2+bt+c$

import numpy as np
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import PolynomialFeatures
x = [1,2,3,4,5,6,7,8,9,10]
y = [45000, 50000, 60000, 80000, 110000, 150000, 200000, 300000, 500000, 1000000]

# 生成多项式特征
'''
x = np.array([[1],[2],[3],[4]])
poly = PolynomialFeatures(2)
poly.fit_transform(x)
结果：
这里的1是偏置项
array([[ 1.,  1.,  1.],
       [ 1.,  2.,  4.],
       [ 1.,  3.,  9.],
       [ 1.,  4., 16.]])
'''
x = np.array(x).reshape((-1,1))
y = np.array(y).reshape((-1,1))
# feature的个数越大，拟合程度越高
poly = PolynomialFeatures(5)
x_trans = poly.fit_transform(x)
print(x_trans)
model = LinearRegression()
model.fit(x_trans, y)

# 图像预测及展示
plt.scatter(x, y)

plt.plot(x, model.predict(x_trans))
plt.legend()
plt.show()

标准方程法

这种情况下，若矩阵不可逆就没有办法求出：

梯度下降法VS标准方程法

岭回归

因为根据矩阵秩的特点：r(AB) <= min(r(A), r(B))，所以一定不满秩。发现加入正则项后，求导的时候，真好可以凑成上述的公式，因此可以避免不满秩。

这里以longley数据集为例：
其2-7列是特征，第一列是结果。

链接：https://pan.baidu.com/s/1uT620nZtGVayteIoxiRpRA
提取码：0525

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from sklearn import linear_model
# 读取2-7列作为特征输入，第1列作为结果输出
x_data = pd.read_csv(r'longley.csv', usecols=range(2,8))
y_data = pd.read_csv(r'longley.csv', usecols=[1])

# 定义了正则项前的系数，用alpha表示。生成0.001-1之间的50个样本测试点
alpha_to_test = np.linspace(0.001, 1, 50)
# Ridge表示岭回归，CV表示交叉验证
model = linear_model.RidgeCV(alphas=alpha_to_test, store_cv_values=True)
model.fit(x_data ,y_data)

# 最好的岭系数
print(model.alpha_)
# 交叉验证的损失，因为有16个样本，每次训练都会取一个样本作为测试集，又因为初始化设置的岭系数有50，因为有50个损失
print(model.cv_values_.shape)

# 先删除一个维度，然后在0维度上求均值，算出每一个岭系数的loss值
# model.alpha_表示最好的岭系数，min(model.cv_values_.squeeze(1).mean(axis=0))表示最小的loss值。
plt.plot(alpha_to_test, model.cv_values_.squeeze(1).mean(axis=0))
plt.scatter(model.alpha_, min(model.cv_values_.squeeze(1).mean(axis=0)), color='r')
plt.show()

LASSO回归

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from sklearn import linear_model
# 读取2-7列作为特征输入，第1列作为结果输出
x_data = pd.read_csv(r'longley.csv', usecols=range(2,8))
y_data = pd.read_csv(r'longley.csv', usecols=[1])

# 定义了正则项前的系数，用alpha表示。生成0.001-1之间的50个样本测试点
alpha_to_test = np.linspace(0.001, 1, 50)
# Ridge表示岭回归，CV表示交叉验证
model = linear_model.LassoCV()
model.fit(x_data ,y_data)

# 最好的岭系数
print(model.alpha_)
# 交叉验证的损失，因为有16个样本，每次训练都会取一个样本作为测试集，又因为初始化设置的岭系数有50，因为有50个损失
print(model.coef_)

弹性网

弹性网是对岭回归和LASSO回归的一种中和的算法。
LASSO回归和弹性网都能使得某些特征的系数为0，即可以排除一些共线的特征（成比例的特征）

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from sklearn import linear_model
# 读取2-7列作为特征输入，第1列作为结果输出
x_data = pd.read_csv(r'longley.csv', usecols=range(2,8))
y_data = pd.read_csv(r'longley.csv', usecols=[1])

# Ridge表示岭回归，CV表示交叉验证
model = linear_model.ElasticNetCV()
model.fit(x_data ,y_data)

# 最好的岭系数
print(model.alpha_)
# 交叉验证的损失，因为有16个样本，每次训练都会取一个样本作为测试集，又因为初始化设置的岭系数有50，因为有50个损失
print(model.coef_)