线性模型对数几率回归C语言,[线性模型] 对数几率回归（Logistic Regression）

公式推导

对数几率回归用于处理二分类问题，其数学基础为对数几率函数，是一种 Sigmoid 函数

\[y = \frac{1}{1+e^{-z}} \tag 1

\]

其函数图像如下

取 $z = \boldsymbol{w}^T\boldsymbol{x}+b$，并对式 $(1)$ 进行一定变换，得

$$

\ln\frac{y}{1-y}= \boldsymbol{w}^T\boldsymbol{x}+b \tag 2\\

$$

可以理解为，使用线性回归的预测结果逼近真实标记的对数几率. 当 $y>0.5$ 时，左式大于 $0$；当 $y<0.5$ 时，左式小于$0$.

按照线性回归模型的求解过程，推导出求解 \(\boldsymbol w\) 和 \(b\) 的迭代公式. 将式 \((2)\) 中的 \(y\) 看作类后验概率估计\(p(y=1|\boldsymbol x)\), 代入式 \((2)\) 可以推出

\[\begin{align}

p(y=1|\boldsymbol x) = \frac{e^{\boldsymbol w^T\boldsymbol x+b}}{1+e^{\boldsymbol w^T \boldsymbol x+b}} \tag 3\\

p(y=0|\boldsymbol x) = \frac{1}{1+e^{\boldsymbol w^T \boldsymbol x+b}} \tag 4\\

\end{align}

\]

令 \(\hat{\boldsymbol w} = (\boldsymbol w;b), \hat{\boldsymbol x_i}=(\boldsymbol x_i;1)\) 采用极大似然估计法估计 \(\hat{\boldsymbol w}\) ,最大化函数

\[L(\hat{\boldsymbol w}) = \sum\limits_{i=1}^m\ln p(y_i|\hat{\boldsymbol x_i},\hat{\boldsymbol w}) \tag 5\\

\]

令\(p_1(\hat{\boldsymbol x_i},\hat{\boldsymbol w}) = p(y=1|\hat{\boldsymbol x_i},\hat{\boldsymbol w}), p_0 = 1-p_1(\hat{\boldsymbol x_i},\hat{\boldsymbol w})\)，代入式 \((5)\)，得

\[L(\hat{\boldsymbol w}) = \sum\limits_{i=1}^m\ln (y_ip_1(\hat{\boldsymbol x_i},\hat{\boldsymbol w})+(1-y_i)p_0(\hat{\boldsymbol x_i},\hat{\boldsymbol w})) \tag 6\\

\]

代入式 \((3)(4)\)，最大化式 \((6)\) 等价于最小化

\[L(\hat{\boldsymbol w}) = \sum\limits_{i=1}^m (-y_i\hat{\boldsymbol w}^T\hat{\boldsymbol x_i}+\ln(1+e^{\hat{\boldsymbol w}^T\hat{\boldsymbol x_i}}))\tag 7\\

\]

使用牛顿迭代法得到迭代更新公式

\[\begin{align}

\hat{\boldsymbol w} &\leftarrow \hat{\boldsymbol w}-\left(\frac{\partial^2L(\hat{\boldsymbol w})}{\partial \hat{\boldsymbol w} \partial\hat{\boldsymbol w}^T} \right )^{-1}\frac{\partial L(\hat{\boldsymbol w})}{\partial\hat{\boldsymbol w}} \tag 8\\

\frac{\partial L(\hat{\boldsymbol w})}{\partial\hat{\boldsymbol w}} &=-\sum\limits_{i=1}^m \hat{\boldsymbol x_i}(y_i-p_1(\hat{\boldsymbol x_i},\hat{\boldsymbol w})) \tag 9\\

\frac{\partial^2L(\hat{\boldsymbol w})}{\partial \hat{\boldsymbol w} \partial\hat{\boldsymbol w}^T} &=\sum\limits_{i=1}^m\hat{\boldsymbol x_i}\hat{\boldsymbol x_i}^Tp_1(\hat{\boldsymbol x_i},\hat{\boldsymbol w})(1-p_1(\hat{\boldsymbol x_i},\hat{\boldsymbol w})) \tag {10}

\end{align}

\]

其中式 \((9)\) 可以向量化为

\[\frac{\partial L(\hat{\boldsymbol w})}{\partial\hat{\boldsymbol w}} = \boldsymbol X^T(p_1({\boldsymbol X},\hat{\boldsymbol w})-\boldsymbol y) \tag{11}

\]

MATLAB 实现

% 生成随机训练样本，直线 y=0.7x+200 上方的为正例，下方为反例

% 输出训练样本在坐标轴上的分布

x = zeros(100, 2);

y = zeros(100, 1);

kb = [0.7,200];

figure;

hold on;

for i = 1:100

x(i,1) = randi(1000,1);

x(i,2) = randi(1000,1);

if kb(1)*x(i,1)+kb(2)>x(i,2)

plot(x(i,1), x(i,2), 'r*');

y(i) = 1;

else

plot(x(i,1),x(i,2), 'b*');

y(i) = 0;

end

end

% 牛顿迭代法解权重

function w = cal(X,y,eps)

[m,n] = size(X);

X = [X ones(m,1)];

n = n + 1;

w = zeros(n, 1);

w(3) = 1000;

prew = zeros(n, 1);

while (true)

flag = 0;

sum1 = X'*(1-1./(1+exp(X*prew))-y);

sum2 = 0;

for i = 1:m

sum2 = sum2 + X(i,:)*X(i,:)'*(1-1./(1+exp(X(i,:)*prew)))*(1./(1+exp(X(i,:)*prew)));

end

w = prew - sum1./sum2;

for i = 1:n

if abs(w(i)-prew(i))>eps

flag = 1;

end

end

if flag==0

break;

end

prew = w;

end

end

% 测试并输出训练效果图像

eps = 0.0001;

w = cal(x, y, eps);

figure;

hold on;

for i=1:100

yy = w(1)*x(i,1)+w(2)*x(i,2)+w(3);

if yy>0

plot(x(i,1), x(i,2), 'r*');

else

plot(x(i,1), x(i,2), 'b*');

end

end

训练集的真实分布(红色正例、蓝色反例)：

模型的分类效果：