[S] Labels wrong with lrm
Dear Jan, Thank you very much for your excellent description of the problem and the self-contained test code. This is a problem that I've been meaning to either document better or solve for some time. The root of the problem is with the builtin S-Plus terms.inner function:
attr(terms.inner(asthma ~ pol(age,kx) + smok),'variables')
expression(age, kx, smok) You can see that terms.inner inappropriately includes kx as an independent variable as it does not know that the first argument to pol is the special variable. When a constant replaces kx, all is well. As I am relying on the C code called by terms.inner to do the job, I don't have a ready solution. I would be happy if someone comes up with a solution. The all.vars function in the R language has the same limitation: all.vars(asthma ~ pol(age,kx) + smok) [1] "asthma" "age" "kx" "smok" Frank Harrell
Jan Brogger wrote:
I get funny behaviour from lrm when using transformations that require a parameter, and the parameter is a variable. It only happens when there are other variables in the model. The transformation must be the last term in the formula, or some labels will be wrong. For example, if the correct labels are "age" "age^2" and "smok", the labels will be "age", "age^2" and "kx" where "kx" is actually the parameter for the transformation, and it should be "smok". I have tried including the variable that contains the parameter for the transformation with datadist, but it doesn't help. Also tried reinstalling newest Hmisc, Design and starting new project directory. Platform: S+2000 Win98. It happens to both rcs and pol. With glm, everything is OK. The value of coefficients are the same. It isn't a major issue, but it might puzzle some and I thought I'd let you know. Solutions are: specify the parameters transformations directly, not in variables, or keep your transformations last in the formula. Any better solutions ? It'd be nice to be able to keep knot locations in a variable, so you can change it only one place in your code. This code reproduces the problem: library(Hmisc,T) library(Design,T) age <- rnorm(500,41.5,12) #age smok <- sample(0:2,500,rep=T) #smoking habit asthma <- sample(0:1,500,rep=T) #asthma mydat <- data.frame(age,smok,asthma) dd <- datadist(mydat) options(datadist="dd") kx <- 2 #this model gives the wrong labels lrm(asthma~pol(age,kx)+smok,data=mydat) #this model gives the right label lrm(asthma~smok+pol(age,kx),data=mydat) #If you don't use a variable, but a value, then everything is ok lrm(asthma~smok+pol(age,2),data=mydat) lrm(asthma~pol(age,2)+smok,data=mydat) #If you only have one term, then everything is OK lrm(asthma~pol(age,2),data=mydat) #It doesn't happen with GLM: glm(asthma~pol(age,2)+smok,data=mydat,family=binomial) #but it does with rcs: kx <- c(15,40,70) lrm(asthma~rcs(age,kx)+smok,data=mydat) lrm(asthma~smok+rcs(age,kx),data=mydat) Jan Brogger, University of Bergen, Norway --------------------------------------------------------------------- This message was distributed by s-news at lists.biostat.wustl.edu. To unsubscribe send e-mail to s-news-request at lists.biostat.wustl.edu with the BODY of the message: unsubscribe s-news
Frank E Harrell Jr Prof. of Biostatistics & Statistics Div. of Biostatistics & Epidem. Dept. of Health Evaluation Sciences U. Virginia School of Medicine http://hesweb1.med.virginia.edu/biostat -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._