Tutorials › Marginal Structural Models › Weighted Cox Score Equations

Fitting MSMs for Survival Outcomes · Tutorial 168 of 501

Weighted Cox Score Equations

Learn to derive and interpret the weighted Cox score equations, check their curvature, and understand what the treatment-effect equation compares across weighted risk sets.

Advanced 11 min read

What You'll Learn

  • Differentiate the weighted Cox partial log-likelihood to obtain its score equation.
  • Distinguish the role of weights for event contributions from their role in risk-set averages.
  • Interpret the treatment score as a weighted event-versus-risk-set contrast.
  • Derive the score derivative and identify conditions for concavity and a unique solution.
  • Recognize numerical and modeling diagnostics that can signal an unstable or invalid fit.

Why the Score Equation Matters

The weighted Cox partial likelihood turns a marginal structural Cox model into an estimating problem: choose the coefficient vector whose weighted score is zero. The score makes the mechanics visible. At each event time, it compares the covariate of the person who failed with the covariate average among those still at risk, using the weights in both parts of that comparison.

This equation is useful for more than computation. Its treatment component reveals what information identifies the treatment coefficient: weighted differences between treatment among people who fail and treatment composition in the corresponding weighted risk sets. A score root is not automatically a causal result, however. Its interpretation still depends on the causal assumptions, weights, and marginal structural model being appropriate.

Definition of the Weighted Score

Suppose there are \(n\) independent individuals. For individual \(i\), let \(N_i(t)\) be the event counting process, \(Y_i(t)\) the at-risk indicator, \(Z_i(t)\) the covariate vector in the Cox marginal structural model, and \(w_i(t)\) the nonnegative weight at time \(t\). We consider distinct event times and assume that covariates and weights used at \(t\) are available immediately before \(t\). At each event time, the weighted risk-set sum must be finite and positive.

For a coefficient vector \(\beta\), define the weighted risk-set sums

$$ S_w^{(0)}(\beta,t) =\sum_{j=1}^n Y_j(t)w_j(t)\exp\{\beta^\mathsf{T}Z_j(t)\}, \qquad S_w^{(1)}(\beta,t) =\sum_{j=1}^n Y_j(t)w_j(t)\exp\{\beta^\mathsf{T}Z_j(t)\}Z_j(t). $$

The weighted partial log-likelihood is

$$ \ell_w(\beta) =\sum_{i=1}^n\int w_i(t) \left[ \beta^\mathsf{T}Z_i(t)-\log S_w^{(0)}(\beta,t) \right]\,dN_i(t). $$

Thus a failure contributes its own weight to the event contribution, while the risk-set denominator uses the weights of all individuals currently at risk. The score is the gradient of this log-likelihood. A weighted Cox fit solves \(U_w(\beta)=0\), where the equation below follows by differentiation.

Definition: The weighted Cox score is $$ U_w(\beta) =\sum_{i=1}^n\int w_i(t) \left[ Z_i(t)-\frac{S_w^{(1)}(\beta,t)}{S_w^{(0)}(\beta,t)} \right]\,dN_i(t). $$ A coefficient estimate \(\hat\beta\) is a solution of \(U_w(\hat\beta)=0\), when such a solution exists.

Deriving the Score and Its Curvature

At a fixed time \(t\), the derivative of \(\beta^\mathsf{T}Z_i(t)\) is \(Z_i(t)\). For the denominator, differentiating the finite sum gives

$$ \nabla_\beta S_w^{(0)}(\beta,t)=S_w^{(1)}(\beta,t), \qquad \nabla_\beta\log S_w^{(0)}(\beta,t) =\frac{S_w^{(1)}(\beta,t)}{S_w^{(0)}(\beta,t)}. $$

Substituting these derivatives into the log-likelihood yields the displayed score. Because the observed event times form a finite sum in a finite data set, this differentiation requires no interchange of an infinite sum or integral with a limit.

The derivative of the score describes how sharply the equation changes around a candidate solution. Define the weighted risk-set probabilities

$$ \pi_j(\beta,t) =\frac{Y_j(t)w_j(t)\exp\{\beta^\mathsf{T}Z_j(t)\}} {S_w^{(0)}(\beta,t)}. $$

At a time with a nonempty weighted risk set these probabilities sum to one. Their mean covariate is \(S_w^{(1)}/S_w^{(0)}\); let \(\operatorname{Cov}_{\pi(\beta,t)}(Z,t)\) denote the covariance matrix of \(Z_j(t)\) under these probabilities.

Result: The derivative of the weighted score is $$ \nabla_\beta U_w(\beta) =-\sum_{i=1}^n\int w_i(t)\, \operatorname{Cov}_{\pi(\beta,t)}(Z,t)\,dN_i(t). $$ In particular, it is negative semidefinite when the weights are nonnegative.

To verify the result, write \(\bar Z_w=S_w^{(1)}/S_w^{(0)}=\sum_j\pi_jZ_j\). Differentiating this ratio gives

$$ \nabla_\beta\bar Z_w =\frac{S_w^{(2)}}{S_w^{(0)}}- \left(\frac{S_w^{(1)}}{S_w^{(0)}}\right) \left(\frac{S_w^{(1)}}{S_w^{(0)}}\right)^\mathsf{T} =\operatorname{Cov}_{\pi(\beta,t)}(Z,t), $$

where \(S_w^{(2)}=\sum_jY_jw_j\exp\{\beta^\mathsf{T}Z_j\}Z_jZ_j^\mathsf{T}\). The event covariate term \(Z_i(t)\) does not depend on \(\beta\). Differentiating the remaining negative risk-set mean in the score therefore gives the stated negative covariance sum.

For any vector \(a\), each covariance matrix satisfies \(a^\mathsf{T}\operatorname{Cov}_{\pi}(Z,t)a =\sum_j\pi_j\{a^\mathsf{T}Z_j-\sum_k\pi_k a^\mathsf{T}Z_k\}^2\geq0\). Multiplying by nonnegative event weights and summing shows \(a^\mathsf{T}\nabla_\beta U_w(\beta)a\leq0\). Hence the weighted partial log-likelihood is concave. If the summed risk-set covariance is positive definite in every nonzero direction throughout the parameter region, the Hessian is negative definite there, so the log-likelihood is strictly concave and has at most one finite maximizer. Concavity does not guarantee that a finite maximizer exists. \(\square\)

Worked Score Calculations

Worked Example: One Weighted Event Contribution

At a particular event time, suppose the weighted risk set has total weight \(3\) for untreated individuals and total weight \(4\) for treated individuals. The model has a single binary covariate \(A\), with coefficient \(\beta\). Evaluate the score contribution at \(\beta=\log 2\) when a treated individual with event weight \(3\) fails.

The untreated risk-set contribution to \(S_w^{(0)}\) is \(3\exp(0)=3\). The treated contribution is \(4\exp(\log 2)=8\). Thus \(S_w^{(0)}=11\), \(S_w^{(1)}=8\), and the risk-set mean of treatment is \(8/11\). The event subject has \(A=1\), so the score contribution is

$$ 3\left(1-\frac{8}{11}\right)=\frac{9}{11}\approx0.818. $$

This is a positive contribution: at that event time, the failing individual's treatment value exceeds the weighted risk-set mean. It is only one contribution to the full score; the fitted coefficient must balance contributions over all observed events.

Worked Example: Treatment-Effect Score Across Two Events

Consider a model with treatment \(A\) as its only covariate. At \(\beta=0\), suppose the first event's weighted risk set has treated proportion \(0.4\), and an untreated person with event weight \(2\) fails. At the second event, the weighted risk set has treated proportion \(0.7\), and a treated person with event weight \(1\) fails.

For a binary treatment, the weighted risk-set mean is its weighted treated proportion. The first event contributes \(2(0-0.4)=-0.8\). The second contributes \(1(1-0.7)=0.3\). Therefore the total treatment score at zero, from these two events, is \(-0.8+0.3=-0.5\). The negative sign says that, at this coefficient value, the weighted event treatment values are below their corresponding weighted risk-set means in aggregate. It does not by itself establish the sign of the eventual fitted coefficient; the score must be considered as a function of \(\beta\).

Worked Example: A Newton Update and Its Arithmetic

Suppose two distinct event times each have weighted treated proportion \(0.5\) at \(\beta=0\). The total event weight among treated failures is \(3.25\), and the total event weight among untreated failures is \(4.25\). For a binary treatment, an event with weight \(q\) contributes \(q(A-p)\), where \(p\) is the weighted treated proportion. The total score is

$$ U_w(0)=3.25(1-0.5)+4.25(0-0.5) =1.625-2.125=-0.5. $$

The derivative contribution at each event is \(-q p(1-p)\). Since \(p(1-p)=0.25\), the total derivative is

$$ U_w'(0)=-(3.25+4.25)(0.25)=-7.5(0.25)=-1.875=-\frac{15}{8}. $$

A Newton step for solving \(U_w(\beta)=0\), starting at zero, is \(\beta_{\mathrm{new}}=0-U_w(0)/U_w'(0)\). Substitution gives

$$ \beta_{\mathrm{new}} =0-\frac{-1/2}{-15/8} =-\frac{1}{2}\frac{8}{15} =-\frac{4}{15} \approx-0.267. $$

Both the score and its derivative are negative at zero, so their ratio is positive and the Newton subtraction moves the estimate in the negative direction. This is a single local update, not a claim that the score is exactly zero at \(-4/15\); further iterations may be needed.

Interpreting the Treatment Equation

When treatment \(A_i(t)\) is binary and is the only covariate, the treatment component of the score is

$$ U_{w,A}(\beta) =\sum_{i=1}^n\int w_i(t) \left[ A_i(t)-\bar A_w(\beta,t) \right]\,dN_i(t), \qquad \bar A_w(\beta,t) =\frac{\sum_jY_j(t)w_j(t)e^{\beta A_j(t)}A_j(t)} {\sum_jY_j(t)w_j(t)e^{\beta A_j(t)}}. $$

At each failure, the equation contrasts the treatment of the person who failed with the treatment average in that event's risk set, tilted by the weights and by the current model coefficient. The estimating equation sets the sum of these weighted contrasts to zero. With additional covariates, the treatment component still uses the treatment entry minus its weighted risk-set mean, where the risk-set probabilities depend on the full linear predictor.

In a marginal structural analysis, causal interpretation requires more than solving this equation. The treatment and censoring weights must be appropriate for the target population; the relevant positivity and sequential exchangeability conditions must hold; and the specified marginal structural Cox model must adequately represent the target hazard relationship. Under those conditions, weighting is intended to make the measured treatment histories behave as if assigned according to the target intervention. An incorrectly specified model or unstable weights can undermine that interpretation even when the numerical score is zero.

Diagnostics, Alternatives, and Limitations

Several checks follow directly from the score equation and its derivation:

  • Check risk-set construction. At every event, confirm that the event subject is in the risk set just before failure and that treatment history and weight correspond to the same time interval.
  • Check both weight placements. Event contributions use the failing individual's weight; the denominator uses weights for all individuals at risk. Omitting either part changes the estimating equation.
  • Inspect score and curvature. A numerical routine should report a small score at the proposed solution. Very flat curvature, a nearly singular Hessian, or large changes across iterations can indicate weak information or near-separation.
  • Inspect weights and weighted risk sets. A few very large weights can make the fit depend heavily on a small number of individuals. Empty or nearly empty weighted treatment groups at event times provide little or no information about the treatment contrast.
  • Check model and causal assumptions. A converged fit cannot diagnose unmeasured confounding, failure of positivity, incorrect weight models, or a misspecified marginal structural Cox model by itself.

The ordinary unweighted Cox score is the special case \(w_i(t)=1\) for everyone. It is not generally interchangeable with a weighted score: it targets the observed risk sets rather than the reweighted target population. Directly fitting an ordinary Cox model after treating weights as if they were ordinary frequency counts can also give a different equation, because the risk-set denominator and event contribution need not be weighted in the same way. Weighted Cox estimation is useful when the causal target calls for reweighting, but it remains a model-based estimating procedure; it does not remove the need to justify the weights or the proportional-hazards form.

A common mistake is to weight only the failures and leave the risk-set mean unweighted, or to weight the risk-set denominator but forget the event weight. The correction is to differentiate the weighted partial log-likelihood actually being used: its event term and its denominator jointly determine the score. Another common mistake is to treat a zero score as proof of a causal effect. It is only the estimating equation's balance condition; causal interpretation depends on the assumptions stated above.

Check Your Understanding

Try each question, then compare your reasoning with the solutions.

  1. At one event time, the weighted risk-set treated proportion is \(0.6\). An untreated individual with event weight \(2\) fails. What is this event's treatment-score contribution?
  2. Why does the derivative of the treatment score have a nonpositive sign when weights are nonnegative?
  3. At \(\beta=0\), suppose \(U_w(0)=0.3\) and \(U_w'(0)=-0.6\). What is one Newton update from zero?
  4. Does a converged weighted score equation alone establish that the treatment coefficient has a causal interpretation? Explain briefly.

Solutions. (1) The contribution is \(2(0-0.6)=-1.2\). (2) The derivative is the negative sum of event weights times weighted risk-set variances; variances and weights are nonnegative. (3) The update is \(0-0.3/(-0.6)=0.5\). (4) No. Convergence establishes only a numerical score root; causal interpretation also requires appropriate weights, positivity, sequential exchangeability, and an adequate marginal structural model.

Clinical Biostats