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
|
#
# Data courtesy of Bercedis Peterson, Duke University.
# v4 of survreg fails due to 2 groups that have only 1 subject; the coef
# for them easily gets out of hand. In fact, this data set is my toughest
# test of the minimizer.
#
# A shrinkage model for this coefficient is therefore interesting
peterson <- data.frame(
scan('data.peterson', what=list(grp=0, time=0, status=0)))
fitp <- survreg(Surv(time, status) ~ factor(grp), peterson)
summary(fitp)
# Now a shrinkage model. Give the group coefficients
# about 1/2 the scale parameter of the original model, i.e., .18.
#
ffit <- survreg(Surv(time, status) ~ frailty(grp, theta=.1), peterson)
ffit
#
# Try 3 degrees of freedom Gaussian fit, since there are 6 groups.
# Compare them to the unconstrained ones. The frailty coefs are
# on a "sum to 0" constraint rather than "first coef=0", so
# some conversion is neccessary
#
ffit3 <- survreg(Surv(time, status) ~ frailty(grp, df=3, dist='gauss'),
peterson)
print(ffit3)
temp <- mean(c(0, fitp$coef[-1]))
temp2 <- c(fitp$coef[1] + temp, c(0,fitp$coef[-1]) - temp)
xx <- rbind(c(nrow(peterson), table(peterson$grp)),
temp2,
c(ffit3$coef, ffit3$frail))
dimnames(xx) <- list(c("N", "factor fit", "frailty fit"),
c("Intercept", paste("grp", 1:6)))
signif(xx,2)
#
# All but the first coef are shrunk towards zero.
#
rm(ffit, ffit3, temp, temp2, xx, fitp)
|
