[R] survreg 3-way interaction

tdenes at cogpsyphy.hu tdenes at cogpsyphy.hu
Fri Jan 28 22:30:20 CET 2011


The fire devoured my laptop, but the PC is still working... :)
Sorry for the lazyness, you are totally right.

Here goes a reproducible example, which resembles the main features of our
dataset:

# three experimental factors, 30 subjects
expgrid <- expand.grid(F1=0:1,F2=0:1,F3=0:1,id=1:30)

# set up the dataset (make it unbalanced)
inpdata <- expgrid[sample(nrow(expgrid),nrow(expgrid)*3,replace=T),]
inpdata <- inpdata[order(inpdata$id),]

# define coefficients
# fixed effects
coefs <- c(5,rnorm(7))
# random effects
subj.coefs <- rnorm(30)

# compute fixed effects
modmat <- model.matrix(logtimes~F1*F2*F3,
	data=data.frame(logtimes=0,inpdata))
fix.eff <- rep(modmat%*%coefs)

# compute logtimes (fixed effects + random effects + noise)
logtimes <- fix.eff+rep(subj.coefs,table(inpdata$id))+
	rnorm(nrow(inpdata),0,0.5)

# random censoring
cens <- runif(nrow(inpdata))<0.2
logtimes[cens] <- logtimes[cens]-rnorm(sum(cens))

# finalize dataset
inpdata$status <- !cens
inpdata$times <- exp(logtimes)
inpdata$id <- factor(inpdata$id)

# run survreg with fixed df in frailty
# (method="aic" leads to error, don't know why)
library(survival)
survreg3 <- survreg(Surv(times,status) ~
	F1*F2*F3 + frailty.gamma(id,df=29), data=inpdata,
	dist="lognormal")
# Error in survpenal.fit(X, Y, weights, offset, init = init, controlvals =
control,  :
#  Invalid pcols or pattr arg

# note, that without the frailty term or without the tree-way interaction
term, it runs smoothly


# the same with gamlss
library(gamlss.cens)
gamlss3 <- gamlss(Surv(times,status) ~
	F1*F2*F3 + random(id,df=29), data=inpdata,
	family=cens(LOGNO))






>  > I was wondering why survreg (in survival package) can not handle
>> three-way interactions. I have an AFT .....
>
>   You have given us no data to diagnose your problem. What do you mean
> by "cannot handle" -- does the package print a message "no 3 way
> interactions", gives wrong answers, your laptop catches on fire when you
> run it, ....?
>   Also, make sure you read help(frailty).  The following is a key line
> "Use of gamma/ml or gaussian/reml with survreg does not lead to valid
> results."
>
> Terry Therneau
>
>
>



More information about the R-help mailing list