Predicting absolute risk using cox regression I am trying to use R to predict the absolute risk of developing adverse events in a cohort, and to compare that with the observed outcome. Should I use survreg or coxph to do this? Anyone kind enough to explain how to do this with R code? 
The mean follow up period of my cohort is only up to 6 years, so am I able to predict the absolute risk for each individual at the end of the follow up period (including both censored and non censored data)?
 A: So by "absolute risk" I'm going to assume you mean either the cumulative probability of an event at time t, or the hazard at time t.
In short, no, a Cox proportional hazards model doesn't really give you back that information - the model itself doesn't calculate the underlying hazard, just the relative difference in the hazard between covariate values. This rather nicely frees you from having to specify the underlying hazard of your outcome, which in many cases is unknown, not particularly of interest in the first place, or difficult to specify using a parametric model.
If you do want to estimate the underlying hazard function, you need to use parametric survival models, such as those used in survreg. There are a large number of tutorials online, including code, for survreg and parametric models.
This is one of my favorites, as it includes some theoretical treatment and a good bit of code. The 'Survival' package documentation is also a good place to start.
A: I agree with the above answers that Cox was not primarily meant for absolute risk calculations, but I believe it was
more for historical reasons (people focused on relative risks and not absolute), and don't see a legitimate reason why not to use it for predictions. There are works in prediction modelling that successfully use Cox already, but I think they prefer glmnet package for Cox family, with flexibility of regularization terms and hdnom package for survival probabilities (hdnom:::glmnet_survcurve).
Back to the question, going in cycles for survfit basic Cox model, I came up with this code to back out estimated individual probabilities of an event:
bh=basehaz(cox_model_fit, newdata = data, centered = FALSE) #baseline cumulative hazard function for all covariates set to 0 (centered = false
lps <- predict(cox_model_fit, data, type = "lp") #linear predictors, betas x covariates for each observation 

time_of_interest = 6 #if you are interested in survival at t=6

i_time_of_interest = match(1, round(bh$time,1) == time_of_interest, nomatch = -100) #annoyingly, no time argument in basehaz function, so find that time in the list (it builds at all censored and event times, so either choose a point like that or need to interpolate between available points)

event_prob_t = 1-exp(-bh[i_time_of_interest,1])^exp(lps) #calculate risk of event from the baseline hazard and linear predictors. returns a vector of probabilities for each person in the "data"

You can then check event_prob_t separately for cases and controls / censored and test the difference etc.
Or, use this function to calculate concordance statistic for survival data- this gives probability that two random observations from cases and controls have probabilities in the right order (lower chances of event for controls than for the cases). For binary outcome c-index is the same as ROC_AUC, so very common measure of how good a model is for classification/binary models.
 survConcordance(Surv(data$survtimed, data$outcome) ~ lps)$concordance
