Comparing groups

Once you have survival curves for two or more groups, the natural question is whether the differences between them are real or could be explained by chance. The log-rank test and its relatives answer this. They compare the observed number of events in each group against the number expected if the groups shared a common survival experience, accumulated across every event time. This guide page explains the test, its weighted variants, and how to read the results.

We continue with the lung dataset, comparing survival between the two sexes.

import greenwood as gw

# Load the bundled lung dataset as a Polars DataFrame
lung = gw.load_dataset("lung", backend="polars")

# Build a right-censored response; status == 2 marks a death
y = gw.Surv.right(lung["time"], event=(lung["status"] == 2))

# Display the response summary
y

Surv(type=right, n=228, events=165)

The response y is a Surv object pairing each subject’s follow-up time with an event indicator, where a / marks a censored observation. The grouping variable lung["sex"] splits these same subjects into the two groups we are about to compare.

The log-rank test

The log-rank test is the standard method for comparing survival across groups. At each event time it forms a small contingency table of who was at risk and who had an event in each group, computes the events expected under the null hypothesis of no difference, and combines the observed-minus-expected discrepancies into a single chi-square statistic.

# Load data and build the response
lung = gw.load_dataset("lung", backend="polars")
y = gw.Surv.right(lung["time"], event=(lung["status"] == 2))

# Run the log-rank test comparing survival between the two sexes
result = gw.logrank_test(y, group=lung["sex"])

The call returns a TestResult object, which we store in result. Displaying it prints a readable summary of the whole test at once.

# Print the test statistic, degrees of freedom, and p-value
result

TestResult(method='Log-rank test', statistic=10.3267, df=1, p_value=0.001311)

The result carries the test statistic, its degrees of freedom (one less than the number of groups), and the p-value. A small p-value is evidence that the survival curves differ. You can also inspect the weighted observed and expected event counts per group, which show the direction of the difference. These live on the result object as the observed and expected attributes, one entry per group.

# Inspect the per-group event counts to see the direction of the difference
print("observed:", result.observed)
print("expected:", result.expected)

observed: {1: 112.0, 2: 53.0}
expected: {1: 91.58173902957279, 2: 73.41826097042721}

Here one group has more observed events than expected and the other fewer, which is the signature of a survival difference. The test works for any number of groups; with three or more it produces an overall test on the corresponding degrees of freedom.

What the log-rank test does and does not tell you

The test tells you whether the curves differ overall. It does not tell you by how much, at what times, or in which direction beyond the observed-versus-expected counts. Pair it with a Kaplan-Meier plot and with effect measures such as a hazard ratio from Cox regression or a restricted-mean-survival difference.

Weighted tests: the G-rho family

The standard log-rank test weights every event time equally, which makes it most sensitive to differences in the tails of the curves where events accumulate. Sometimes you care more about early differences. The Fleming-Harrington, or G-rho, family generalizes the log-rank test by weighting each event time by a power of the pooled survival probability.

The rho parameter controls this. With rho=0 you get the ordinary log-rank test. With rho=1 you get the Peto-Peto test, a Wilcoxon-type test that gives more weight to early event times and is therefore more sensitive to early separation of the curves.

This returns the same kind of TestResult as before, so you can read it the same way. Only the weighting has changed, so compare its statistic and p-value against the unweighted test above to see how much the emphasis on early event times matters here.

# Load data and build the response
lung = gw.load_dataset("lung", backend="polars")
y = gw.Surv.right(lung["time"], event=(lung["status"] == 2))

# Peto-Peto test: rho=1 weights early event times more heavily than the standard log-rank
gw.logrank_test(y, group=lung["sex"], rho=1)

TestResult(method='G-rho test (rho=1, gamma=0.0)', statistic=12.7142, df=1, p_value=0.0003629)

Choosing between them is a matter of what alternative you want to detect, and the choice should be made before looking at the results rather than by picking whichever gives the smaller p-value.

Choose the weighting in advance

Running several weighted tests and reporting the most significant one inflates the false positive rate. Decide on the test, and its rho, based on subject-matter reasoning before you see the data.

Robust handling of edge cases

