# Linear predictor from coefficients of Cox PH model

I need to calculate the linear predictor of a Cox PH model by hand.

I can get continuous and binary variables to match the output of predict.coxph (specifying 'lp') but I can't seem to figure out how to calculate it for categorical variables with more than 2 levels.

My aim is to assess calibration of a published model in my own data-I only have coefficients so need to be able to do this by hand.

This post on stackoverflow describes how to calculate for continuous variables...

Any advice would be appreciated! Thanks

R example...

URL   <- "http://socserv.mcmaster.ca/jfox/Books/Companion/data/Rossi.txt"

summary(Rossi[,c("week", "arrest", "fin")])
#      week           arrest        fin
# Min.   : 1.00   Min.   :0.0000   no :216
# 1st Qu.:50.00   1st Qu.:0.0000   yes:216
# Median :52.00   Median :0.0000
# Mean   :45.85   Mean   :0.2639
# 3rd Qu.:52.00   3rd Qu.:1.0000
# Max.   :52.00   Max.   :1.0000

library(survival)

#for binary variable
fitCPH <- coxph(Surv(week, arrest) ~ fin, data=Rossi)    #Cox-PH model
(coefCPH <- coef(fitCPH))                               # estimated coefficients
#   finyes
#-0.3690692

#[1]  0.1845346  0.1845346  0.1845346 -0.1845346  0.1845346  0.1845346

head(((as.numeric(Rossi$$fin) - 1) - mean(as.numeric(Rossi$$fin) - 1)) * coef(fitCPH))
#[1]  0.1845346  0.1845346  0.1845346 -0.1845346  0.1845346  0.1845346

#for categorical variable
set.seed(170981)
Rossi$categorical.example <- as.factor(sample(1:3,nrow(Rossi),replace = TRUE)) summary(Rossi[,c("week", "arrest", "categorical.example")]) # week arrest categorical.example # Min. : 1.00 Min. :0.0000 1:156 # 1st Qu.:50.00 1st Qu.:0.0000 2:138 # Median :52.00 Median :0.0000 3:138 # Mean :45.85 Mean :0.2639 # 3rd Qu.:52.00 3rd Qu.:1.0000 # Max. :52.00 Max. :1.0000 fitCPH2 <- coxph(Surv(week, arrest) ~ categorical.example, data=Rossi) #Cox-PH model (coefCPH2 <- coef(fitCPH2)) #categorical.example2 categorical.example3 # -0.181790 -0.103019 head(predict(fitCPH2,type="lp")) #[1] 0.09098066 -0.01203832 -0.01203832 0.09098066 -0.09080938 -0.09080938 #How to calculate by hand??  • It seems that neither riskRegression:::coxLP.coxph nor predict(fitCPH,type="lp") center the categorical and binary variables anymore. Any confirmation? Commented Jan 12, 2023 at 23:51 ## 2 Answers You need to know the how the dummy variables are generated when you fit the model. Otherwise, it is impossible to get linear predictor. For example, suppose categorical variable X has 3 levels. Generally, 2 dummy variables will be generated. One coding scheme is as following:  dummy variables ------------------------------------------- X level X1 x2 ------------------------------------------- 1 0 0 2 1 0 3 0 1 --------------------------------------------  Then you can get the linear predictor for subjects in different X level by apply different regression coefficients. For your example, calculate a shift first. Total number of objects = 432. Among them 138 with categorical.example = 2, and 138 with categorical.example =3. So averages are 138/432 and 138/432. The shift is $$s = -(-0.181790 \times 138/432 + -0.103019 \times 138/432) = 0.09098066$$. Then for each individual, the linear predictor is $$X_1 \beta_1 + X_2\beta_2 +s$$ For categorical.example = 1, we have $$X_1 = 0$$ and $$X_2=0$$, so the linear predictor = $$s = 0.09098066$$ For categorical.example = 2, we have $$X_1 = 1$$ and $$X_2=0$$, so the linear predictor = $$\beta_1 + s = -0.181790 + 0.09098066 = -0.9080938$$ For categorical.example = 2, we have $$X_1 = 1$$ and $$X_2=0$$, so the linear predictor = $$\beta_2 + s = -0.103019 + 0.09098066 = -0.01203832$$ • Thank you-could you explain this using the R example added to original question? Thanks in advance for your help! Commented Dec 5, 2018 at 10:34 Thanks @user158565 for your help, I figured it out from the dummy variable clue you posted! Code below if anyone interested! dv1 <- as.numeric(Rossi$categorical.example)-1    #make 0,1,2 rather than 1,2,3
dv1[dv1==2] <- 0

dv2 <- as.numeric(Rossi\$categorical.example)-2    #make -1,0,1 rather than 1,2,3
dv2[dv2==-1] <- 0

meandv1  <- mean(dv1)
meandv2  <- mean(dv2)