library(mfp) # CRAN - mfp

data(GBSG)

# htreatは層化変数 / tumsize,posnodal,esmはFP変換

# S={-2,-1,-0.5,0,0.5,1,2,3}が冪乗の候補(ただし，0は対数変換)

# CoxPHを当てはめ（変数選択も同時に行う）

f <- mfp(Surv(rfst,cens)~strata(htreat)+age+fp(tumsize)+

fp(posnodal)+fp(prm)+fp(esm)+menostat+tumgrad,family=cox,

data=GBSG,select=0.05,verbose=T)

f$fptable # 変数選択の結果

f$fit # 最終モデルの結果

plot(survfit(f$fit),col=c("red","green"))

legend("topright",legend=c("htreat = no","htreat = yes"),lty=1,col=c("red","green"))