简化表示为
代数法表示：
矩阵表示：
注：
m为样本个数，n为特征个数加一
设定初始参数,不断迭代，使得
最小化，θ的迭代公式如下：
当J为凸函数时，梯度下降法相当于让参数\(\theta\)不断向J的最小值位置移动
梯度下降法的缺陷：如果函数为非凸函数，有可能找到的并非全局最优值，而是局部最优值。
令：
利用矩阵的链式求导法则，和两个矩阵的求导公式：
整理得到：
- 最小二乘法需要计算\(X^TX\)的逆矩阵，有可能它的逆矩阵不存在，这样就没有办法直接用最小二乘法了，此时梯度下降法仍然可以使用。当然，我们可以通过对样本数据进行整理，去掉冗余特征。让\(X^TX\)的行列式不为0，然后继续使用最小二乘法。
- 当样本特征n非常的大的时候，计算\(X^TX\)的逆矩阵是一个非常耗时的工作（nxn的矩阵求逆），甚至不可行。此时以梯度下降为代表的迭代法仍然可以使用。那这个n到底多大就不适合最小二乘法呢？如果你没有很多的分布式大数据计算资源，建议超过10000个特征就用迭代法吧。或者通过主成分分析降低特征的维度后再用最小二乘法。
- 如果拟合函数不是线性的，这时无法使用最小二乘法，需要通过一些技巧转化为线性才能使用，此时梯度下降仍然可以用
- 特殊情况下，当样本量m很少，小于特征数n的时候，这时拟合方程是欠定的，常用的优化方法都无法去拟合数据。当样本量m等于特征数n的时候，用方程组求解就可以了。当m大于n时，拟合方程是超定的，也就是我们常用与最小二乘法的场景了。
将\(J(\theta)\)泰勒展开到二阶：
两边求导，得：
重复迭代式可以利用Hessian矩阵：
牛顿法的收敛速度非常快，但海森矩阵的计算较为复杂，尤其当参数的维度很多时，会耗费大量计算成本。可以用其他矩阵替代海森矩阵，用拟牛顿法进行估计
拟牛顿法的思路是用一个矩阵替代计算复杂的海森矩阵H，因此要找到符合H性质的矩阵。
要求得海森矩阵符合的条件，同样对泰勒公式求导\(f‘(x) = f‘(x_0) + f‘‘(x_0)x -f‘‘(x_0)x_0\)
令\(x = x_1\)，即迭代后的值，代入可得：
更一般的
x_k为第k个迭代值
即找到矩阵G，使得它符合上式。 常用的拟牛顿法的算法包括DFP，BFGS等。
回到开始的线性模型
, 如果这里不仅仅是x的一次方，比如增加二次方，那么模型就变成了多项式回归。这里写一个只有两个特征的p次方多项式回归的模型：
令\(*x*_0=1,x_1=x_1,x_2=x_2,x_3=x^2_1,x_4=x^2_2,x_5=x_1x_2\)
这样我们就得到了下式：
可以发现，又重新回到了线性回归，这是一个五元线性回归，可以用线性回归的方法来完成算法。对于每个二元样本特征(x1,x2), 得到一个五元样本特征(1,x1,x2,x21,x22,x1x2)，通过这个改进的五元样本特征，重新把不是线性回归的函数变回线性回归。
输出Y不满足和X的线性关系，但是lnY 和X满足线性关系，模型函数如下：
这样对与每个样本的输入y，我们用 lny去对应， 从而仍然可以用线性回归的算法去处理这个问题。
Lasso 回归
Lasso回归可以使得一些特征的系数变小，甚至还是一些绝对值较小的系数直接变为0。增强模型的泛化能力。Lasso回归的求解办法一般有坐标轴下降法（coordinate descent）和最小角回归法（ Least Angle Regression）。
Ridge回归
Ridge回归在不抛弃任何一个特征的情况下，缩小了回归系数，使得模型相对而言比较的稳定，但和Lasso回归比，这会使得模型的特征留的特别多，模型解释性差。
Ridge回归的求解比较简单，一般用最小二乘法。这里给出用最小二乘法的矩阵推导形式，和普通线性回归类似。
#生成数据
import numpy as np
#生成随机数
np.random.seed(1234)
x = np.random.rand(500,3)
#构建映射关系，模拟真实的数据待预测值,映射关系为y = 4.2 + 5.7*x1 + 10.8*x2
y = x.dot(np.array([4.2,5.7,10.8]))
import numpy as np
from sklearn.linear_model import LinearRegression
import matplotlib.pyplot as plt
%matplotlib inline
# 调用模型
lr = LinearRegression(fit_intercept=True)
# 训练模型
lr.fit(x,y)
print("估计的参数值为：%s" %(lr.coef_))
# 计算R平方(拟合优度，取值[0,1], 越接近1拟合越好)
print(‘R2:%s‘ %(lr.score(x,y)))
# 任意设定变量，预测目标值
x_test = np.array([2,4,5]).reshape(1,-1)
y_hat = lr.predict(x_test)
print("预测值为: %s" %(y_hat))
class LR_GD():
def __init__(self):
self.w = None
def fit(self,X,y,alpha=0.02,loss = 1e-10): # 设定步长为0.002,判断是否收敛的条件为1e-10
y = y.reshape(-1,1) #重塑y值的维度以便矩阵运算
[m,d] = np.shape(X) #自变量的维度
self.w = np.zeros((d)) #将参数的初始值定为0
tol = 1e5
while tol > loss:
h_f = X.dot(self.w).reshape(-1,1)
theta = self.w + alpha*np.mean(X*(y - h_f),axis=0) #计算迭代的参数值
tol = np.sum(np.abs(theta - self.w))
self.w = theta
def predict(self, X):
# 用已经拟合的参数值预测新自变量
y_pred = X.dot(self.w)
return y_pred
if __name__ == "__main__":
lr_gd = LR_GD()
lr_gd.fit(x,y)
print("估计的参数值为：%s" %(lr_gd.w))
x_test = np.array([2,4,5]).reshape(1,-1)
print("预测值为：%s" %(lr_gd.predict(x_test)))
class LR_LS:
def __init__(self, fit_intercept=True):
self.w = None
self.fit_intercept = fit_intercept # bias 截距
def fit(self, X, y):
if self.fit_intercept:
X = np.c_[np.ones(X.shape[0]), X]
pseudo_inverse = np.dot(np.linalg.inv(np.dot(X.T, X)), X.T)
self.w = np.dot(pseudo_inverse, y)
def predict(self, X):
if self.fit_intercept:
X = np.c_[np.ones(X.shape[0]), X]
return np.dot(X, self.w)
def first_derivativ(X, y, w):
return np.matmul(X.T, np.matmul(X, w) - y)
def second_derivative(X):
return np.matmul(X.T, X)
def get_error(X, y, w):
tmp = np.matmul(X, w) - y
return np.matmul(tmp.T, tmp) / 2
def get_min_m(X, y, sigma, delta, d, w, g):
m = 0
while True:
w_new = w + pow(sigma, m) * d
left = get_error(X, y , w_new)
right = get_error(X, y , w) + delta * pow(sigma, m) * g.T * d
if left <= right:
break
else:
m += 1
return m
class LR_NEWTON():
def __init__(self):
self.w = None
def fit(self, X, y, max_iter=50, sigma=0.5, delta=0.25):
‘‘‘
sigma : (0, 1)
delta : (0. 0.5)
‘‘‘
print(type(X))
n = X.shape[1]
self.w = np.mat(np.ones((n, 1)))
it = 0
while it <= max_iter:
g = first_derivativ(X, y, self.w)
G = second_derivative(X)
try:
d = -np.linalg.inv(G) * g
except Exception:
print(‘不可逆‘)
break
m = get_min_m(X, y, sigma, delta, d, self.w, g) # 得到最小的m
self.w = self.w + pow(sigma, m) * d
if it % 10 == 9:
print ("\t---- itration: ", it, " , error: ", get_error(X, y , self.w)[0, 0])
it += 1
def predict(self, X):
y_pred = X.dot(self.w)
return y_pred
if __name__ == "__main__":
lr_nt = LR_NEWTON()
lr_nt.fit(x,y.reshape(-1,1))
print("估计的参数值为：%s" %(lr_nt.w))
x_test = np.array([2,4,5,2,4,5,2,4,5]).reshape(1,-1)
print("预测值为：%s" %(lr_nt.predict(x_test)))
class RidgeRegression:
def __init__(self, fit_intercept=True):
self.w = None
self.fit_intercept = fit_intercept
def fit(self, X, y, alpha=1):
if self.fit_intercept:
X = np.c_[np.ones(X.shape[0]), X]
A = alpha * np.eye(X.shape[1])
# A = 0
pseudo_inverse = np.dot(np.linalg.inv(np.dot(X.T, X) + A), X.T)
self.w = np.dot(pseudo_inverse, y)
def predict(self, X):
if self.fit_intercept:
X = np.c_[np.ones(X.shape[0]), X]
return np.dot(X, self.w)
if __name__ == "__main__":
lr_ridge = RidgeRegression()
lr_ridge.fit(x,y.reshape(-1,1))
print("估计的参数值为：%s" %(lr_ridge.w))
x_test = np.array([2,4,5,2,4,5,2,4,5]).reshape(1,-1)
print("预测值为：%s" %(lr_ridge.predict(x_test)))
原文：https://www.cnblogs.com/54hys/p/12741850.html