Greenwood validates input data upfront, so weighted tests handle edge cases robustly. Times can be zero (events at baseline are valid in survival analysis), and weighted tests like Fleming-Harrington work correctly without silent data loss or dimension mismatches. This robustness extends to all censoring patterns and follows the validation rules documented in the Validation and reproducibility section.

Comparing more than two groups

The log-rank test extends to any number of groups, giving one overall test. We switch to the veteran dataset, which records a cell type with four levels.

# Load the veteran dataset and build the response
vet = gw.load_dataset("veteran", backend="polars")
y_vet = gw.Surv.right(vet["time"], event=vet["status"])

# Overall log-rank test across all four cell types
gw.logrank_test(y_vet, group=vet["celltype"])

TestResult(method='Log-rank test', statistic=25.4037, df=3, p_value=1.271e-05)

A small p-value here says the four cell types do not all share the same survival, but it does not say which ones differ. For that, test each pair and correct for the fact that you are running several tests at once. The pairwise_logrank_test() function does both.

# Test every pair of cell types; p_adjusted applies Holm multiple-comparison correction
gw.pairwise_logrank_test(y_vet, group=vet["celltype"], format="polars")

	group1 str	group2 str	statistic f64	p_value f64	p_adjusted f64
PolarsRows6Columns5
0	adeno	large	17.6693215293	2.62831687848e-05	0.000157699012709
1	adeno	smallcell	0.096843191971	0.755651328687	0.755651328687
2	adeno	squamous	12.0454836411	0.000519180149284	0.00259590074642
3	large	smallcell	9.37090414818	0.00220456751896	0.00661370255689
4	large	squamous	0.822593978656	0.364422837484	0.728845674968
5	smallcell	squamous	11.57367392	0.000668921225041	0.00267568490017

Each row is one pair, with its own statistic, the raw p_value, and a p_adjusted value that accounts for the multiple comparisons. The default correction is Holm’s but you can pass correction="bh" for Benjamini-Hochberg, "bonferroni", or "none". Read significance from the adjusted column, not the raw one.

Restricted mean survival time (RMST) comparisons

While the log-rank test tells you whether groups differ, it does not tell you by how much. The restricted mean survival time (RMST) is an intuitive effect measure: the expected survival time over a fixed time window tau, computed as the area under the survival curve. For example, 1-year RMST is the average survival time during the first year. Unlike hazard ratios, RMST is on the survival time scale, making it easy to interpret and explain to non-technical audiences.

Two-group RMST comparison

Compare RMST between two groups with confidence intervals and hypothesis tests using rmst_test(). Here we compare 1-year RMST between the two sexes in the lung dataset.

# Load data and build the response
lung = gw.load_dataset("lung", backend="polars")
y = gw.Surv.right(lung["time"], event=(lung["status"] == 2))

# Compare 1-year RMST between groups
result = gw.rmst_test(y, tau=365, group=lung["sex"])

# Display the RMST per group, their difference, and the p-value
result

RMSTResult(method='RMST difference (tau=365)', estimate=-55.9703, se=14.9581, 95% CI=[-85.2877, -26.6529], p_value=0.0001827)

The result shows the RMST for each group, their difference, standard error, 95% confidence interval, and a two-sided p-value testing equality. A positive difference means group 1 has higher RMST (longer average survival).

You can extract specific components:

# Extract the key effect-size components from the result object
print(f"RMST difference: {result.estimate:.1f} days")
print(f"95% CI: [{result.lower_ci:.1f}, {result.upper_ci:.1f}]")
print(f"P-value: {result.p_value:.4f}")

RMST difference: -56.0 days
95% CI: [-85.3, -26.7]
P-value: 0.0002

Alternative estimands

By default, rmst_test() returns the RMST difference. You can also compute the ratio or percentage difference using the estimand parameter.

Ratio of RMSTs (group 1 / group 2):

# Ratio estimand: group 1 RMST divided by group 2 RMST
result_ratio = gw.rmst_test(y, tau=365, group=lung["sex"], estimand="ratio")

# A value < 1 means group 1 has shorter average survival
result_ratio.estimate

0.8118425788462916

A ratio of 0.8 means group 1 has 80% of the RMST of group 2 (20% shorter average survival).

Percentage difference ((group 1 - group 2) / group 2 × 100):

