Matlab 中的精度精度

Question

我正在尝试实现一个Logistic Regression算法，部分Matlab代码如下。

for i = 1 : MAX_ITR
    % Calculate mu
    mu = 1.0 ./ (1.0 + exp(-(X * w)));

    % Calculate h to check for convergence
    h(i) = sum(-y .* log(mu) - (1 - y) .* log(1 - mu)) + (lambda / 2) * norm(w([2:end]))^2;

    % Calculate gradient and hessian.
    G = lambda .* w; 
    G(1) = 0;       % Set extra term for gradient to 0
    L = lambda .* eye(D + 1); 
    L(1) = 0;       % Set extra term for Hessian to 0
    grad = (X' * (mu - y)) + G;
    S = diag((mu .* (1 - mu)));
    H = (X' * S * X) + L;

    % Update w
    w = w - H\grad;
end

显然mu 的值不能为 1，因为指数不能为 0。但是，在某些值中，指数的计算结果非常小，例如 1.6629e-05。这会导致该实例的 mu 值非常接近 1，即 0.999983371689452。

我已经逐次运行代码迭代，前 4 次迭代很好，因为 mu 不包含任何这样的“接近 1”值。但是，第五次迭代确实如此，因此，我的h 得到了一个NaN，算法不会收敛。

我已经运行了digits 命令，它的值是32，所以我不知道是什么导致了这个问题。

预计到达时间： 在@rayryeng 的建议之后更新了代码： D是X中的特征个数，m是训练样本个数X

   for i = 1 : MAX_ITR
        % Initialize arrays
        grad = zeros(D+1,1);
        h(i) = 0;
        H = zeros(D+1,D+1);

        for j = 1 : m
            % Calculate mu
            mu = sigmoid(X(j,:) * w);

            % Calculate h (to check for convergence)
            h(i) = h(i) - (1/m)*(y(j) * log(mu) + (1 - y(j)) * log(1 - mu)) + (lambda / (2 * m)) * norm(w(2:end))^2;

            % Calculate gradient and Hessian
            G = lambda * w; 
            G(1) = 0;       % Set extra term for gradient to 0
            L = lambda * eye(numFeatures + 1); 
            L(1) = 0;       % Set extra term for Hessian to 0
            grad = grad - (1/m) * ((X(j,:)' * (mu - y(j))) + G);
            S = diag((mu .* (1 - mu)));
            H = H - (1/m) * ((X(j,:)' * S * X(j,:)) + L);
        end
        %fprintf('h(%d) = %0.5f\n', i, h(i));

        % Update w
        w = w - H\grad;
    end

我不再收到NaN 错误，并且h 值似乎在几次迭代后收敛。

Answer 1

在我的实现中，

• X是输入数据m x n（一列1s被添加到 初始数据X = [ones(m, 1) X];)

• y是输出数据。

• theta 的初始值为 0

sigmoid函数定义如下：

function g = sigmoid(z)
    g = 1./(1+exp(-z));
end

代价函数定义如下：

function [J, grad] = costFunction(theta, X, y)
    m = length(y);
    J = (1/m)*sum(-y .* log(sigmoid(X*theta)) - (1-y) .* log(1-sigmoid(X*theta)));
    grad = (X'*(sigmoid(X*theta)-y))/m;
end

在主代码中，我使用fminunc 函数来寻找最优theta。

options = optimset('GradObj', 'on', 'MaxIter', 400);
[theta, cost] = fminunc(@(t)(costFunction(t, X, y)), initial_theta, options);

最后，它与您的实现不同，但找出问题的解决方案可能会很有用。