{party} を使った Survival CART 解析（Ref: http://rcommanderdeigakutoukeikaiseki.com/survival_CART.html ）。終端ノードの生存曲線の色や太さを変えられるようにしたテンプレート

SurvivalCART_with_party.R

# {party} を使った Survival CART 解析
# Ref: http://rcommanderdeigakutoukeikaiseki.com/survival_CART.html

# {party} の読み込み
if(!require("party")){
  install.packages("party", dep = TRUE, repos = "https://cloud.r-project.org/")
}
library("party")

# party::node_surv() を一部書き換えて、終端ノードの生存曲線の色と太さを指定できるようにする
# {party} のソース R/Plot.R より
# sv_col : 生存曲線の色　（デフォルトは黒 = 1）
# sv_lwd : 生存曲線の太さ（デフォルトは 1）
node_surv2 <- function(ctreeobj,
                       ylines = 2,
                       id = TRUE,
                       sv_col = 1, sv_lwd = 1, 
                       ...)
{
  survobj <- response(ctreeobj)[[1]]
  if (!("Surv" %in% class(survobj))) 
    stop(sQuote("ctreeobj"), " is not a survival tree")
  
  ### panel function for Kaplan-Meier curves in nodes
  rval <- function(node) {
    km <- party:::mysurvfit(survobj, weights = node$weights, ...)
    
    a <- party:::dostep(km$time, km$surv)
    
    yscale <- c(0,1)
    xscale <- c(0, max(survobj[,1]))
    
    top_vp <- viewport(layout = grid.layout(nrow = 2, ncol = 3,
                                            widths = unit(c(ylines, 1, 1), 
                                                          c("lines", "null", "lines")),  
                                            heights = unit(c(1, 1), c("lines", "null"))),
                       width = unit(1, "npc"), 
                       height = unit(1, "npc") - unit(2, "lines"),
                       name = paste("node_surv", node$nodeID, sep = ""))
    
    pushViewport(top_vp)
    grid.rect(gp = gpar(fill = "white", col = 0))
    
    ## main title
    top <- viewport(layout.pos.col=2, layout.pos.row=1)
    pushViewport(top)
    mainlab <- paste(ifelse(id, paste("Node", node$nodeID, "(n = "), "n = "),
                     sum(node$weights), ifelse(id, ")", ""), sep = "")
    grid.text(mainlab)
    popViewport()
    
    plot <- viewport(layout.pos.col=2, layout.pos.row=2,
                     xscale=xscale, yscale=yscale,
                     name = paste("node_surv", node$nodeID, "plot", 
                                  sep = ""))
    
    pushViewport(plot)
    # grid.lines(a$x/max(survobj[,1]), a$y)
    grid.lines(a$x/max(survobj[,1]), a$y, gp = gpar(col = sv_col, lwd = sv_lwd))
    grid.xaxis()
    grid.yaxis()
    grid.rect(gp = gpar(fill = "transparent"))
    upViewport(2)
  }
  
  return(rval)
}
class(node_surv2) <- "grapcon_generator"


# Survival CART 解析
#   data       : 解析するデータ
#       EZR / RCommander でアクティブデータセットを使う場合は、 data = get(activeDataSet())
#   time       : 観察期間
#   event      : イベント変数
#   event_code : イベント：JMP式なら0
#   var        : 説明変数（単変量）
#   plot       : グラフが不要なら FALSE
#   line_color : 生存曲線の色
#   line_width : 生存曲線の太さ
generate_tree <- function(data,
                          time       = "Time",
                          event      = "Censor",
                          event_code = 1,
                          var        = "Variable",
                          plot       = TRUE,
                          line_color = "red",
                          line_width = 1)
{
  # 作業用データを作る
  data_temp <- data.frame(cbind(
    time  = data[[time]],
    event = data[[event]] == event_code,
    var   = data[[var]]
  ))
  
  model <- ctree(
    Surv(time, event) ~ var,
    data = data_temp,
    controls = ctree_control(
      testtype  = "Univariate",
      # 少数すぎる終端Nodeができないよう調整
      minbucket = 20
    )
  )
  
  # グラフが必要なときは出力
  if(plot)
    plot(model,
         terminal_panel = node_surv2,
         tp_args = list(sv_col = line_color, sv_lwd = line_width))
  
  # Node番号をつけて返す
  return(
    transform(cbind(
      subset(data, select = names(data) %in% c(time, event, var)),
      Node = model@where
    ))
  )
}