# Percentage-difference estimand: (group 1 - group 2) / group 2 × 100
result_pct = gw.rmst_test(y, tau=365, group=lung["sex"], estimand="percentage_difference")

# Negative means group 1 has shorter average survival as a percentage
result_pct.estimate

-18.81574211537084

All three estimands use the delta method to construct standard errors and confidence intervals, and support flexible confidence levels via the conf_level= parameter.

Pairwise RMST comparisons

When you have three or more groups, use pairwise_rmst_test() to compare all pairs with multiple-comparison adjustment. Here we test which cell types in the veteran dataset differ in 1-year RMST.

# Load the veteran dataset and build the response
vet = gw.load_dataset("veteran", backend="polars")
y_vet = gw.Surv.right(vet["time"], event=vet["status"])

# Pairwise RMST differences with Holm adjustment
pairs = gw.pairwise_rmst_test(y_vet, tau=365, group=vet["celltype"])

# Display the full pairwise comparison table
pairs

	group1 str	group2 str	rmst1 f64	rmst2 f64	estimate f64	se f64	lower_ci f64	upper_ci f64	statistic f64	p_value f64	p_adjusted f64
PolarsRows6Columns11
0	adeno	large	65.5555555556	162.234567901	-96.6790123457	22.7443674012	-141.257153303	-52.1008713881	-4.25067932821	2.13123118245e-05	0.000127873870947
1	adeno	smallcell	65.5555555556	77.6352880658	-12.0797325103	17.1396326857	-45.6727952825	21.5133302619	-0.704783628202	0.480944898525	0.961889797051
2	adeno	squamous	65.5555555556	170.642372598	-105.086817043	26.1899192619	-156.418115554	-53.7555185313	-4.01249106543	6.00813536891e-05	0.000300406768445
3	large	smallcell	162.234567901	77.6352880658	84.5992798354	24.7760388636	36.0391359831	133.159423688	3.41456034603	0.000638850399879	0.00255540159952
4	large	squamous	162.234567901	170.642372598	-8.40780469693	31.7174551961	-70.5728745625	53.7572651687	-0.265084466738	0.790944393171	0.961889797051
5	smallcell	squamous	77.6352880658	170.642372598	-93.0070845323	27.9724458036	-147.832070867	-38.1820981978	-3.32495360561	0.000884332493635	0.0026529974809

Each row compares one pair of groups, showing the RMST estimates, difference, and both raw and adjusted p-values. The adjustment method (default: Holm) controls for multiple testing. Pass correction="bh" for Benjamini-Hochberg or "bonferroni" for Bonferroni correction.

Choosing the restriction time tau

The choice of tau is critical: RMST at tau=365 only includes follow-up through 1 year, while tau=1825 extends to 5 years. Choose tau based on clinical relevance before looking at the data. Common choices are clinically meaningful milestones (1 year, 5 years, etc.).

# Extend the horizon to 5 years and compare again
result_5y = gw.rmst_test(y, tau=1825, group=lung["sex"])

# Display the 5-year RMST difference
result_5y.estimate

-172.70616141861058

Stratified RMST comparisons

When your two groups differ in the distribution of a confounding variable, a plain RMST comparison mixes the groups’ true survival difference with the compositional imbalance. Stratification solves this: RMST is estimated separately within each stratum and the per-stratum differences are combined with inverse-variance weights before any test is computed.

Pass strata= to rmst_test(). Here we compare 1-year RMST between the two sexes while controlling for ECOG performance status (ph.ecog), which reflects disease severity and is correlated with both sex and prognosis.

import polars as pl

# Load the lung dataset
lung = gw.load_dataset("lung", backend="polars")

# Drop the small number of rows with missing ECOG score
lung_clean = lung.filter(pl.col("ph.ecog").is_not_null())

# Build the response from the cleaned frame
y = gw.Surv.right(lung_clean["time"], event=(lung_clean["status"] == 2))

# Stratified RMST: compare sexes while controlling for ECOG performance status
result_strat = gw.rmst_test(
    y,
    tau=365,
    group=lung_clean["sex"],
    strata=lung_clean["ph.ecog"],
)

# Display the pooled stratified result
result_strat

RMSTResult(method='Stratified RMST difference (tau=365)', estimate=-53.3006, se=14.3155, 95% CI=[-81.3585, -25.2428], p_value=0.0001967)

