1.6.8 Cox proportional hazards (survival)
Source: "glmnet/examples/06-cox.rkt"
The Cox family fits a survival model. Each observation carries a follow-up time and a 0/1 event indicator: 1 means the event (death, failure, churn) was observed at that time, 0 means the subject was still event-free when observation stopped, so the time is only a lower bound. That second case is right-censoring, and handling it is the reason not to regress the times directly.
The model is h(t | x) = h₀(t) exp(xβ): every subject’s hazard is the same baseline h₀(t) scaled by exp(xβ). The baseline absorbs any intercept, so Cox fits have none, and cox-fit has no #:intercept? keyword. The penalty is the same as for every other family.
1.6.8.1 The data
In this small cohort x₁ is a risk factor: subjects with higher x₁ have the event sooner, and the three lowest-risk subjects are censored. x₂ is noise.
> (define X '((0.5 1.0) (1.0 2.0) (1.5 1.0) (2.0 2.0) (2.5 1.0) (3.0 2.0) (3.5 1.0) (4.0 2.0))) > (define times '(12.0 10.0 11.0 8.0 6.0 7.0 4.0 3.0)) > (define statuses '(0 0 0 1 1 1 1 1))
1.6.8.2 Fitting
> (define fit (cox-fit X times statuses #:lambda 0.1)) > fit #<glmnet:cox λ=0.1 dev=0.6516 nz=1/2>
> (coef fit) '#(1.9229529670643517 0.0)
A Cox fit has no intercept, so coef holds only the coefficients. β₁ is positive, so risk rises with x₁, and the noise coefficient is exactly 0.0. Exponentiated, a coefficient is a hazard ratio: one more unit of x₁ multiplies the hazard at every time by
> (exp (vector-ref (cox-result-coefficients fit) 0)) 6.8411303037210045
1.6.8.3 Ranking by risk
cox-linear-predictor returns xβ for each row and cox-relative-risk returns exp(xβ), the hazard relative to the baseline. They are predict with #:type 'link and #:type 'response:
> (cox-linear-predictor fit '((1.0 1.0) (2.0 1.0) (3.0 1.0))) '(1.9229529670643517 3.8459059341287034 5.768858901193055)
> (cox-relative-risk fit '((1.0 1.0) (2.0 1.0) (3.0 1.0))) '(6.8411303037210045 46.801063832489845 320.17217603082736)
Each step of one unit in x₁ multiplies the relative risk by the hazard ratio above. Because the baseline hazard is not estimated, these numbers rank subjects against each other; they do not predict a survival time or a survival curve.
1.6.8.4 Varying λ
cox-path fits a list of λ in one call, largest first (Regularization paths):
> (define (round3 x) (/ (round (* 1000 x)) 1000))
> (define (rounded v) (for/list ([b (in-vector v)]) (round3 b))) > (define path (cox-path X times statuses #:lambda '(0.01 0.1 0.5 1.0)))
> (for ([lam (in-vector (glmnet-path-lambda path))] [beta (in-vector (glmnet-path-coefficients path))] [dev (in-vector (glmnet-path-dev-ratio path))]) (printf "λ = ~a: β = ~a, dev-ratio = ~a\n" lam (rounded beta) (round3 dev)))
λ = 1.0: β = (0.0 0.0), dev-ratio = 0.0
λ = 0.5: β = (0.389 0.0), dev-ratio = 0.237
λ = 0.1: β = (1.923 0.0), dev-ratio = 0.652
λ = 0.01: β = (4.201 -1.381), dev-ratio = 0.77
With eight subjects and five events, a small penalty lets the noise predictor in with a sizeable coefficient; the lasso removes it by λ = 0.1, and by λ = 1 it has removed x₁ as well.
1.6.8.5 Censoring and events
The partial likelihood is built from the event times, so at least one observation must be an event:
> (cox-fit X times '(0 0 0 0 0 0 0 0) #:lambda 0.1) cox-fit: at least one observation must be an event (status =
1)
1.6.8.6 When to use it
Use the Cox family for time-to-event outcomes with censoring: patient survival, machine failure, customer churn. The lasso Cox model is a standard way to pick a few prognostic factors from many candidates.