forked from r-forge/surrosurv
-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathconvergence.R
110 lines (100 loc) · 3.76 KB
/
convergence.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
97
98
99
100
101
102
103
104
105
106
107
108
109
110
# ---------------------------------------------------------------------------- #
convals <- function(x, ...) UseMethod('convals')
# ---------------------------------------------------------------------------- #
# ---------------------------------------------------------------------------- #
convals.surrosurv <- function(x, ...) {
models <- c('Clayton', 'Plackett', 'Hougaard', 'Poisson')
models <-
models[sapply(models, function(y)
any(grepl(y, names(x))))]
copulas <- function(cop) {
f <- function() {
res <- t(sapply(x[[cop]], function(y) {
vals <- c(maxSgrad = tryCatch(
max(abs(with(
attributes(y$optimxRES), solve(hessian, grad)
))),
error = function(x)
return(Inf)
),
minHev = tryCatch(
min(eigen(attr(
y$optimxRES, 'hessian'
))$values),
error = function(x)
return(-1)
),
minREev = tryCatch(
min(eigen(y$VarCor2)$values),
error = function(x)
return(-1)
))
return(vals)
}))
rownames(res) <- paste(cop, rownames(res))
return(res)
}
return(f)
}
Clayton = copulas('Clayton')
Plackett = copulas('Plackett')
Hougaard = copulas('Hougaard')
Poisson = function() {
which_models <- x[grepl('Poisson', names(x))]
convres <- t(sapply(which_models,
function(y)
tryCatch({
vals <- c(maxSgrad = tryCatch(
max(abs(with(
y$optinfo$derivs, solve(Hessian, gradient)
))),
error = function(x)
return(Inf)
),
minHev = tryCatch(
min(eigen(y$optinfo$derivs$Hessian)$values),
error = function(x)
return(-1)
),
minREev = tryCatch(
min(eigen(y$VarCor$trialref)$values),
error = function(x)
return(-1)
))
return(vals)
}, error = function(x)
return(c(
Inf, -1, -1
)))))
names(convres) <- names(which_models)
return(convres)
}
fitRES <-
do.call(rbind, lapply(models, function(x)
eval(call(paste(
x
)))))
fitRES[grepl('unadj|PoissonI', rownames(fitRES)), 'minREev'] <- NA
class(fitRES) <- 'conv'
return(fitRES)
}
# ---------------------------------------------------------------------------- #
# ---------------------------------------------------------------------------- #
convergence <- function(x, kkttol = 1e-2, kkt2tol = 1e-8, ...)
UseMethod('convergence')
# ---------------------------------------------------------------------------- #
# ---------------------------------------------------------------------------- #
convergence.surrosurv <- function(x, kkttol = 1e-2, kkt2tol = 1e-8, ...) {
ConvRES <- convals(x)
checkConvRES <- cbind(ConvRES[, 1, drop = FALSE] <= kkttol,
ConvRES[, 2:3] > kkt2tol)
class(checkConvRES) <- 'conv'
return(checkConvRES)
}
# ---------------------------------------------------------------------------- #
# ---------------------------------------------------------------------------- #
print.conv <- function(x, na.print = '---', ...) {
class(x) <- 'matrix'
print(x, na.print = na.print, ...)
}
# ---------------------------------------------------------------------------- #