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
>
>
>
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
