生活中我们会遇到很多分类问题,如一封邮件是否是垃圾邮件等。那么能否建立一种模型来给出垃圾邮件的概率?通过前面线性回归一文我们知道线性回归的结果是一个连续值,并且它的范围是无法限定的,而分类问题需要我们给出一个0.0~1.0之间的概率值,那么有没有什么办法可以将线性回归的结果值映射为一个0.0~1.0之间的概率值呢?有!它就是逻辑函数。
逻辑函数是一种S函数,函数形式为:
$$ p(t)=\frac{1}{1+e^{-t}} $$
逻辑函数的函数图像为(来源于维基百科):
利用逻辑函数我们可以构造二项逻辑回归模型(binomial logistic regression model)。二项逻辑回归模型中的二项意为二分类,也就是分类类别 $y\in{0,1}$ 。
设特征数量为 $n$ 即,特征为 $x$ ,特征的参数为 $\theta$ ,我们定义如下:
\begin{align} \\ x & = \begin{bmatrix} 1 \ x_1 \ \cdots \ x_n \end{bmatrix} \\ \\ \theta & = \begin{bmatrix} \theta_0 \\ \theta_1 \\ \vdots \\ \theta_n \end{bmatrix} \\ \\ x\theta & = \theta_0 + \theta_1x_1 + \theta_2x_2 + \cdots + \theta_nx_n \\ \\ h_\theta(x) & = \frac{1}{1+e^{x\theta}} \end{align}
$h_\theta(x)$ 就是我们的概率函数,我们可以假设它为分类类别为1的概率,那么类别1和类别0的概率分别为:
\begin{align} \\ P(y=1 \mid x) & =h_\theta(x) \\ \\ P(y=0 \mid x) & =1-h_\theta(x) \\ \end{align}
综合起来可以写成: \begin{align} \\ P(y \mid x) & =h_\theta(x)^y(1-h_\theta(x))^{1-y} \\ \end{align}
设训练样本数为 $m$,训练样本集为 $X$ ,训练输出集为 $Y$ ,如下: \begin{align} X & = \begin{bmatrix} x^{0} \\ x^{1} \\ \cdots \\ x^{m-1} \end{bmatrix} \\ \\ Y & = \begin{bmatrix} y^{0} \\ y^{1} \\ \vdots \\ y^{m-1} \end{bmatrix} \\ \end{align}
我们的目标是已知 $X$ 和 $Y$ 的情况下得到最优的 $\theta$。
哪个 $\theta$ 是最优的?我们需要先定义似然函数:
\begin{align} \\ L(\theta) &= \prod_{i=0}^{m-1} P(y^i \mid x^i) \\ \\ L(\theta) &= \prod_{i=0}^{m-1}h_\theta(x^i)^{y^i} (1-h_\theta(x^i))^{1-y^i} \\ \end{align}
我们在似然函数中引入自然对数以方便后续的求导,则:
\begin{align} \\ L(\theta) &= \log(\prod_{i=0}^{m-1}h_\theta(x^i)^{y^i} (1-h_\theta(x^i))^{1-y^i}) \\ L(\theta) &= \sum_{i=0}^{m-1}\log(h_\theta(x^i)^{y^i} (1-h_\theta(x^i))^{1-y^i}) \\ L(\theta) &= \sum_{i=0}^{m-1}(y^i\log(h_\theta(x^i)) + (1-y^i)\log(1-h_\theta(x^i))) \\ \end{align}
很明显似然函数最大值对应的 $\theta$ 就是我们求解的目标,所以问题变为: $$ \max_\theta L_\theta $$
使用梯度上升法可以帮助我们找到似然函数的最大值,参数 $\theta_j$的梯度为:
$$ \frac{\partial L(\theta)}{\partial \theta_j}=\frac{\partial}{\partial \theta_j} \sum_{i=0}^{m-1}(y^i\log(h_\theta(x^i)) + (1-y^i)\log(1-h_\theta(x^i))) $$
假设样本数为1,则得到如下: \begin{align} \\ \frac{\partial L(\theta)}{\partial \theta_j} &= \frac{\partial}{\partial \theta_j} (y\log(h_\theta(x)) + (1-y)\log(1-h_\theta(x))) \\ \\ \frac{\partial L(\theta)}{\partial \theta_j} &= \frac{y}{h_\theta(x)} \frac{\partial}{\partial \theta_j} h_\theta(x) + \frac{1-y}{1-h_\theta(x)}\frac{\partial}{\partial \theta_j}(1-h_\theta(x)) \\ \\ \frac{\partial L(\theta)}{\partial \theta_j} &= (\frac{y}{h_\theta(x)} - \frac{1-y}{1-h_\theta(x)})\frac{\partial}{\partial \theta_j} h_\theta(x) \\ \end{align}
又因为:
\begin{align} \\ \frac{\partial}{\partial \theta_j} h_\theta(x) &= \frac{\partial}{\partial \theta_j} \frac{1}{1+e^{x\theta}} \\ \\ \frac{\partial}{\partial \theta_j} h_\theta(x) &=\frac{-1}{(1+e^{x\theta})^2} \frac{\partial}{\partial \theta_j}e^{x\theta} \\ \\ \frac{\partial}{\partial \theta_j} h_\theta(x) &=\frac{-e^{x\theta}}{(1+e^{x\theta})^2} \frac{\partial}{\partial \theta_j} x\theta \\ \\ \frac{\partial}{\partial \theta_j} h_\theta(x) &=-h_\theta(x)(1-h_\theta(x)) \frac{\partial}{\partial \theta_j} x\theta \\ \\ \frac{\partial}{\partial \theta_j} h_\theta(x) &=-h_\theta(x)(1-h_\theta(x)) x_j \\ \end{align}
所以:
\begin{align} \\ \frac{\partial L(\theta)}{\partial \theta_j} &= (\frac{y}{h_\theta(x)} - \frac{1-y}{1-h_\theta(x)})(-h_\theta(x)(1-h_\theta(x)) x_j) \\ \\ \frac{\partial L(\theta)}{\partial \theta_j} &= -(y(1-h_\theta(x)) - (1-y)h_\theta(x))x_j \\ \\ \frac{\partial L(\theta)}{\partial \theta_j} &= -(y-h_\theta(x))x_j \\ \end{align}
多个样本的正确公式为: \begin{align} \\ \frac{\partial L(\theta)}{\partial \theta_j} &= -\sum_{i=0}^{m-1} (y^i-h_\theta(x^i))x_j^i \\ \end{align}
梯度上升法下的 $\theta_j$ 的更新公式为:($\alpha$ 为学习速度)
\begin{align} \\ \theta_j &:= \theta_j + \alpha \frac{\partial L(\theta)}{\partial \theta_j} \\ \\ \theta_j &:= \theta_j - \alpha \sum_{i=0}^{m-1} (y^i-h_\theta(x^i))x_j^i \\ \end{align}
关于随机梯度,批量梯度以及学习速度在上一文线性回归中我们已经介绍过,所以这边不再解释。
我们定义如下的接口:
- typedef LMatrix<float> LRegressionMatrix; /// @brief 回归中的ZERO分类 #ifndef REGRESSION_ZERO #define REGRESSION_ZERO 0.0f #endif /// @brief 回归中的ONE分类 #ifndef REGRESSION_ONE #define REGRESSION_ONE 1.0f #endif /// @brief 逻辑回归(分类) class LLogisticRegression { public: /// @brief 构造函数 LLogisticRegression(); /// @brief 析构函数 ~LLogisticRegression(); /// @brief 训练模型 /// 如果一次训练的样本数量为1, 则为随机梯度下降 /// 如果一次训练的样本数量为M(样本总数), 则为梯度下降 /// 如果一次训练的样本数量为m(1 < m < M), 则为批量梯度下降 /// @param[in] xMatrix 样本矩阵, 每一行代表一个样本, 每一列代表样本的一个特征 /// @param[in] yVector(列向量) 样本标记向量, 每一行代表一个样本, 值只能为REGRESSION_ONE或REGRESSION_ZERO /// @param[in] alpha 学习速度, 该值必须大于0.0f /// @return 成功返回true, 失败返回false(参数错误的情况下会返回失败) bool TrainModel(IN const LRegressionMatrix& xMatrix, IN const LRegressionMatrix& yVector, IN float alpha); /// @brief 使用训练好的模型预测数据 /// @param[in] xMatrix 需要预测的样本矩阵 /// @param[out] yVector 存储预测的结果向量(列向量), 值为REGRESSION_ONE标记的概率 /// @return 成功返回true, 失败返回false(模型未训练或参数错误的情况下会返回失败) bool Predict(IN const LRegressionMatrix& xMatrix, OUT LRegressionMatrix& yVector) const; /// @brief 计算似然值, 似然值为0.0~1.0之间的数, 似然值值越大模型越好 /// @param[in] xMatrix 样本矩阵, 每一行代表一个样本, 每一列代表样本的一个特征 /// @param[in] yVector(列向量) 样本输出向量, 每一行代表一个样本, 值只能为REGRESSION_ONE或REGRESSION_ZERO /// @return 成功返回似然值, 失败返回-1.0f(参数错误的情况下会返回失败) float LikelihoodValue(IN const LRegressionMatrix& xMatrix, IN const LRegressionMatrix& yVector) const; private: CLogisticRegression* m_pLogisticRegression; ///< 逻辑回归实现类 };
LMatrix是我们自定义的矩阵类,用于方便机器学习的一些矩阵计算,关于它的详细代码可以查看链接:猛戳我
我们为LLogisticRegression设计了三个方法TrainModel,Predict以及LikelihoodValue,用于训练模型,预测新数据以及计算似然值。
我们看一下TrainModel的实现:
- LRegressionMatrix X; Regression::SamplexAddConstant(xMatrix, X); const LRegressionMatrix& Y = yVector; LRegressionMatrix& W = m_wVector; LRegressionMatrix XT = X.T(); /* 如果h(x) = 1/(1 + e^(X * W)) 则 wj = wj - α * ∑((y - h(x)) * xj) 如果h(x) = 1/(1 + e^(-X * W)) 则 wj = wj + α * ∑((y - h(x)) * xj) */ LRegressionMatrix XW(X.RowLen, 1); LRegressionMatrix DW; LRegressionMatrix::MUL(X, W, XW); for (unsigned int m = 0; m < XW.RowLen; m++) { this->Sigmoid(XW[m][0], XW[m][0]); } LRegressionMatrix::SUB(Y, XW, XW); LRegressionMatrix::MUL(XT, XW, DW); LRegressionMatrix::SCALARMUL(DW, -1.0f * alpha, DW); LRegressionMatrix::ADD(W, DW, W);
逻辑函数定义如下:
- /// @brief S型函数 /// @param[in] input 输入值 /// @param[out] output 存储输出值 void Sigmoid(float input, float& output) const { output = 1.0f/(1.0f + exp(input)); }
我们看一下Predict的实现:
- LRegressionMatrix X; Regression::SamplexAddConstant(xMatrix, X); yVector.Reset(X.RowLen, 1, 0.0f); LRegressionMatrix::MUL(X, m_wVector, yVector); for (unsigned int m = 0; m < yVector.RowLen; m++) { this->Sigmoid(yVector[m][0], yVector[m][0]); }
最后就是LikelihoodValue的实现:
- LRegressionMatrix predictY; this->Predict(xMatrix, predictY); float likelihood = 1.0f; for (unsigned int i = 0; i < yVector.RowLen; i++) { if (yVector[i][0] == REGRESSION_ONE) likelihood *= predictY[i][0]; else if (yVector[i][0] == REGRESSION_ZERO) likelihood *= (1.0f - predictY[i][0]); else return -1.0f; } return likelihood;
以上完整的代码可以在链接:猛戳我查看,我们的逻辑回归被定义在文件LRegression.h和LRegression.cpp中。
机器学习基础1-线性回归及C++实现
机器学习基础2-二项逻辑回归及C++实现
来源: http://www.coderjie.com/blog/604c0804dbeb11e7841d00163e0c0e36