  Treatment specific survival curves are an important tool to illustrate the
treatment effect in studies with time-to-event outcomes. In non-randomized
studies, unadjusted estimates can lead to biased depictions due to confounding.
Multiple methods to adjust survival curves for confounders exist. However, it
is currently unclear which method is the most appropriate in which situation.
Our goal is to compare forms of Inverse Probability of Treatment Weighting, the
G-Formula, Propensity Score Matching, Empirical Likelihood Estimation and
augmented estimators as well as their pseudo-values based counterparts in
different scenarios with a focus on their bias and goodness-of-fit. We provide
a short review of all methods and illustrate their usage by contrasting the
survival of smokers and non-smokers, using data from the German Epidemiological
Trial on Ankle-Brachial-Index. Subsequently, we compare the methods using a
Monte-Carlo simulation. We consider scenarios in which correctly or incorrectly
specified models for describing the treatment assignment and the time-to-event
outcome are used with varying sample sizes. The bias and goodness-of-fit is
determined by taking the entire survival curve into account. When used
properly, all methods showed no systematic bias in medium to large samples. Cox
regression based methods, however, showed systematic bias in small samples. The
goodness-of-fit varied greatly between different methods and scenarios. Methods
utilizing an outcome model were more efficient than other techniques, while
augmented estimators using an additional treatment assignment model were
unbiased when either model was correct with a goodness-of-fit comparable to
other methods. These doubly-robust methods have important advantages in every
considered scenario.
