# Calculating the Log Hazard Ratio with SE in R

In R, is there a predefined function that will give me the log hazard ratio and its standard error for a black male (as shown in the example below) given the output of coxph regression?

library(survival)
library(KMsurv)

#Kidney transplant data from Klein and Moeshberger. Massage data to make
#results look like those in book
data(kidtran)
data2 <- kidtran
data2$Gender <- "male" data2[data2$gender==2,7] <- "female"
data2$Race <- "white" data2[data2$race==2,8] <- "black"
data2$Gender <- as.factor(data2$Gender)
data2$Race <- as.factor(data2$Race)
data2$Race <- relevel(data2$Race,ref="white")

fit2 <- coxph(Surv(time,delta) ~ Gender * Race, data=data2)
summary(fit2)

#Relative log risk for a black male (reference white female) from
#page 252 in Klein and Moeshberger
(coef(fit2)[3] + coef(fit2)[2] + coef(fit2)[1])
sqrt(sum(diag(fit2$var)) + 2*fit2$var[2,1] + 2*fit2$var[3,1] + 2*fit2$var[3,2])

#Let's use predict
black.male <- data.frame(
Gender="male",
Race="black"
)

white.female <- data.frame(
Gender="female",
Race="white"
)

bm <- predict(fit2,newdata=black.male,se.fit=TRUE)

#bm in terms of original coefficients
coef(fit2)[3]*(1-fit2$means[3]) + coef(fit2)[2]*(1-fit2$means[2]) + coef(fit2)[1]*(1-fit2$means[1]) wf <- predict(fit2,newdata=white.female,se.fit=TRUE) #Relative log risk and se for a black male (reference white female) bm$fit - wf$fit sqrt(bm$se.fit^2 + wfse.fit^2)  - ## 1 Answer I do not think there exists a function like logHazardRatio("black males", "white females") which would return an estimate and a standard error. However, you can easily recover them as below. Answer prior to your edit Your Cox model is as follows, $$h_j(t) = h_0(t) \exp(\beta_1 \textrm{gender}_j + \beta_2 \textrm{race}_j + \beta_3 \textrm{gender}_j * \textrm{race}_j),$$ with $$\textrm{gender} = \left\{ \begin{array}{l} 0 \textrm{ if female} \\ 1 \textrm{ if male,} \\ \end{array} \right.$$ $$\textrm{race} = \left\{ \begin{array}{l} 0 \textrm{ if white} \\ 1 \textrm{ if black,} \\ \end{array} \right.$$ and we have $$\frac{\textrm{hazard in black males}}{\textrm{hazard in white females}} = \frac{h_0(t) \exp(\beta_1 + \beta_2 + \beta_3)}{h_0(t)}$$ so that $$\log \left( \frac{\textrm{hazard in black males}}{\textrm{hazard in white females}} \right) = \beta_1 + \beta_2 + \beta_3.$$ Your estimate coef(fit2)[3] + coef(fit2)[2] + coef(fit2)[1] is thus correct :-) To get a standard error, note that \begin{align*} \textrm{Var}(\hat{\beta}_1 + \hat{\beta}_2 + \hat{\beta}_3) & = \textrm{Var}(\hat{\beta}_1) + \textrm{Var}(\hat{\beta}_2) + \textrm{Var}(\hat{\beta}_3) \\ & + 2 \textrm{Cov}(\hat{\beta}_1, \hat{\beta}_2) + 2 \textrm{Cov}(\hat{\beta}_1, \hat{\beta}_3) + 2 \textrm{Cov}(\hat{\beta}_2, \hat{\beta}_3). \end{align*} Estimates for all of these terms can be found in the estimated covariance matrix of(\hat{\beta}_1, \hat{\beta}_2, \hat{\beta}_3)'$: > fit2$var
[,1]        [,2]        [,3]
[1,]  0.03941978  0.02576211 -0.03945987
[2,]  0.02576211  0.09731841 -0.09729282
[3,] -0.03945987 -0.09729282  0.18242564


For example, $\widehat{\textrm{Cov}}(\hat{\beta}_1, \hat{\beta}_3) = -0.03945987.$

-
Yes, but is there a function in R that does all of this for you? I thought that predict might do it, but I have been unable to get it to give me answers like the ones you provided. – Tom Nov 23 '12 at 20:36
I do not think there exists a function like 'logHazardRatio("black males", "white females")' which would return an estimate and a standard error. But it is not too tedious to go through calculations like in my answer... – ocram Nov 23 '12 at 20:38
I could not figure out how to use predict to do this calculation. Apparently predict gives the risk relative to a person who has mean values for all covariates, so perhaps I can just calculate the log hazard ratio for black male and subtract out the log hazard ratio for white female. This seems to work, but the SE is different than what you calculated above. – Tom Nov 24 '12 at 18:26
@Tom: 1) Can you say what 'bm' (in your edited question) is supposed to represent in terms of beta coefficients? 2) You computed Var('bm' - 'wf') assuming 'bm' and 'wf' being uncorrelated... which is probably false because your result differs from what we know to be true... – ocram Nov 24 '12 at 19:36
I added what bm is in terms of the original coefficients. I don't understand completely why R centers on the mean values of the factors for predicted values. Yes, it is unfortunate that the SE is wrong because otherwise this approach is giving me what I wanted. – Tom Nov 24 '12 at 21:42