-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathComparisionFLX.R
96 lines (66 loc) · 3 KB
/
ComparisionFLX.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
### This function compares feature selection and C-Index recovery and prediction ###
### Let's see what comes up
ComparisionFLX = function(){
gr.km <- kmeans(Y, F, nstart =10)
gr.km.rand <- adjustedRandIndex(c.true,as.factor(gr.km$cluster))
#################Flexmix Clustering ################################################################
##############################################################################
data <- data.frame(y =time, x = Y)
## The cross validation folds for choosing lambda
fo <- sample(rep(seq(10), length = nrow(data)))
gr.flx <- flexmix(y ~ ., data = data, k = F, cluster = gr.km$cluster, model = FLXMRglmnet(foldid = fo, adaptive= FALSE, family = c("gaussian")), control = list(iter.max = 500))
recovRandIndex.flx <<- as.numeric(adjustedRandIndex(c.true,as.factor(clusters(gr.flx))))
linear.recov.flx <- as.numeric(unlist(predict(gr.flx, newdata = data, aggregate = TRUE)))
recovCIndex.flx <<- as.numeric(survConcordance(smod ~ exp(-linear.recov.flx))[1])
############### Penalized FlexMix #######################################################
################################################################################
## Fit a AFT model with FLXmix clustering
linear.flx <- c(0)
beta.flx <- matrix(0, nrow = D, ncol = F)
for ( q in 1:F){
ind <- which(clusters(gr.flx) == q)
L= length(ind)
time.tmp <- time[ind]
censoring.tmp <- censoring[ind]
Y.tmp <- Y[ind,]
reg <- cv.glmnet(x = Y.tmp, y = time.tmp, family = "gaussian")
coeff.pred <- coef(object =reg, newx = Y.tmp, s= "lambda.min")
rel.coeff <- coeff.pred[2:(D+1)]
beta.flx[1:D,q] <- ((rel.coeff != 0)+0)
linear.flx[ind] <- predict(object = reg, newx = Y.tmp, s = "lambda.min")
}
recovCIndex.flx.aft <<- as.numeric(survConcordance(smod ~ exp(-linear.flx))[1])
### Use flexmix clustering
data.new <- data.frame(x = Y.new)
linear.pred.flx <- as.numeric(unlist(predict(gr.flx, newdata = data.new, aggregate = TRUE)))
predCIndex.flx <<- as.numeric(survConcordance(smod.new ~ exp(-linear.pred.flx))[1])
### Use FLexmix to make predictions on future probabilities ####
# cl.new <- posterior(gr.flx, newdata = data.new)
# cl.pred.flx <- c(0)
#
# for ( i in 1: dim(Y.new)[1]){
# cl.pred.flx[i] <- which.max(cl.new[i,])
#
# }
# predRandIndex.flx <<- adjustedRandIndex(c.true.new,as.factor(cl.pred.flx))
#
#
# ### Use PAFT now to build the corresponding cluster specific PAFT models
# pre.flx <- c(0)
#
# for ( q in 1:F){
# ind <- which(clusters(gr.flx) == q)
# ind.new <- which(cl.pred.flx == q)
#
# time.tmp <- time[ind]
#
# Y.tmp <- Y[ind,]
# Y.tmp.new <- Y.new[ind.new,]
#
# reg <- cv.glmnet(x = Y.tmp, y = time.tmp, family = "gaussian")
# pre.flx[ind.new] <- predict(object = reg, newx = Y.tmp.new, s = "lambda.min")
# }
# predCIndex.flx.aft <<- as.numeric(survConcordance(smod.new ~ exp(-pre.flx))[1])
#
#
}