plot_coxphSGD_estimates

library(coxphSGD)
library(survival)
library(reshape2)
set.seed(456)
x <- matrix(sample(0:1, size = 20000, replace = TRUE), ncol = 2)
head(x)
dCox <- dataCox(10^4, lambda = 3, rho = 2, x,
                beta = c(2,2), cens.rate = 5)
head(dCox)

batch_id <- sample(1:90, size = 10^4, replace = TRUE)
dCox_split <- split(dCox, batch_id)
results <-
  coxphSGD(formula     = Surv(time, status) ~ x.1+x.2,
           data        = dCox_split,
           epsilon     = 1e-5,
           learn.rates = function(x){1/(100*sqrt(x))},
           beta.zero   = c(0,0),
           max.iter    = 10*90)
           
coeff_by_iteration <-
  as.data.frame(
    do.call(
      rbind,
      results$coefficients
      )
    )
head(coeff_by_iteration)


coxph_loglik <- function(beta, formula, data) {
  coxph(formula, init=beta, control=list('iter.max'=0), data =data)$loglik[2]
}
coxph_loglik <- function(beta, formula, data) {
  coxph(formula, init=beta, control=list('iter.max'=0), data =data)$loglik[2]
}
calculate_outer_cox_3 <- function(dCox){
  ## contours
  outer_res <- outer(seq(0,4, length = 25),
                     seq(0,4, length = 25),
                     Vectorize( function(beta1,beta2){
                       coxph_loglik(beta=c(beta1,beta2), Surv(time, status)~x.1+x.2-1, dCox)
                     } )
  )
  outer_res_melted <- melt(outer_res)
  outer_res_melted$Var1 <- as.factor(outer_res_melted$Var1)
  levels(outer_res_melted$Var1) <- as.character(seq(0,4, length = 25))
  outer_res_melted$Var2 <- as.factor(outer_res_melted$Var2)
  levels(outer_res_melted$Var2) <- as.character(seq(0,4, length = 25))
  outer_res_melted$Var1 <- as.numeric(as.character(outer_res_melted$Var1))
  outer_res_melted$Var2 <- as.numeric(as.character(outer_res_melted$Var2))
  return(outer_res_melted)
}
calculate_outer_cox_3(dCox) -> outerCox
save(outerCox, file = 'dev/outerCox.rda')
#d2ggplot <- coeff_by_iteration
beta.zero <- c(0,0)
solution <- c(2,2)
library(ggplot2)
ggplot() +
stat_contour(aes(x=outerCox$Var1,
                 y=outerCox$Var2,
                 z=outerCox$value),
             bins = 40, alpha = 0.25) +
  geom_path(aes(coeff_by_iteration[['x.1']],
                coeff_by_iteration[['x.2']]),
                #group = d2ggplot$version,
                #colour = d2ggplot$version),
            size = 1) +
  theme_bw(base_size = 20) +
  theme(panel.border = element_blank(),
        legend.key = element_blank(),
        legend.position = "top") +
  scale_colour_brewer(palette="Dark2",
                      name = 'Algorithm \n & Steps') +
  geom_point(aes(x = beta.zero[1], y = beta.zero[2]),
             col = "black",
             size = 4, shape = 17) +
  geom_point(aes(x = solution[1], y = solution[2]),
             col = "black", size = 4, shape = 15) +
  geom_point(aes(x = summary(coxph(Surv(time, status) ~ x.1+x.2, data = dCox))$coeff[1,1],
                 y = summary(coxph(Surv(time, status) ~ x.1+x.2, data = dCox))$coeff[2,1]),
             col = "black", size = 4, shape = 13) +
  xlab("X1") +
  ylab("X2") -> p