I try to understand the "score residual" from cox ph model in R.
There is one reference site I have used for it.
Can someone help me to get the score residual manually (I mean, using formula) in R?
library(survival)
cph1 <- coxph(Surv(futime, fustat)~rx+age , data=ovarian)
residual <- residuals(cph1, type="score")
The score residuals in a Cox model are based on the "score function" that is solved to get the coefficient estimates and the martingale residuals that underlie the Cox model and its extensions. The "score residual" is actually a matrix of residuals, with an entry for each individual and each covariate coefficient estimated in the model.
For each covariate, you need its risk-weighted average $\bar x(t)$ among all cases at risk at each event time (used in solving the score function) and the (change in the) martingale residual for each individual at each event time. The martingale residuals in turn depend on the overall cumulative hazard estimated from the Cox model, $\hat \Lambda(t)$, for individual $i$: $\widehat M_i= \delta_i -r \hat \Lambda(t_i)$, where $\delta_i$ is the event indicator and $r$ is the relative risk.
Then the score residual for an individual and a covariate $x$ can be written:
$$\int_0^\infty (1- \bar x(t)) d \hat M_i(t), $$
which represents a sum over all event times. $d \hat M_i(t)$ is the change in the individual's martingale residual at each event time for which the individual is at risk.
Therneau and Grambsch work through a simple example in their Appendix E. Even that simple example with a single covariate and only 6 observations requires a fair amount of effort to work through in detail. For R coxph models the calculations are coded in C; you could examine the source code for more details.
