The R package coxphMIC implements a sparse estimation method for Cox proportional hazards models by minimizing the approximated information criteria (BIC by default). To load the pacakge,
require(coxphMIC)
## Loading required package: coxphMIC
pbc
DataTo illsutrate, we consider the well-known pbc
data from the survival package.
library(survival); data(pbc);
dat <- pbc; dim(dat);
## [1] 418 20
A couple preprocessing step is made, including listwise deletion for missing values:
dat$status <- ifelse(pbc$status==2, 1, 0)
dat$sex <- ifelse(pbc$sex=="f", 1, 0)
dat <- na.omit(dat);
dim(dat); head(dat)
## [1] 276 20
## id time status trt age sex ascites hepato spiders edema bili chol
## 1 1 400 1 1 58.76523 1 1 1 1 1.0 14.5 261
## 2 2 4500 0 1 56.44627 1 0 1 1 0.0 1.1 302
## 3 3 1012 1 1 70.07255 0 0 0 0 0.5 1.4 176
## 4 4 1925 1 1 54.74059 1 0 1 1 0.5 1.8 244
## 5 5 1504 0 2 38.10541 1 0 1 1 0.0 3.4 279
## 7 7 1832 0 2 55.53457 1 0 1 0 0.0 1.0 322
## albumin copper alk.phos ast trig platelet protime stage
## 1 2.60 156 1718.0 137.95 172 190 12.2 4
## 2 4.14 54 7394.8 113.52 88 221 10.6 3
## 3 3.48 210 516.0 96.10 55 151 12.0 4
## 4 2.54 64 6121.8 60.63 92 183 10.3 4
## 5 3.53 143 671.0 113.15 72 136 10.9 3
## 7 4.09 52 824.0 60.45 213 204 9.7 3
The function coxphMIC()
does the sparse estimation and returns a S3 class object.
fit.mic <- coxphMIC(formula=Surv(time, status)~.-id, data=dat, method="BIC", scale.x=TRUE)
names(fit.mic)
## [1] "opt.global" "opt.local" "min.Q" "gamma" "beta"
## [6] "VCOV.gamma" "se.gamma" "se.beta" "BIC" "result"
## [11] "call"
The print.coxphMIC()
provides a summary table that tabulates β estiamtes and the reparameterized γ estiamtes.
print(fit.mic)
##
## Call:
## coxphMIC(formula = Surv(time, status) ~ . - id, data = dat, method = "BIC",
## scale.x = TRUE)
##
## Table of Estimated Coefficients via MIC:
##
## beta0 gamma se.gamma LB UB z.stat p.value
## trt -0.0622 0.0000 0.1071 -0.21000 0.21000 0.0000 1.0000
## age 0.3041 0.3309 0.1219 0.09191 0.56986 2.7138 0.0067
## sex -0.1204 0.0000 0.1086 -0.21285 0.21280 -0.0002 0.9998
## ascites 0.0224 0.0000 0.0991 -0.19423 0.19423 0.0000 1.0000
## hepato 0.0128 0.0000 0.1259 -0.24684 0.24684 0.0000 1.0000
## spiders 0.0460 0.0000 0.1118 -0.21908 0.21907 -0.0001 1.0000
## edema 0.2733 0.2224 0.1066 0.01345 0.43139 2.0861 0.0370
## bili 0.3681 0.3909 0.1142 0.16713 0.61472 3.4237 0.0006
## chol 0.1155 0.0000 0.1181 -0.23144 0.23148 0.0002 0.9999
## albumin -0.2999 -0.2901 0.1248 -0.53478 -0.04543 -2.3239 0.0201
## copper 0.2198 0.2518 0.1050 0.04605 0.45757 2.3986 0.0165
## alk.phos 0.0022 0.0000 0.0837 -0.16411 0.16412 0.0000 1.0000
## ast 0.2308 0.2484 0.1128 0.02732 0.46938 2.2023 0.0276
## trig -0.0637 0.0000 0.0858 -0.16820 0.16820 0.0000 1.0000
## platelet 0.0840 0.0000 0.1129 -0.22134 0.22134 0.0000 1.0000
## protime 0.2344 0.2293 0.1046 0.02424 0.43428 2.1917 0.0284
## stage 0.3881 0.3692 0.1476 0.07983 0.65856 2.5007 0.0124
## beta.MIC se.beta.MIC
## trt 0.0000 NA
## age 0.3309 0.1074
## sex 0.0000 NA
## ascites 0.0000 NA
## hepato 0.0000 NA
## spiders 0.0000 NA
## edema 0.2224 0.0939
## bili 0.3909 0.0890
## chol 0.0000 NA
## albumin -0.2901 0.1103
## copper 0.2518 0.0868
## alk.phos 0.0000 NA
## ast 0.2484 0.1025
## trig 0.0000 NA
## platelet 0.0000 NA
## protime 0.2293 0.1022
## stage 0.3692 0.1243
The plot.coxphMIC()
gives a bar-plot fo the fitting info.
plot(fit.mic)