Survival Analysis With Genomic and Proteomic Features

Why Time-to-Event Outcomes Need Their Own Analysis

A multi-omics study may ask whether a molecular profile is associated with time until relapse, progression, or death. A simple comparison of average observed follow-up times is generally inadequate: some participants have not experienced the event when observation ends, and others may leave follow-up earlier. Their event times are not known; treating their last observed times as if they were event times would distort the analysis.

Survival analysis uses the information that a participant remained event-free up to a known time, even when the event itself was not observed. In an integrated genomic and proteomic analysis, it can estimate how features are associated with the instantaneous event rate while adjusting for other measured features and clinical covariates. These are adjusted associations, not automatically causal effects.

Definition: Let \(T\) be the time from a specified origin to the event of interest, and let \(C\) be the time from that origin to censoring. The observed time is \(Y=\min(T,C)\), and the event indicator is \(\delta=I(T\le C)\), where \(\delta=1\) means the event was observed and \(\delta=0\) means it was right-censored. At time \(t\), a participant is in the risk set if they are still under observation and event-free immediately before \(t\).

For a person with features \(z=(z_1,\ldots,z_p)\), the hazard \(h(t\mid z)\) is the instantaneous event rate at time \(t\), conditional on being event-free just before \(t\) and on the modeled features. It is a rate, not the probability of experiencing the event at exactly time \(t\). The hazard ratio compares these rates for two feature profiles at the same time.

The Cox Model and Its Assumptions

The Cox proportional-hazards model is useful when the goal is to estimate feature associations without specifying a particular parametric shape for the baseline hazard. For example, \(z\) might include standardized expression of one gene, a protein abundance measurement, age, and treatment group. These predictors must refer to the same defined time origin and be available at or before that origin for the intended baseline analysis.

Definition: The Cox proportional-hazards model is

h(t\mid z)=h_0(t)\exp(\beta_1z_1+\cdots+\beta_pz_p).

Here \(h_0(t)\) is an unspecified baseline hazard and \(\beta_j\) is the log hazard ratio for a one-unit increase in \(z_j\), holding the other included predictors fixed. For a one-unit increase, the hazard ratio is \(\exp(\beta_j)\). The model assumes that this ratio does not change with time.

The word adjusted means that a coefficient compares participants differing in one modeled predictor while the other modeled predictors are held fixed. It does not mean that all confounding has been removed. A feature measured on a different scale has a different one-unit interpretation; standardizing a feature makes its coefficient refer to a one-standard-deviation increase, not to a raw assay unit.

Interpretation depends on assumptions. Participants should be independent across individuals, the event and time origin should be defined consistently, and the log hazard should be adequately represented by the included predictors. Right censoring must be non-informative conditional on the modeled predictors: among participants with the same modeled information, censoring should not depend on their otherwise unobserved future event time. Proportional hazards requires each modeled feature's hazard ratio to remain constant over time. The data also need enough events and sufficient overlap in feature values to estimate the coefficients stably.

These assumptions require scientific scrutiny as well as statistical checks. For instance, if high-risk patients are more likely to be lost to follow-up because their health deteriorates, censoring may be informative even after adjustment for the recorded features. A competing event, such as death before relapse when relapse is the endpoint, is not automatically ordinary censoring: it may change the scientific quantity being estimated and call for a competing-risks analysis.

How the Partial Likelihood Is Constructed

The Cox model estimates \(\beta\) using the partial likelihood. At an observed event time, it compares the feature profile of the person who had the event with the profiles of everyone in that time's risk set. The unknown baseline hazard cancels from this comparison. For simplicity, the expression below assumes distinct event times; tied event times require a tie-handling method, such as the Efron approximation or an exact partial likelihood.

Definition: If \(D\) is the set of observed event times, \(i(t)\) is the participant who experiences the event at time \(t\), and \(R(t)\) is the risk set just before \(t\), the partial likelihood for distinct event times is

L(\beta)=\prod_{t\in D} \frac{\exp\{\beta^\mathsf{T}z_{i(t)}\}} {\sum_{j\in R(t)}\exp\{\beta^\mathsf{T}z_j\}}.

The estimate maximizes this likelihood, or equivalently its logarithm. Each factor is the model-based relative weight of the observed event participant among those at risk at that time.

This construction explains both the strength and the boundary of the method. Estimation of feature effects does not require choosing a baseline hazard distribution, but the model still requires proportional hazards and a defensible censoring assumption. The partial likelihood also does not make correlated genomic and proteomic measurements independent or solve confounding by itself.

Worked Example: A Partial-Likelihood Contribution

Suppose three hypothetical participants are at risk at month 4. Participant A experiences the event at month 4; B and C remain event-free at that time. A has a binary genomic feature value \(x=1\), while B and C each have \(x=0\). Consider the contribution to the partial likelihood at the candidate coefficient \(\beta=\log(2)\).

Model and assumptions for this calculation. Use a Cox model with this single feature, assume distinct event times, and suppose the event and risk set have been recorded correctly. This calculation evaluates the likelihood contribution at a candidate coefficient; it does not claim that the candidate is the fitted estimate.

At month 4, the numerator is \(\exp\{\log(2)\times1\}=2\). The denominator includes all three risk-set members: A contributes 2, while B and C each contribute \(\exp\{\log(2)\times0\}=1\). Therefore, the event-time contribution is

\frac{2}{2+1+1}=\frac{1}{2}=0.5.

If B has an event later, at month 7, and C remains at risk, the risk set at that later event contains B and C, both with \(x=0\). That contribution is \(1/(1+1)=0.5\). The partial likelihood for these two events at the candidate coefficient is thus \(0.5\times0.5=0.25\). The example shows how the risk sets, rather than merely the overall numbers of events in each feature group, determine the Cox likelihood.

Interpreting an Adjusted Feature Association

A coefficient is expressed on the log-hazard scale. Exponentiating it gives a hazard ratio, which is usually easier to describe. The reference is a one-unit increase in the feature as it was entered into the model. Always specify whether that means one assay unit, one standard deviation, or a different defined increment.

Worked Example: A Genomic Feature Adjusted for Protein Abundance

Consider a hypothetical cohort in which the Cox model includes a standardized genomic expression feature, a standardized protein abundance feature, and age. The estimated coefficient for genomic expression is \(\hat{\beta}=0.405\), with \(SE=0.180\). Assume the model is otherwise suitable, its proportional-hazards assumption is plausible, and censoring is non-informative conditional on these predictors.

For a one-standard-deviation increase in genomic expression, holding protein abundance and age fixed, the estimated hazard ratio is

\exp(0.405)\approx1.499.

The estimated hazard is therefore about \(49.9\%\) higher at any given time for participants differing by one standard deviation in genomic expression but having the same modeled protein abundance and age. This is an adjusted association with the hazard, not a statement that the probability of an event is 49.9% higher, and not evidence by itself that changing expression would cause the event rate to change.

A large-sample Wald calculation gives \(z=\hat{\beta}/SE=0.405/0.180=2.25\), with a two-sided \(p\)-value of approximately 0.0244. A 95% confidence interval on the log-hazard scale is \(0.405\pm1.96(0.180)\), or \((0.0522,0.7578)\). Exponentiating the endpoints gives an approximate 95% confidence interval for the hazard ratio of \((1.054,2.134)\). Under the model and the procedure's repeated-sampling assumptions, this interval is compatible with hazard ratios from about 1.054 to 2.134 per standard deviation increase, conditional on the other included predictors. It is not a 95% probability statement about this particular interval containing the parameter.

Checking Proportional Hazards and Censoring

A fitted coefficient can be misleading if its association changes over follow-up. Proportional hazards can be assessed by examining scaled Schoenfeld residuals against time, testing their association with time, and using plots of estimated log cumulative hazards across groups as a visual check. A systematic time trend in residuals is a warning that a constant hazard ratio may not describe the data. These checks may have low power when there are few events, so a non-significant test does not establish that proportional hazards holds.

A direct diagnostic is to add a feature-by-time term. For example, replacing a constant coefficient \(\beta\) by \(\beta+\gamma\log(t/12)\), with time measured in months, permits the feature's association to vary with time. A nonzero \(\gamma\) indicates time variation under this specified model. The model should be chosen and interpreted with care; adding a time interaction changes the assumption and the meaning of the coefficient.

Worked Example: A Time-Varying Feature Association

In a hypothetical diagnostic analysis, a one-standard-deviation protein feature is modeled with coefficient \(\beta+\gamma\log(t/12)\), where \(t\) is months since enrollment. The estimated time-interaction coefficient is \(\hat{\gamma}=0.320\), with \(SE=0.140\). Suppose participants are independent, event times and censoring are recorded correctly, and the specified log-time interaction is a reasonable diagnostic form.

The Wald statistic for \(H_0:\gamma=0\) is \(z=0.320/0.140=2.286\), giving a two-sided \(p\)-value of approximately 0.0223. This provides evidence against a constant coefficient in this specified model, so a single proportional-hazards ratio for this feature may be inadequate. For a time increase by a factor of \(e\), the feature's hazard ratio is multiplied by \(\exp(0.320)\approx1.377\), because \(\log(et/12)-\log(t/12)=1\). This describes a change in the feature hazard ratio across time, not a 37.7% increase in an individual participant's hazard merely from the passage of time.

Censoring also needs inspection. Suppose, among 20 hypothetical high-feature participants and 20 low-feature participants who were eligible for the same follow-up schedule, 8 high-feature and 2 low-feature participants leave follow-up early. The observed early-loss proportions are \(8/20=0.40\) and \(2/20=0.10\), respectively, a fourfold descriptive ratio. This imbalance is a reason to investigate whether loss is related to prognosis after conditioning on modeled variables; it does not, by itself, prove informative censoring. Compare reasons and timing of loss, follow-up procedures, and relevant clinical predictors across feature levels. Statistical checks cannot establish independence from unobserved future event times, so study-design knowledge remains essential.

Other diagnostics help distinguish a proportional-hazards failure from other model problems. Martingale residuals can flag an inadequate functional form for a continuous predictor; influential-case diagnostics can reveal whether a small number of events dominate a coefficient. Large standard errors, unstable estimates under modest changes in feature selection, or very sparse event counts relative to the number of predictors signal that the model may be overfit. With many omics features, prespecification, dimension reduction, or penalization may be necessary, but ordinary coefficient confidence intervals and \(p\)-values after data-driven selection can be too optimistic.

What the Cox Model Does—and Does Not—Provide

The Cox model is not the only way to analyze time-to-event data. A Kaplan–Meier curve can describe event-free survival over time for a group, while a log-rank test compares groups without adjusting for multiple covariates. These methods are useful descriptive or unadjusted alternatives, but they do not estimate a multivariable adjusted feature coefficient. Parametric survival models specify a distributional form for event times and can support direct estimates of survival times or absolute risks when their assumptions are credible. Accelerated failure-time models express associations in terms of time ratios rather than hazard ratios. Choice among these approaches depends on the scientific estimand, assumptions, and data quality.

Cox hazard ratios are conditional on remaining event-free at each time, and are not generally interchangeable with risk ratios or differences in survival probability. A hazard ratio can also be hard to interpret clinically when hazards cross or when proportional hazards is implausible. If absolute prognosis is the goal, report survival probabilities at clinically meaningful times as well, with appropriate uncertainty and a stated baseline estimation method.

In multi-omics applications, correlated features can make individual adjusted coefficients imprecise or unstable: two layers may measure overlapping biology, so their conditional associations can differ substantially from their marginal associations. Cohort and platform effects, batch-related confounding, missing assays, and feature-selection leakage can also compromise interpretation. The earlier material on diagnosing cohort and platform confounding is relevant to this risk. A credible report states the time origin, event definition, censoring rules, feature scaling, covariates, number of events, handling of ties and missing profiles, diagnostics, and validation strategy.

Key takeaway: A Cox coefficient describes an adjusted log hazard ratio under a model with proportional hazards. Interpret its exponentiated value in the feature's actual units, and accompany it with evidence about model form, censoring, uncertainty, and the cohort context.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions below.

A protein coefficient is \(-0.223\) per one-standard-deviation increase. What is the hazard ratio, and how should it be interpreted?
Why does a right-censored participant still contribute information to a Cox analysis?
A scaled Schoenfeld residual plot shows a consistent trend with time. What assumption is in doubt, and what is one diagnostic follow-up?
Why is a competing event not always appropriately treated as ordinary right censoring?

Solutions. (1) \(\exp(-0.223)\approx0.800\). Conditional on other modeled predictors, the estimated hazard is about 20.0% lower per one-standard-deviation increase, if proportional hazards and the other assumptions are adequate. (2) The participant contributes to each risk set up to the censoring time, documenting event-free follow-up through that time, even though their later event time is unknown. (3) Proportional hazards is in doubt; one can test or model a feature-by-time interaction and examine whether it accounts for the trend. (4) A competing event may prevent the event of interest from occurring and may have a distinct scientific interpretation; censoring it as though its occurrence were independent can target the wrong quantity or bias inference.