The result carries a stratified=True flag and the method string includes "Stratified" to make it clear which analysis was run. The estimate is the pooled RMST difference, and se reflects the combined uncertainty across all strata.

# Display the pooled estimate and its uncertainty
print(f"Pooled RMST difference: {result_strat.estimate:.1f} days")
print(f"95% CI: [{result_strat.lower_ci:.1f}, {result_strat.upper_ci:.1f}]")
print(f"P-value: {result_strat.p_value:.4f}")
print(f"Stratified: {result_strat.stratified}")

Pooled RMST difference: -53.3 days
95% CI: [-81.4, -25.2]
P-value: 0.0002
Stratified: True

How the pooling works. Within each stratum s, the RMST difference d_s = \text{RMST}_{1s} - \text{RMST}_{2s} and its variance \sigma_s^2 = \text{SE}_{1s}^2 + \text{SE}_{2s}^2 are computed from independent Kaplan-Meier curves. The combined estimate is:

\hat{\Delta} = \frac{\sum_s w_s\, d_s}{\sum_s w_s}, \quad w_s = \frac{1}{\sigma_s^2}, \quad \text{SE}(\hat{\Delta}) = \frac{1}{\sqrt{\sum_s w_s}}

Strata in which one group is absent are silently skipped, so unbalanced designs are handled without any special handling.

Stratified RMST vs. stratified log-rank

Both methods control for a confounder, but they answer different questions:

Stratified log-rank: is there a statistically significant difference in the survival curves? No magnitude.
Stratified RMST: by how many days does survival differ on average within the follow-up window, after adjusting for the stratification variable?

Use the log-rank test when you need a significance test under proportional hazards. Use stratified RMST when you want an interpretable effect size that does not require that assumption.

RMST vs. log-rank: complementary tools

log-rank test: answers whether curves differ overall (yes/no).
RMST comparison: answers by how much they differ and how long the difference persists.

Use both: run the log-rank test for significance, then estimate RMST difference for the magnitude and direction of the effect. They often agree, but RMST is easier to communicate to people unfamiliar with hazards.

Stratified tests

Sometimes you want to compare groups while controlling for a nuisance variable, so that the comparison is made within each level of that variable and then combined. A stratified log-rank does this. Here we compare the two treatments while stratifying by cell type, which removes any cell-type imbalance from the comparison.

# Load the veteran dataset and build the response
vet = gw.load_dataset("veteran", backend="polars")
y_vet = gw.Surv.right(vet["time"], event=vet["status"])

# Stratified log-rank: compare treatments while controlling for cell type
gw.logrank_test(y_vet, group=vet["trt"], strata=vet["celltype"])

TestResult(method='Stratified log-rank test', statistic=0.7017, df=1, p_value=0.4022)

The statistic is pooled across the strata, matching R’s survdiff with a strata() term. Stratification is the test-level counterpart to a stratified Cox model, where each stratum has its own baseline hazard.

Testing for trend in ordered groups

When your groups are naturally ordered as disease stages (I, II, III, IV), dose levels (low, medium, high), or educational attainment (primary, secondary, tertiary), a general log-rank test asks the right question but it doesn’t use the information about ordering. A trend test is more powerful because it looks specifically for a linear relationship across the ordered levels.

The trend test compares a linear contrast across groups to see if survival changes systematically with group ordering. It produces a one-degree-of-freedom chi-square test, which is more powerful than the multi-degree-of-freedom log-rank test when groups are truly ordered.

Default scoring

By default, groups are sorted and assigned scores 0, 1, 2, …, (k - 1). Here we test whether survival differs across cell types in the veteran dataset, treating them as an ordered grouping.

# Load the veteran dataset and build the response
vet = gw.load_dataset("veteran", backend="polars")
y = gw.Surv.right(vet["time"], event=vet["status"])

# Test linear trend across ordered cell types (default integer scores)
gw.trend_test(y, group=vet["celltype"])

TestResult(method='Linear trend test', statistic=2.3681, df=1, p_value=0.1238)

The result is a TestResult with df=1, indicating a one-degree-of-freedom test. The observed and expected event counts per group reveal which groups contribute most to the trend, and the p-value tells you if the trend is statistically significant.

