# Biostatistics III in R

#### Code by Andreas Karlsson and Mark Clements

# Cox model for cause-specific mortality for melanoma (all stages)

Use Cox regression to model the cause-specific survival of patients with skin melanoma (including all stages).

You may have to install the required packages the first time you use them. You can install a package by `install.packages("package_of_interest")` for each package you require.

    library(biostat3)

Load the melanoma data and explore it.

    data(melanoma)
    melanoma <- transform(melanoma,
                          death_cancer = ifelse( status == "Dead: cancer", 1, 0))

## (a)

First fit the model with sex as the only explanatory variable. Does there appear to be a difference in survival between males and females?

    # Cox regression with time-since-entry as the timescale
    # Note that R uses the Efron method for approximating the likelihood in the presence of ties,
    # whereas Stata (and some other software) uses the Breslow method
    cox1 <- coxph(Surv(surv_mm, death_cancer) ~ sex, data=melanoma)
    summary(cox1)

## (b)

Is the effect of sex confounded by other factors (e.g. age, stage, subsite, period)? After controlling for potential confounders, does there still appear to a difference in survival between males and females?

    cox2 <- coxph(Surv(surv_mm, death_cancer) ~ sex + agegrp + stage + subsite + year8594, data=melanoma)
    summary(cox2)

## (c)

Consider the hypothesis that there exists a class of melanomas where female sex hormones play a large role in the etiology. These hormone related cancers are diagnosed primarily in women and are, on average, less aggressive (i.e., prognosis is good). If such a hypothesis were true we might expect the effect of sex to be modified by age at diagnosis (e.g., pre versus post menopausal). Test whether this is the case.

    cox3 <- coxph(Surv(surv_mm, death_cancer) ~ agegrp + agegrp * sex, data=melanoma)
    summary(cox3)

## (d)

Decide on a ‘most appropriate’ model for these data. Be sure to evaluate the proportional hazards assumption.

    cox4 <- coxph(Surv(surv_mm, death_cancer) ~ sex + year8594 + agegrp + subsite + stage, data=melanoma)
    summary(cox4)

    ## Test proportional hazards assumption
    cox.zph(cox4, transform="log")
    cox.zph(cox4, transform="identity") # Stata default

    ## Stratified Cox model; separate baseline hazard functions are fit for each strata.
    cox5 <- coxph(Surv(surv_mm, death_cancer) ~ sex + year8594 + agegrp + subsite + strata(stage), data=melanoma)
    summary(cox5)

    ## Test proportional hazards assumption
    cox.zph(cox5, transform="log")
    cox.zph(cox5, transform="identity")

    cox5 <- coxph(Surv(surv_mm, death_cancer) ~ sex * agegrp + year8594 + agegrp + subsite + strata(stage), data=melanoma)
    summary(cox5)
    anova(cox5)