普通线性模型:
假设主要有以下几点:
1.响应变量Y和误差正态性:响应变量Y和误差项服从正态分布,并且误差是一个白噪声过程,因而具有零均值,同方差的特性。
2.预测量X和未知参数的非随机性:预测量X具有非随机性,可测不存在测量误差;未知参数认为是未知但不具有随机性的常数。
3.主要研究响应变量的均值E[Y]
4.连接方式:响应变量E[Y]与预测变量的线性组合。f(x)=x
广义线性模型是在普通线性模型的基础上,将上述四点模型假设进行推广而得出的应用范围更广,更实用的回归模型。
1.响应变量的分布推广至指数分布族:比如正态分布,泊松分布,二项分布,负二项分布,伽马分布,逆高斯分布
2.预测变量X和未知参数的非随机性
3.研究对象:广义线性模型的主要研究对象任然是响应变量的均值E[Y]
4.连接方式;理论上可以是任意的,而不再局限于f(x)=x.
因变量y不再只是正态分布,而是扩大为指数族中的任一分布。
解释变量x的线性组合不再直接用于解释因变量y的均值u,而是通过一个联系函数g来解释g(u).
线性回归是一种parametric learning algorithm,而局部加权线性回归是一种non-parametric learning algorithm。Parametric learning algorithm有固定的(指的是:值的大小是固定的)、有限的参数、通过训练样本,找到合适的参数后,对于之后未知的输入,我们可以直接利用这组参数得出其相应的预测输出。而non-parametric learning algorithm需要的计算量与输入的训练集大小成正比,对于每次新的输入,需要计算相应的参数后,才能求出相应的预测输出。
局部加权回归的缺点是:需要拿着训练数据去做预测。
当训练样本的features大于样本数时,那么矩阵X就不是满秩的。无法计算X.T*X。为了解决这个问题,引入了ridge regression的概念。ridge regession加入了lambda*I到矩阵X.T*X
当两个或更多的特征有关联时,如果用正规的最小二乘法回归,可能会有一个非常大的正的权值和一个负的非常大的权值。
自适应学习率的学习:
基于Armijo准则计算搜索方向上的最大步长,其基本思想是沿着搜索方向移动一个比较大的步长估计值,然后以迭代形式不断缩减步长,直到该步长使得函数值f(xk+a*dk)相对于当前函数值f(xk)的减小程度大于预设的期望值为止。
# -*- coding:utf8 -*-
import math
import matplotlib.pyplot as plt
def f(w, x):
N = len(w)
i = 0
y = 0
while i < N - 1:
y += w[i] * x[i]
i += 1
y += w[N - 1] # 常数项
return y
def gradient(data, w, j):
M = len(data) # 样本数
N = len(data[0])
i = 0
g = 0 # 当前维度的梯度
while i < M:
y = f(w, data[i])
if (j != N - 1):
g += (data[i][N - 1] - y) * data[i][j]
else:
g += data[i][N - 1] - y
i += 1
return g / M
def gradientStochastic(data, w, j):
N = len(data) # 维度
y = data[N - 1] - f(w, data)
if (j != N - 1):
return y * data[j]
return y # 常数项
def isSame(a, b):
n = len(a)
i = 0
while i < n:
if abs(a[i] - b[i]) > 0.01:
return False
i += 1
return True
def fw(w, data):
M = len(data) # 样本数
N = len(data[0])
i = 0
s = 0
while i < M:
y = data[i][N - 1] - f(w, data[i])
s += y ** 2
i += 1
return s / 2
def fwStochastic(w, data):
y = data[len(data) - 1] - f(w, data)
y **= 2
return y / 2
def numberProduct(n, vec, w):
N = len(vec)
i = 0
while i < N:
w[i] += vec[i] * n
i += 1
def assign(a):
L = []
for x in a:
L.append(x)
return L
# a = b
def assign2(a, b):
i = 0
while i < len(a):
a[i] = b[i]
i += 1
def dotProduct(a, b):
N = len(a)
i = 0
dp = 0
while i < N:
dp += a[i] * b[i]
i += 1
return dp
# w当前值;g当前梯度方向;a当前学习率;data数据
def calcAlpha(w, g, a, data):
c1 = 0.3
now = fw(w, data)
wNext = assign(w)
numberProduct(a, g, wNext)
next = fw(wNext, data)
# 寻找足够大的a,使得h(a)>0
count = 30
while next < now:
a *= 2
wNext = assign(w)
numberProduct(a, g, wNext)
next = fw(wNext, data)
count -= 1
if count == 0:
break
# 寻找合适的学习率a
count = 50
while next > now - c1*a*dotProduct(g, g):
a /= 2
wNext = assign(w)
numberProduct(a, g, wNext)
next = fw(wNext, data)
count -= 1
if count == 0:
break
return a
# w当前值;g当前梯度方向;a当前学习率;data数据
def calcAlphaStochastic(w, g, a, data):
c1 = 0.01 # 因为是每个样本都下降,所以参数运行度大些,即:激进一些
now = fwStochastic(w, data)
wNext = assign(w)
numberProduct(a, g, wNext)
next = fwStochastic(wNext, data)
# 寻找足够大的a,使得h(a)>0
count = 30
while next < now:
if a < 1e-10:
a = 0.01
else:
a *= 2
wNext = assign(w)
numberProduct(a, g, wNext)
next = fwStochastic(wNext, data)
count -= 1
if count == 0:
break
# 寻找合适的学习率a
count = 50
while next > now - c1*a*dotProduct(g, g):
a /= 2
wNext = assign(w)
numberProduct(a, g, wNext)
next = fwStochastic(wNext, data)
count -= 1
if count == 0:
break
return a
def normalize(g):
s = 0
for x in g:
s += x * x
s = math.sqrt(s)
i = 0
N = len(g)
while i < N:
g[i] /= s
i += 1
def calcCoefficient(data, listA, listW, listLostFunction):
M = len(data) # 样本数目
N = len(data[0]) # 维度
w = [0 for i in range(N)]
wNew = [0 for i in range(N)]
g = [0 for i in range(N)]
times = 0
alpha = 100.0 # 学习率随意初始化
same = False
while times < 10000:
i = 0
while i < M:
j = 0
while j < N:
g[j] = gradientStochastic(data[i], w, j)
j += 1
normalize(g) # 正则化梯度
alpha = calcAlphaStochastic(w, g, alpha, data[i])
#alpha = 0.01
numberProduct(alpha, g, wNew)
print "times,i, alpha,fw,w,g:\t", times, i, alpha, fw(w, data), w, g
if isSame(w, wNew):
if times > 5: #防止训练次数过少
same = True
break
assign2(w, wNew) # 更新权值
listA.append(alpha)
listW.append(assign(w))
listLostFunction.append(fw(w, data))
i += 1
if same:
break
times += 1
return w
if __name__ == "__main__":
fileData = open("d8.txt")
data = []
for line in fileData:
d = map(float, line.split('\t'))
data.append(d)
fileData.close()
listA = [] # 每一步的学习率
listW = [] # 每一步的权值
listLostFunction = [] # 每一步的损失函数值
w = calcCoefficient(data, listA, listW, listLostFunction)
# 绘制学习率
plt.plot(listA, 'r-', linewidth=2)
plt.plot(listA, 'go')
plt.xlabel('Times')
plt.ylabel('Ratio/Step')
plt.show()
# 绘制损失
listLostFunction = listLostFunction[0:100]
listLostFunction[0] /= 2
plt.plot(listLostFunction, 'r-', linewidth=2)
plt.plot(listLostFunction, 'gv', alpha = 0.75)
plt.xlabel('Times')
plt.ylabel('Loss Value')
plt.grid(True)
plt.show()
# 绘制权值
X = []
Y = []
for d in data:
X.append(d[0])
Y.append(d[1])
plt.plot(X, Y, 'cp', label=u'Original Data', alpha=0.75)
x = [min(X), max(X)]
y = [w[0] * x[0] + w[1], w[0] * x[1] + w[1]]
plt.plot(x, y, 'r-', linewidth=3, label='Regression Curve')
plt.legend(loc='upper left')
plt.grid(True)
plt.show()
以上代码为:分别用了梯度下降算法和随机梯度下降算法求线性回归的方程。学习率使用Armijo准则计算。
采用多项式插值法拟合简单函数,然后根据简单函数估计函数的极值点,这样选择合适的步长的效率会高很多。
现在拥有的数据为:Xk的函数值f(Xk)及其导数f'(Xk),再加上第一次尝试的步长a0。如果a0满足条件,显然算法退出;若a0不满足条件,则可以根据上述信息可以构造一个二次函数:
显然,导数为0处的最优值为:
若a1满足Armijo准则,则输出该学习率;否则继续迭代。
#include "stdafx.h"
#include <math.h>
double f(double x) {
return x * x;
}
double g(double x) {//一阶导数
return 2 * x;
}
double g2(double x) {//二阶导数
return 2;
}
//x为当前值 d为x处的导数 a为输入的学习率
//使用Armijo求学习率
double GetA_Armijo(double x, double d, double a) {
double c1 = 0.3;
double now = f(x);
double next = f(x-a*d);
int count = 30;
while (next < now) {
a *= 2;
next = f(x-a*d);
count --;
if (count == 0)
break;
}
count = 50;
while(next > now - c1*a*d*d) {
a /= 2;
next = f(x-a*d);
count --;
if(count == 0) {
break;
}
}
return a;
}
//使用二次插值法求学习率
// x为当前值 d为梯度的方向 a为前一步的学习率
double GetA_Quad(double x, double d, double a) {
double c1 = 0.03;
double now = f(x);
double next = f(x-a*d);
int count = 30;
while (next < now) {
a *= 2;
next = f(x-a*d);
count -= 1;
if (count == 0) {
break;
}
}
while (next < now-c1*a*d*d) {
a = d*a*a/(now + d*a-next);//二次插值求学习率
}
}
int _tmain(int argc, _TCHAR* argv[])
{
double x = 1.5;
double d; //一阶导 梯度方向
double a=0.01; //学习率
for (int i = 0; i < 100; i++)
{
d=g(x);
//a = GetA_Armijo(x,d,10);
}
return 0;
}
逻辑回归和广义线性模型有什么关系呢?
Softmax Regression: