机器学习之回归模型

基本形式

线性模型(linear model)就是试图经过属性的线性组合来进行预测的函数,基本形式以下:
python

f(x)=wTx+b

许多非线性模型可在线性模型的基础上经过引入层结构或者高维映射(好比核方法)来解决。线性模型有很好的解释性。

线性回归

线性回归要求均方偏差最小:
git

(w,b)=argmini=1m(f(xi)yi)2

均方偏差有很好的几何意义,它对应了经常使用的欧式距离(Euclidean distance)。基于均方偏差最小化来进行模型求解称为最小二乘法(least square method),线性回归中,最小二乘发就是试图找到一条直线,使得全部样本到直线的欧式距离之和最小。

咱们把上式写成矩阵的形式:
github

w=argmin(yXw)T(yXw)

这里咱们把b融合到w中,X中最后再加一列1。为了求最小值,咱们对w求导并令其为0:
2XT(Xwy)=0

XTX 为满秩矩阵(full-rank matrix)时是可逆的。此时:
w=(XTX)1XTy

xi=(xi,1) ,能够获得线性回归模型:
f(xi)=xTi(XTX)1XTy

最小二乘法的代码以下:web

def standRegres(xArr,yArr):
    xMat = mat(xArr); yMat = mat(yArr).T
    xTx = xMat.T*xMat
    if linalg.det(xTx) == 0.0:
        print "This matrix is singular, cannot do inverse"
        return
    ws = xTx.I * (xMat.T*yMat)
    return ws

然而,现实状况是不少时候, XTX 并不满秩,这样 w 的解可能有多个解,它们都能使得均方偏差最小化。碰到数据特征比样本数还要多的状况,一种方案是引入正则化(regularization)项。另外一种方案是能够“缩减系数”,好比用岭回归(ridge regression)、LASSO法(该方法效果好可是计算复杂),前向逐步回归(能够获得与LASSO差很少的效果,但更容易实现)、LAR、PCA回归等。算法

局部加权线性回归

咱们在介绍岭回归以前先看一下局部加权线性回归。这首针对线性回归容易出现“欠拟合”现象提出的,由于它求的是具备最小均方偏差的无偏估计。若是模型欠拟合将不能取得最好的效果。因此有些方法在估计中引入一些误差,从而下降均方偏差。其中一个方法就是局部加权线性回归(Locally Wegihted Linear Regression,LWLR),咱们给待预测点附近的每一个点赋予必定的权重。回归系数 w 形式以下:
dom

w=(XTWX)1XTWy

w 就是用来给每一个数据点赋权重的矩阵。
LWLWR使用核方法来赋权重,最经常使用的就是高斯核:
w(i,i)=exp(|xix|2k2)

这样构建了一个只含对角元素的权重矩阵 w ,而且点x与x(i)越近,w(i,i)将会越大。k决定了对附近的点赋多大权重。
函数代码以下:

def lwlr(testPoint,xArr,yArr,k=1.0):
    xMat = mat(xArr); yMat = mat(yArr).T
    m = shape(xMat)[0]
    weights = mat(eye((m)))
    for j in range(m): #next 2 lines create weights matrix
        diffMat = testPoint - xMat[j,:] 
        weights[j,j] = exp(diffMat*diffMat.T/(-2.0*k**2))
    xTx = xMat.T * (weights * xMat)
    if linalg.det(xTx) == 0.0:
        print "This matrix is singular, cannot do inverse"
        return
    ws = xTx.I * (xMat.T * (weights * yMat))
    return testPoint * ws

文末给出的下载连接中regreesion.py文件中,有个线性加权回归例子。机器学习

岭回归

就是在 XTX 加上一个 λI 使得矩阵非奇异,这样 XTX+λI 就能够求逆了。因而:
ide

w=(XTX+λI)1XTy

岭回归最早用来处理特征数多于样本数的状况,如今也用于在估计中加入误差,从而获得更好的估计。这里经过引入 λ 来限制全部 w 之和,经过引入该惩罚项,可以减小不重要的参数。这种方法叫作“缩减”(shrinkage)。咱们经过预测偏差最小化来获得 λ ,即经过选取不一样的 λ 来重复训练测试过程,最终获得一个使预测偏差最小的 λ
相关代码以下:

def ridgeRegres(xMat,yMat,lam=0.2):
    xTx = xMat.T*xMat
    denom = xTx + eye(shape(xMat)[1])*lam
    if linalg.det(denom) == 0.0:
        print "This matrix is singular, cannot do inverse"
        return
    ws = denom.I * (xMat.T*yMat)
    return ws

def ridgeTest(xArr,yArr):
    xMat = mat(xArr); yMat=mat(yArr).T
    yMean = mean(yMat,0)
    yMat = yMat - yMean     #to eliminate X0 take mean off of Y
    #regularize X's
    xMeans = mean(xMat,0)   #calc mean then subtract it off
    xVar = var(xMat,0)      #calc variance of Xi then divide by it
    xMat = (xMat - xMeans)/xVar
    numTestPts = 30
    wMat = zeros((numTestPts,shape(xMat)[1]))
    for i in range(numTestPts):
        ws = ridgeRegres(xMat,yMat,exp(i-10))
        wMat[i,:]=ws.T
    return wMat

ridgeRegres用于计算回归系数,ridgeTest用于在一组 λ 上测试结果。
一样在使用特征时,都要先进行归一化处理。svg

LASSO

LASSO是least absolute shrinkage and selection operator的简称。咱们注意到增长以下约束,最小二乘法回归会获得与岭回归同样的公式:
函数

k=1nw2kλ

上式限定了全部回归系数的平方和不大于 λ 。在用最小二乘法回归时,当两个或更多特征相关时,可能会获得一个很大的正系数和一个很大的负系数,正由于上面约束,岭回归避免了这个问题。
LASSO的约束为:

k=1n|wk|λ

虽然用绝对值取代了平方和,但结果却差异很大。在 λ 足够小的时候,一些系数会所以被迫缩减为0,这样结果就有很好的稀疏性,这个特性能够帮助咱们更好地理解数据。两个看起来相差无几,但细微的变化却极大增长了计算复杂度。为了在此约束下求解,须要使用二次规划算法。

前向逐步回归

前向逐步回归能够获得与LASSO差很少的结果,但更加简单,是一种贪心算法,每一步都尽量减小偏差。一开始,全部权重设为1,而后每一步所作的决策是对某个权重增长或减小一个很小的值。
代码以下:

def stageWise(xArr,yArr,eps=0.01,numIt=100):
    xMat = mat(xArr); yMat=mat(yArr).T
    yMean = mean(yMat,0)
    yMat = yMat - yMean     #can also regularize ys but will get smaller coef
    xMat = regularize(xMat)
    m,n=shape(xMat)
    returnMat = zeros((numIt,n)) #testing code remove
    ws = zeros((n,1)); wsTest = ws.copy(); wsMax = ws.copy()
    for i in range(numIt):#could change this to while loop
        #print ws.T
        lowestError = inf; 
        for j in range(n):
            for sign in [-1,1]:
                wsTest = ws.copy()
                wsTest[j] += eps*sign
                yTest = xMat*wsTest
                rssE = rssError(yMat.A,yTest.A)
                if rssE < lowestError:
                    lowestError = rssE
                    wsMax = wsTest
        ws = wsMax.copy()
        returnMat[i,:]=ws.T
    return returnMat

逐步线性回归主要优势在于它能够帮助人们理解现有模型并作出改进。当构建一个模型后,运行该算法找出重要的特征,这样就能够及时中止对那些不重要特征的收集。当应用缩减方法时,模型增长了误差(bias),也减小了方差。

对数概率回归(逻辑回归)

上面是回归任务,若是咱们须要的是分类任务呢?咱们只须要找一个单调可微函数将分类任务的真实标记与线性回归的预测值就能够了。

二分类最理想的是“单位阶跃函数”(unit-step function)。可是此函数不连续,因而咱们但愿找一个单调可微的替代函数:

y=11+ez

logistic function是一种Sigmoid函数,咱们将z代入线性回归模型:
y=11+ewTx+b

变换一下:
lny1y=wTx+b

若是y是正样本的可能性,那么1-y就是负样本的可能性,二者比值称为概率,左边的意义就是对数概率。
实际上上式是在用线性回归的模型预测结果去逼近真实标记的对数概率。虽然名字叫“回归”,倒是一种分类方法。这种方法优势不少,它直接对分类可能性进行建模,无需事先假设数据分布。此外,对率函数是任意阶可导的凸函数,利用现有的数值优化算法能够直接获得最优解。
下面来看如下具体推导:
假设y=1是二分类中正样本,
lnp(y=1|x)p(y=0|x)=wTx+b

p(y=1|x)=ewTx+b1+ewTx+b

p(y=0|x)=11+ewTx+b

因而,咱们采用极大似然来估计w和b:
l(w,b)=i=1mlnp(yi|xi;w,b)

即令每一个样本属于其真实标记的几率越大越好。
咱们能够再转换成最小值问题,用经典的数值优化算法如梯度降低法(gradient descent method),牛顿法(Newton method)解决。

在此处咱们简便起见,直接采用梯度上升算法。咱们知道梯度算子老是指向函数值最快的方向。咱们每次向增加最快的方向移动一个值,称为步长,记为 a ,因此梯度上升算法的公式为:

w=w+af(w)w

该公式将一直迭代执行,直到达到某个中止条件(好比迭代次数)为止。
实现代码为:(固然以前数据先标准化)

def gradAscent(dataMatIn, classLabels):
    dataMatrix = mat(dataMatIn)             
    labelMat = mat(classLabels).transpose() 
    m,n = shape(dataMatrix)
    alpha = 0.001
    maxCycles = 500
    weights = ones((n,1))
    for k in range(maxCycles):              
        h = sigmoid(dataMatrix*weights)
        error = (labelMat - h)              #vector subtraction
        weights = weights + alpha * dataMatrix.transpose()* error
    return weights

代码中求导很差操做,这里咱们用差分能够起到一样的效果。
梯度上升算法在每次更新回归系数时要遍历整个数据集,改进方法是一次仅用一个样本点来更新回归系数,称为“随机梯度上升算法”,这是一种在线算法。上述代码仍有改进之处,好比 a 在每次迭代时候都会调整,这会缓解数据波动或者高频波动。 a 会随着获得次数不断减少,但永远不会减到0,这样屡次迭代以后新数据仍具备必定的影响。
改进后的随机梯度上升代码以下:

def stocGradAscent(dataMatrix, classLabels, numIter=150):
    m,n = shape(dataMatrix)
    weights = ones(n)   #initialize to all ones
    for j in range(numIter):
        dataIndex = range(m)
        for i in range(m):
            alpha = 4/(1.0+j+i)+0.0001    #apha decreases with iteration, does not 
            randIndex = int(random.uniform(0,len(dataIndex)))#go to 0 because of the constant
            h = sigmoid(sum(dataMatrix[randIndex]*weights))
            error = classLabels[randIndex] - h
            weights = weights + alpha * error * dataMatrix[randIndex]
            del(dataIndex[randIndex])
    return weights

这里用个例子来使用Logistic回归来预测患有疝病的马存活问题。可是咱们发现了样本数据缺失问题,一般的解决方法有:

  • 使用可用特征均值填补
  • 使用特殊值来填补缺失值,如-1
  • 忽略有缺失值的样本
  • 使用类似样本的均值填补
  • 使用另外的机器学习算法预测缺失值

这里咱们所有用0来代替缺失值,这样在更新时不会影响系数。即特征对应0时,系数将保持0.5。咱们在数据中还发现了标签缺失,这个在强化学习中能够经过更靠近某些标签就将其标记为那个标签,在本例中因为此类较少,咱们直接丢弃。

线性判别分析

(Linear Discriminant Analysis,LDA)与LDA(Latent Dirichlet Allocation)算法不要搞混了哦,虽然简写都是同样的。
LDA的思想很简单:设法将样例投影到一条直线上,使得同类样例的投影点近可能接近,异类样例的投影尽量原理,在对新样本进行分类时,将其投影到这条直线上,再根据投影点位置来肯定分类。

多分类学习

二分类问题推广到多分类主要用拆解法,主要策略有:
- 一对一
- 一对多
- 多对多
这样那个容易碰到类别不平衡问题,就是某类的数据量特别大或者特别少。

最后就是本文数据集和代码的下载地址啦,请点击这里