Re: [Rd] [R] predict()
On Apr 14, 2011, at 16:52 , Ivo Shterev wrote:
> Dear Dr. Therneau,
>
> Thank you for your response. Just to point out that we didn't experience any
> problems with the lm() function under R version 2.12.2 (2011-02-25):
I couldn't reproduce it from Terry's description either, but there _is_ an
issue which parallels the coxph one: In model.frame.lm, we have
function (formula, ...)
{
dots <- list(...)
nargs <- dots[match(c("data", "na.action", "subset"), names(dots),
0)]
if (length(nargs) || is.null(formula$model)) {
fcall <- formula$call
fcall$method <- "model.frame"
fcall[[1L]] <- as.name("lm")
fcall[names(nargs)] <- nargs
env <- environment(formula$terms)
if (is.null(env))
env <- parent.frame()
eval(fcall, env, parent.frame())
}
else formula$model
}
So under some circumstances, we evaluate formula$call (where "formula" is
actually an "lm" object --- and age-old design infelicity) with modified
arguments. This evaluation takes place in the environment of the model formula,
nested in the parent frame, but neither necessarily contains the arguments to
formula$call:
> t3
function(dat,fm)
{
model.frame(lm(fm,data=dat),subset=1:5)
}
> testdat
otime event x
1 0.84345726 0 -0.4456620
2 0.57661027 0 1.2240818
3 1.32905487 0 0.3598138
4 0.03157736 0 0.4007715
5 0.05621098 0 0.1106827
6 0.31650122 1 -0.5558411
7 0.31422729 1 1.7869131
8 0.14526680 1 0.4978505
9 2.72623646 1 -1.9666172
10 0.02915345 1 0.7013559
> fm2
otime ~ x
> t3(testdat,fm2)
Error in model.frame(formula = fm, data = dat, subset = 1:5, drop.unused.levels
= TRUE) :
object 'fm' not found
It looks more than a bit iffy to try and set this right. One suggestion could
be to explicitly replace the formula part of fcall with the formula contained
in the model object. Or, the original evaluation frame should be kept with the
model object and used rather than parent.frame().
This should probably move to r-devel.
>
>> set.seed(123)
>> testdat=data.frame(y=rexp(10),event=rep(0:1,each=5),x=rnorm(10))
>> testfm=as.formula('y~x')
>>
>> testfun=function(dat,fm)
> + {
> + predict(lm(fm,data=dat),newdata=dat)
> + }
>> testfun(testdat,testfm)
> 1 2 3 4 5 6
> 1.00414971 0.07061335 0.55381658 0.53091759 0.69310319 1.06574974
> 7 8 9 10
> -0.24405980 0.47664172 1.85449993 0.36286396
>> sessionInfo()
> R version 2.12.2 (2011-02-25)
> Platform: x86_64-pc-linux-gnu (64-bit)
>
> locale:
> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
> [3] LC_TIME=en_US.UTF-8LC_COLLATE=en_US.UTF-8
> [5] LC_MONETARY=C LC_MESSAGES=en_US.UTF-8
> [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
> [9] LC_ADDRESS=C LC_TELEPHONE=C
> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
>
> attached base packages:
> [1] stats graphics grDevices utils datasets methods base
>
> loaded via a namespace (and not attached):
> [1] tools_2.12.2
>>
>
>
> Regards
> Ivo
>
>
>
> --- On Thu, 4/14/11, Terry Therneau wrote:
>
>> From: Terry Therneau
>> Subject: Re: predict()
>> To: "Ivo Shterev"
>> Cc: r-h...@r-project.org
>> Date: Thursday, April 14, 2011, 2:59 PM
>> --- begin included message ---
>> I am experimenting with the function predict() in two
>> versions of R and
>> the R extension package "survival".
>>
>> library(survival)
>>
>> set.seed(123)
>> testdat=data.frame(otime=rexp(10),event=rep(0:1,each=5),x=rnorm(10))
>> testfm=as.formula('Surv(otime,event)~x')
>>
>> testfun=function(dat,fm)
>> {
>>
>> predict(coxph(fm,data=dat),type='lp',newdata=dat)
>> }
>>
>> -- end inclusion
>>
>> The question was: this works under survival 2.35-8, but
>> not survival
>> 2.36-5
>>
>> Answer: The predict and underlying model.frame functions
>> for coxph were
>> brought into line with lm and other models. The major
>> advantage is that
>> I now deal with factors and data dependent predictors like
>> ns()
>> correctly.
>> You've shown a disadvantage of which I was not
>> aware. Using your
>> example but replacing coxph() by lm() with otime ~x as the
>> model I get a
>> similar failure. I'd like to ask a wider audience of
>> R-devel since it
>> is bigger than coxph.
>>
>> Terry Therneau
>>
>>
>>
>
> __
> r-h...@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.
--
Peter Dalgaard
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Email: pd@cbs.dk Priv: pda...@gmail.com
Re: [Rd] predict()
On Fri, 2011-04-15 at 09:10 +0200, peter dalgaard wrote: > I couldn't reproduce it from Terry's description either, but there > _is_ an issue which parallels the I'll start with an apology: as a first guess to understand the problem with predict.coxph I tried something parallel to Ivo's example using lm, during a 5 minute slice of time between meetings. I got a similar error message, posted a query based on that, and ran off. The real root of my error message, however, was a simple typographical mistake. My apologies for the premature note, whose error was gently pointed out by all three of you. I should have waited till I had more time. Now for the real problem, which is an oversight in my design. When updating the predict.coxph, residuals.coxph, survfit.coxph and survexp.coxph functions to more modern processing of the newdata argument (proper handling of factors, etc), I made a decision to have all of them call model.frame.coxph and model.matrix.coxph. Model.matrix in particular has some special steps, and it would be better to have only one point of modification. However, this introduces one more level of indirection predict >- model.frame.coxph -> model.frame and the parent.frame() reference in model.frame.coxph now points to predict when we actually want it to point to the parent of predict. In predict.lm the parent.frame() argument lives in the predict.lm code. I see three paths to correction: 1. Make model.frame.coxph smarter. Peter's example seems to show that this is hard. 2. Change the line in predict.coxph (and the others) mf <- model.frame(object) to some form of eval() that causes the parent.frame() operation to reach back properly. I'm not yet sure what variant of eval or do.call or ... this is; the effect I'm looking for is similar to what NextMethod() does. 3. Change my call to mf <- model.frame(object$terms, ... mimic the call in lm This may be the simplest, since model.frame.coxph does not do anything special. Last, I need to check whether any of the above issues affect model.matrix calls. Any comments on which is the better path would be welcomed. Terry T. __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] predict()
An addition to my prior post: my option 3 is not as attractive as I thought. In several cases predict.coxph needs to reconstruct the original data frame, for things that were not saved in the coxph object. The half dozen lines to redo the orignal call with the same data, na.action, etc options is simple, but a nuisance to reproduce in multiple spots. Terry __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] DESCRIPTION file and Rd examples
It would probably be helpful if you could point us to the source version
of your package so that we can take a look what happens.
Which version of R, and the other required packages are you talking about?
Uwe Ligges
On 15.04.2011 05:00, Dario Strbenac wrote:
I have a confusing error from R CMD check that I don't get when running the
example manually by hand.
In the \examples section of an Rd file, I create a GRanges object, then I call
a function with the GRanges object, whose first 2 lines are
require(GenomicRanges)
annoDF<- as.data.frame(anno) # anno is the GRanges object.
and that second line gives:
Error in as.data.frame.default(anno) :
cannot coerce class 'structure("GRanges", package = "GenomicRanges")' into a
data.frame
Calls: annoGR2DF ... annoGR2DF -> .local -> as.data.frame ->
as.data.frame.default
I have GenomicRanges listed in my Imports: field, and IRanges in the Suggests:
of the DESCRIPTION file (it's require()d elsewhere). I'm trying to avoid
putting packages in Depends: , so my package loads fast. Any tips of what I'm
not understanding properly ?
Thanks.
--
Dario Strbenac
Research Assistant
Cancer Epigenetics
Garvan Institute of Medical Research
Darlinghurst NSW 2010
Australia
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] DESCRIPTION file and Rd examples
On Apr 14, 2011, at 11:00 PM, Dario Strbenac wrote:
> I have a confusing error from R CMD check that I don't get when running the
> example manually by hand.
>
> In the \examples section of an Rd file, I create a GRanges object, then I
> call a function with the GRanges object, whose first 2 lines are
>
>require(GenomicRanges)
require() is doesn't guarantee that the package will load, so I think what you
meant to write was more
if (require(GenomicRanges, quietly=TRUE)) {
...
>annoDF <- as.data.frame(anno) # anno is the GRanges object.
>
> and that second line gives:
>
> Error in as.data.frame.default(anno) :
> cannot coerce class 'structure("GRanges", package = "GenomicRanges")' into a
> data.frame
> Calls: annoGR2DF ... annoGR2DF -> .local -> as.data.frame ->
> as.data.frame.default
>
> I have GenomicRanges listed in my Imports: field, and IRanges in the
> Suggests: of the DESCRIPTION file (it's require()d elsewhere). I'm trying to
> avoid putting packages in Depends: , so my package loads fast. Any tips of
> what I'm not understanding properly ?
>
> Thanks.
>
> --
> Dario Strbenac
> Research Assistant
> Cancer Epigenetics
> Garvan Institute of Medical Research
> Darlinghurst NSW 2010
> Australia
>
> __
> R-devel@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
>
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] DESCRIPTION file and Rd examples
On 04/15/2011 11:18 AM, Simon Urbanek wrote:
On Apr 14, 2011, at 11:00 PM, Dario Strbenac wrote:
I have a confusing error from R CMD check that I don't get when running the
example manually by hand.
In the \examples section of an Rd file, I create a GRanges object, then I call
a function with the GRanges object, whose first 2 lines are
require(GenomicRanges)
require() is doesn't guarantee that the package will load, so I think what you
meant to write was more
if (require(GenomicRanges, quietly=TRUE)) {
...
annoDF<- as.data.frame(anno) # anno is the GRanges object.
and that second line gives:
Error in as.data.frame.default(anno) :
cannot coerce class 'structure("GRanges", package = "GenomicRanges")' into a
data.frame
Calls: annoGR2DF ... annoGR2DF -> .local -> as.data.frame ->
as.data.frame.default
Try IRanges::as.data.frame(anno)
I'm guessing that your call finds base::as.data.frame, perhaps because
some earlier example has already require'd GenomicRanges, and that it's
definition of as.data.frame has been masked some time in between.
It's too complicated to debug in detail; R CMD check produces a file
.Rcheck/-Ex.R that contains the compiled example code, and it
is in the evaluation of this file that the error occurs. So you could
dissect it to discover the gory details.
Martin
I have GenomicRanges listed in my Imports: field, and IRanges in the Suggests:
of the DESCRIPTION file (it's require()d elsewhere). I'm trying to avoid
putting packages in Depends: , so my package loads fast. Any tips of what I'm
not understanding properly ?
Thanks.
--
Dario Strbenac
Research Assistant
Cancer Epigenetics
Garvan Institute of Medical Research
Darlinghurst NSW 2010
Australia
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
--
Computational Biology
Fred Hutchinson Cancer Research Center
1100 Fairview Ave. N. PO Box 19024 Seattle, WA 98109
Location: M1-B861
Telephone: 206 667-2793
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