Custom scoring

If the ordering is not uniform (for instance, some groups have very different prognoses) you can specify custom scores. Scores must be numeric and reflect your belief about the ordering.

# Load the veteran dataset and build the response
vet = gw.load_dataset("veteran", backend="polars")

y = gw.Surv.right(vet["time"], event=vet["status"])

# Custom scores: adeno and large get weight 1, smallcell and squamous get weight 3
scores = {"adeno": 1, "large": 1, "smallcell": 3, "squamous": 3}

# Apply the custom scoring to reflect domain knowledge about group ordering
gw.trend_test(y, group=vet["celltype"], scores=scores)

TestResult(method='Linear trend test', statistic=0.0263, df=1, p_value=0.8713)

Different scores can change the test statistic, reflecting your scientific hypothesis about how the groups are ordered.

Weighted trend tests

Like the log-rank test, the trend test can use Fleming-Harrington weights to emphasize early or late differences. With rho=1, the Peto-Peto variant emphasizes early event times, which is useful if you suspect early separation of the curves.

# Load the veteran dataset and build the response
vet = gw.load_dataset("veteran", backend="polars")
y = gw.Surv.right(vet["time"], event=vet["status"])

# Peto-Peto trend test: rho=1 weights early events more heavily
gw.trend_test(y, group=vet["celltype"], rho=1)

TestResult(method='G-rho trend test (rho=1, gamma=0.0)', statistic=0.6858, df=1, p_value=0.4076)

Stratified trend tests

You can also stratify a trend test to control for confounders while still testing for linear trend across ordered groups. The stratified test computes the trend separately in each stratum and then combines the results.

# Load the veteran dataset and build the response
vet = gw.load_dataset("veteran", backend="polars")
y = gw.Surv.right(vet["time"], event=vet["status"])

# Test trend in celltype while controlling for treatment group
gw.trend_test(y, group=vet["celltype"], strata=vet["trt"])

TestResult(method='Stratified linear trend test', statistic=1.7525, df=1, p_value=0.1856)

When to use trend tests vs. pairwise tests

use trend_test() when groups are ordered and you care about detecting a linear trend. It’s more powerful for this specific hypothesis.
use pairwise_logrank_test() when you want to know which specific pairs of unordered groups differ. It tests all pairwise differences, not a linear trend.
use logrank_test() for an overall test of differences without assuming any structure on the groups.

Study design: sample size and power

The log-rank test also answers a planning question asked before any data are collected: how large must a study be? Under proportional hazards, its power depends on the data only through the number of events, so Greenwood provides Schoenfeld’s calculators. To detect a given hazard ratio at a target power, the logrank_n_events() function returns the events required.

gw.logrank_n_events(hazard_ratio=0.5, power=0.9, alpha=0.05)

Running it the other way, logrank_power() gives the power a planned study would have for a given number of events, which is useful for checking whether an existing dataset can support a comparison.

gw.logrank_power(hazard_ratio=0.5, n_events=60)

0.7656462084964221

Events are not the same as subjects. To turn the required events into a sample size, divide by the probability that a subject has the event during the study, which you estimate from the expected survival and the planned follow-up. The logrank_sample_size() function does this.

gw.logrank_sample_size(hazard_ratio=0.5, prob_event=0.4, power=0.9)

What the calculators assume

These formulas assume proportional hazards, a two-group comparison, and (by default) a balanced, two-sided test at alpha=0.05. Use allocation for an unbalanced design and sides=1 for a one-sided test. A balanced design needs the fewest events.

From tests to models

A group comparison test is the right tool when you have one categorical grouping and want a yes-or-no answer about a difference. When you need to adjust for other variables, handle a continuous predictor, or estimate the size of an effect, you move to regression. The log-rank test has a close connection to the Cox model: the score test from a Cox model with a single group indicator is essentially the log-rank test, which is why the two usually agree.

Next steps

You can now describe and formally compare survival across groups.

Visualizing survival shows how to present grouped curves clearly.
Cox regression estimates hazard ratios and adjusts for multiple covariates.
Revisit Kaplan-Meier for restricted-mean-survival summaries, which give an interpretable effect size to accompany a significant test.