在Matlab中实现正则化的逻辑回归代码

Question

我正在尝试使用正则化逻辑回归，这里有Matlab中的简单公式：

成本函数：

J(theta) = 1/m*sum((-y_i)*log(h(x_i)-(1-y_i)*log(1-h(x_i))))+(lambda/2*m)*sum(theta_j)

梯度:

∂J(theta)/∂theta_0 = [(1/m)*(sum((h(x_i)-y_i)*x_j)] if j=0

∂j(theta)/∂theta_n = [(1/m)*(sum((h(x_i)-y_i)*x_j)]+(lambda/m)*(theta_j) if j>1

这不是Matlab代码，只是公式。

到目前为止，我已经完成了以下内容:

function [J, grad] = costFunctionReg(theta, X, y, lambda)

J = 0;
grad = zeros(size(theta));

temp_theta = [];

%cost function

%get the regularization term

for jj = 2:length(theta)

    temp_theta(jj) = theta(jj)^2;
end

theta_reg = lambda/(2*m)*sum(temp_theta);

temp_sum =[];

%for the sum in the cost function

for ii =1:m

   temp_sum(ii) = -y(ii)*log(sigmoid(theta'*X(ii,:)'))-(1-y(ii))*log(1-sigmoid(theta'*X(ii,:)'));

end

tempo = sum(temp_sum);

J = (1/m)*tempo+theta_reg;

%regulatization
%theta 0

reg_theta0 = 0;

for jj=1:m
 reg_theta0(jj) = (sigmoid(theta'*X(m,:)') -y(jj))*X(jj,1)
end    

reg_theta0 = (1/m)*sum(reg_theta0)

grad_temp(1) = reg_theta0

%for the rest of thetas

reg_theta  = [];
thetas_sum = 0;

for ii=2:size(theta)
    for kk =1:m
        reg_theta(kk) = (sigmoid(theta'*X(m,:)') - y(kk))*X(kk,ii)
    end
    thetas_sum(ii) = (1/m)*sum(reg_theta)+(lambda/m)*theta(ii)
    reg_theta = []
end

for i=1:size(theta)

    if i == 1
        grad(i) = grad_temp(i)
    else
        grad(i) = thetas_sum(i)
    end
end
end

成本函数得出的结果是正确的，但我不知道为什么梯度（一步）就不对了。成本函数得出的J = 0.6931是正确的，而梯度grad = 0.3603 -0.1476 0.0320则不正确，成本从2开始是因为参数theta（1）不需要规则化，有没有帮助？我想代码有问题，但四天后我还是找不到错误。谢谢。

- Pedro.Alonso

Answer 1

向量化：

function [J, grad] = costFunctionReg(theta, X, y, lambda)

hx = sigmoid(X * theta);
m = length(X);

J = (sum(-y' * log(hx) - (1 - y')*log(1 - hx)) / m) + lambda * sum(theta(2:end).^2) / (2*m);
grad =((hx - y)' * X / m)' + lambda .* theta .* [0; ones(length(theta)-1, 1)] ./ m ;

end

- Franck Dernoncourt

---

请问您能否解释一下在 J 中添加正则化项 +lambda * sum(theta(2:end).^2) / (2*m) 的意义是什么，因为我们传递的 initial_theta = zeros(size(X, 2), 1);，当乘以该项时会得到值为0，那么这样做有什么意义呢？因为代价函数仍将保持不变... - Incpetor

@Inceptor361 在第一次调用 costFunctionReg 时，theta 的值为0，但在第一次迭代后它们将被更改。 - Franck Dernoncourt

这里只是指出，J = (sum(... 中的第一个 sum 是多余的，因为它的参数已经具有 {1x1} 的维度，也就是一个数字。 - ratijas

我想问一下，在 [ ((hx - y)' * X / m)' + lambda .* theta .* [0; ones(length(theta)-1, 1)] ./ m ] 这个式子中，.* 运算符只需要在 theta 和 [0,ones(... )] 之间使用，因为其他地方使用简单的 * 就足够了，我的理解对吗？ - Charbel

@FranckDernoncourt FYI，有人提出了一个新问题关于这个答案。 - TylerH

Answer 2

我使用了更多变量，这样你可以清楚地看到哪些来自常规公式，哪些来自“添加的正则化成本”。此外，在Matlab/Octave中使用“向量化”而不是循环是一个很好的实践。通过这样做，您将获得更优化的解决方案。

 function [J, grad] = costFunctionReg(theta, X, y, lambda)

    %Hypotheses
    hx = sigmoid(X * theta);

    %%The cost without regularization
    J_partial = (-y' * log(hx) - (1 - y)' * log(1 - hx)) ./ m;


    %%Regularization Cost Added
    J_regularization = (lambda/(2*m)) * sum(theta(2:end).^2);

    %%Cost when we add regularization
    J = J_partial + J_regularization;

    %Grad without regularization
    grad_partial = (1/m) * (X' * (hx -y));

    %%Grad Cost Added
    grad_regularization = (lambda/m) .* theta(2:end);

    grad_regularization = [0; grad_regularization];

    grad = grad_partial + grad_regularization;

- Fernando Cardenas

---

很好，回答不错，尽管应该是“正则化”，:) - LancelotHolmes

你能告诉我为什么要执行 grad_regularization = [0; grad_regularization]; 吗？ - Arnav Das

Answer 3

终于搞定了，重写了4遍之后，这才是正确的代码：

function [J, grad] = costFunctionReg(theta, X, y, lambda)
J = 0;
grad = zeros(size(theta));

temp_theta = [];

for jj = 2:length(theta)

    temp_theta(jj) = theta(jj)^2;
end

theta_reg = lambda/(2*m)*sum(temp_theta);

temp_sum =[];

for ii =1:m

   temp_sum(ii) = -y(ii)*log(sigmoid(theta'*X(ii,:)'))-(1-y(ii))*log(1-sigmoid(theta'*X(ii,:)'));

end

tempo = sum(temp_sum);

J = (1/m)*tempo+theta_reg;

%regulatization
%theta 0

reg_theta0 = 0;

for i=1:m
    reg_theta0(i) = ((sigmoid(theta'*X(i,:)'))-y(i))*X(i,1)
end

theta_temp(1) = (1/m)*sum(reg_theta0)

grad(1) = theta_temp

sum_thetas = []
thetas_sum = []

for j = 2:size(theta)
    for i = 1:m

        sum_thetas(i) = ((sigmoid(theta'*X(i,:)'))-y(i))*X(i,j)
    end

    thetas_sum(j) = (1/m)*sum(sum_thetas)+((lambda/m)*theta(j))
    sum_thetas = []
end

for z=2:size(theta)
    grad(z) = thetas_sum(z)
end


% =============================================================

end

如果对任何人有帮助，或者任何人对如何做得更好有评论，请留言。 :)

- Pedro.Alonso

---

谢谢，你能解释一下吗？

1. 为什么我们在成本J中跳过theta（1）？
2. 为什么我们在grad（1）中忽略了lambda / m * theta？

- user4398985

如果我记得没错的话，通过查看代码，theta(1)没有被跳过，而是单独计算的，我认为这样更容易理解。至于第二个问题，我不确定当时想要实现什么。 - Pedro.Alonso

我认为 grad(1) 被跳过正则化，是因为它对应于您添加到数据中的“1”列的权重。 - Ciprian Tomoiagă

Answer 4

这是一个消除循环的答案。

m = length(y); % number of training examples

predictions = sigmoid(X*theta);
reg_term = (lambda/(2*m)) * sum(theta(2:end).^2);
calcErrors = -y.*log(predictions) - (1 -y).*log(1-predictions);
J = (1/m)*sum(calcErrors)+reg_term;

% prepend a 0 column to our reg_term matrix so we can use simple matrix addition
reg_term = [0 (lambda*theta(2:end)/m)'];
grad = sum(X.*(predictions - y)) / m + reg_term;

- Rick