Gradient Descent, Newton's Method, and IRLS
Linear regression has a formula that gives the ordinary least-squares coefficients in one calculation. Logistic regression does not. Its probabilities change nonlinearly with the coefficients, so fitting starts from a guess and improves that guess repeatedly until the log-loss stops changing.
The log-loss is convex, which means it has one bowl-shaped surface rather than many competing local minima. Gradient descent follows the slope of that bowl. Newton's method also measures its curvature. Iteratively reweighted least squares (IRLS) is the same Newton update written as a weighted linear-regression problem. All three methods begin with the same gradient and Hessian:
$$\nabla J(\beta)=\frac{1}{n}X^\top(p-y),\qquad \nabla^2 J(\beta)=\frac{1}{n}X^\top W X,\qquad W=\operatorname{diag}\big(p_i(1-p_i)\big),$$Here, $n$ is the number of rows, $p$ contains the current predicted probabilities, and $y$ contains the observed labels. The gradient gives the uphill direction and steepness. The Hessian describes how that slope changes in different directions. The running eight-row study-hours example has the final coefficients $\hat\beta_0=-5.7703$ and $\hat\beta_1=2.5646$, with log-loss $0.3129$.
Gradient Descent
Gradient descent uses the gradient alone. Subtracting the gradient moves the coefficients downhill:
$$\beta_{t+1}=\beta_t-\eta\,\nabla J(\beta_t)=\beta_t-\frac{\eta}{n}X^\top(p_t-y),$$The learning rate $\eta$ controls the distance moved. A tiny rate is safe but slow. A rate that is too large can jump across the minimum repeatedly or make the loss grow. One step costs $O(nd)$ for $n$ rows and $d$ coefficient columns. It does not require a $d\times d$ matrix, so gradient-based methods work well with many features or sparse data.
The drawback is the number of steps. Near the minimum, let $L$ and $\mu$ be the largest and smallest curvatures. Their ratio $\kappa=L/\mu$ is the condition number. A large $\kappa$ means the loss surface is a long, narrow valley: a rate small enough for the steep direction makes little progress in the flat direction. With the best fixed rate, the slowest direction keeps roughly $(\kappa-1)/(\kappa+1)$ of its error after each step. The study-hours example has $\kappa\approx94$, so that direction keeps about $97.9\%$ of its error per step. This explains the slow path shown later.
Newton's Method
Newton's method uses both slope and curvature. At the current coefficients, it replaces the loss surface with a simple quadratic bowl and moves to the bottom of that temporary bowl:
$$\beta_{t+1}=\beta_t-\big[\nabla^2 J(\beta_t)\big]^{-1}\nabla J(\beta_t)=\beta_t-\big(X^\top W_t X\big)^{-1}X^\top(p_t-y).$$The inverse in the formula explains the mathematics. A numerical implementation solves a linear system instead of forming the inverse directly. Since every diagonal value in $W$ is nonnegative, the Hessian is positive semidefinite. With independent feature columns and useful weights, it is positive definite and gives a downhill direction.
Once Newton's method is close to the solution, its error shrinks quadratically, so the number of correct digits can roughly double in one step. The study-hours example reaches the minimum in five steps. Each step is more expensive than gradient descent: forming $X^\top W X$ costs $O(nd^2)$ and solving the system costs $O(d^3)$. Newton's method therefore suits a moderate number of dense features, while gradient or quasi-Newton methods suit larger problems.
IRLS: The Weighted Least-Squares View
The Newton update can be rearranged into a weighted least-squares formula:
$$\beta_{t+1}=\big(X^\top W_t X\big)^{-1}X^\top W_t\,\xi_t,\qquad \xi_{t,i}=z_{t,i}+\frac{y_i-p_{t,i}}{p_{t,i}(1-p_{t,i})},$$Here, $z_t=X\beta_t$ is the current linear score, and $\xi_t$ is a temporary target called the working response. One IRLS iteration has four steps:
- The current scores $z_t$ and probabilities $p_t$ are computed.
- Row $i$ receives the weight $w_i=p_i(1-p_i)$.
- The working response $\xi_i$ is built from the current score and probability error.
- A weighted linear regression of $\xi$ on $X$ produces the next coefficients.
The weights and working responses change after every fit, which gives the method its name. The update is still exactly the Newton update, not a different approximation.
At the initial coefficients $\beta^{(0)}=(0,0)$, every score is $0$, every probability is $0.5$, and every weight is $0.25$. The working response is $-2$ for a label of $0$ and $+2$ for a label of $1$. Weighted least squares on these transformed labels produces $\beta^{(1)}=(-3,\,1.3333)$, exactly matching the first Newton step.
At the next iteration, the weights are no longer equal. They range from $0.0806$ for confident rows to $0.2432$ for rows near probability $0.5$. A logistic probability changes most rapidly near $0.5$, so those uncertain rows have more influence on the local quadratic fit.
▶ One IRLS Step as a Weighted Least-Squares Fit
Starting from all-zero coefficients on the study-hours data: labels become working responses, a weighted line is fit, and the result equals the first Newton step. Every value is computed at runtime.
Conditioning and Standardization
The Hessian uses the effective design $W^{1/2}X$. Its condition number measures how sensitive an optimization step is to small numerical changes. A value near $1$ is easy to handle; a very large value means some coefficient directions are much flatter than others.
Two common problems create poor conditioning:
- Different feature scales. A feature measured in millions can dominate one measured in fractions. Standardizing the non-intercept columns usually fixes this source.
- Separation. As probabilities approach $0$ or $1$, the weights $p_i(1-p_i)$ approach $0$. The Hessian then loses information even when the original features were scaled well.
A large condition number at the first iteration often points to scaling or collinearity. A condition number that grows rapidly during fitting can warn that the data is approaching separation.
Damping and Line Search
A Newton direction points downhill when the Hessian is positive definite, but the full step can still be too long when the current coefficients are far from the minimum. Damping multiplies that step by a fraction:
$$\beta_{t+1}=\beta_t-\alpha_t\big[\nabla^2 J(\beta_t)\big]^{-1}\nabla J(\beta_t),$$A backtracking line search tries $\alpha=1$, then smaller values such as $1/2$, $1/4$, and $1/8$ until the loss falls enough. The Armijo condition is a common rule for deciding what counts as enough. Near the minimum, the full step is usually accepted and Newton's fast convergence returns.
If the Hessian is nearly singular, the solver can instead use $\nabla^2J+\lambda I$. The small positive diagonal prevents an almost-flat direction from creating an enormous coefficient change. This is a numerical safeguard. A deliberate regularization penalty is still needed when the model itself requires shrinkage or the data is separated.
Convergence and Complexity
| Method | Uses | Convergence | Cost per step | Best when |
|---|---|---|---|---|
| Gradient descent | gradient | steady but often slow | $O(nd)$ | many or sparse features |
| Newton / IRLS | gradient and full Hessian | very fast near the minimum | $O(nd^2+d^3)$ | moderate, dense feature sets |
| Quasi-Newton (L-BFGS) | gradient and a small curvature history | between the two | $O(nd)$ plus history | many features with smooth loss |
The animation compares the methods on the same study-hours data. Both start from zero coefficients. Newton takes a few curvature-adjusted steps toward $(-5.7703,\,2.5646)$, while fixed-rate gradient descent moves more slowly. The loss plot on the right shows the progress after each iteration.
▶ Gradient Descent Versus Newton: Path and Loss
Left, the coefficient trajectory of each method; right, the log-loss per iteration. Both are computed at runtime from the eight rows, with gradient descent at a fixed rate and Newton using the exact Hessian.
Separation and Divergence
Data is linearly separable when one boundary classifies every training row correctly. Scaling the boundary's coefficients upward then pushes every predicted probability closer to its correct extreme. The log-loss keeps moving toward zero, so unregularized logistic regression has no finite best coefficient vector.
Gradient descent responds by growing the coefficient norm without settling. Newton's method pushes the probabilities toward $0$ and $1$, which makes the weights vanish and the Hessian singular. Enormous coefficients, near-zero weights, or a loss of Hessian rank therefore indicate separation rather than successful convergence. The next page explains how a penalty restores a finite solution.
Edge Cases and the Implementation
- Learning rate too large. Fixed-rate gradient descent diverges for $\eta\ge2/L$. Feature scaling changes $L$, so it also changes the safe range.
- Singular Hessian. Collinear columns or vanishing weights can make $X^\top W X$ singular. A small diagonal term or a rank-aware solver keeps the numerical step defined.
- Weak stopping rule. Reliable fitting checks the gradient, relative loss change, coefficient change, and an iteration limit. A rapidly growing coefficient norm needs a separate separation warning.
- Stochastic fitting. The shuffle seed, batch size, learning-rate schedule, and stopping rule are part of the reproducible model definition.
The C++ gradient and Hessian section computes the shared quantities used by every method. Its solver section then applies damped Newton updates. A dedicated Newton and IRLS implementation chapter is also planned.