[R] Getting an error calling MASS::boxcox in a function
Bert Gunter
bgunter@4567 @end|ng |rom gm@||@com
Sat Jul 8 23:08:07 CEST 2023
Aha. Many thanks, John. Never would have gotten there on my own.
-- Bert
On Sat, Jul 8, 2023 at 2:01 PM John Fox <jfox using mcmaster.ca> wrote:
>
> Hi Bert,
>
> On 2023-07-08 3:42 p.m., Bert Gunter wrote:
> > Caution: This email may have originated from outside the organization. Please exercise additional caution with any links and attachments.
> >
> >
> > Thanks John.
> >
> > ?boxcox says:
> >
> > *************************
> > Arguments
> >
> > object
> >
> > a formula or fitted model object. Currently only lm and aov objects are handled.
> > *************************
> > I read that as saying that
> >
> > boxcox(lm(z+1 ~ 1),...)
> >
> > should run without error. But it didn't. And perhaps here's why:
> > BoxCoxLambda <- function(z){
> > b <- MASS:::boxcox.lm(lm(z+1 ~ 1), lambda = seq(-5, 5, length.out =
> > 61), plotit = FALSE)
> > b$x[which.max(b$y)] # best lambda
> > }
> >
> >> lambdas <- apply(dd,2 , BoxCoxLambda)
> > Error in NextMethod() : 'NextMethod' called from an anonymous function
> >
> > and, indeed, ?UseMethod says:
> > "NextMethod should not be called except in methods called by UseMethod
> > or from internal generics (see InternalGenerics). In particular it
> > will not work inside anonymous calling functions (e.g.,
> > get("print.ts")(AirPassengers))."
> >
> > ****BUT ****
> > BoxCoxLambda <- function(z){
> > b <- MASS:::boxcox(z+1 ~ 1, lambda = seq(-5, 5, length.out = 61),
> > plotit = FALSE)
> > b$x[which.max(b$y)] # best lambda
> > }
> >
> >> lambdas <- apply(dd,2 , BoxCoxLambda)
> >> lambdas
> > [1] 0.1666667 0.1666667
>
> As it turns out, it's the update() step in boxcox.lm() that fails, and
> the update takes place because $y is missing from the lm object, so the
> following works:
>
> BoxCoxLambda <- function(z){
> b <- boxcox(lm(z + 1 ~ 1, y=TRUE),
> lambda = seq(-5, 5, length.out = 101),
> plotit = FALSE)
> b$x[which.max(b$y)]
> }
>
> >
> > The identical lambdas do not seem right to me;
>
> I think that's just an accident of the example (using the BoxCoxLambda()
> above):
>
> > apply(dd, 2, BoxCoxLambda, simplify = TRUE)
> [1] 0.2 0.2
>
> > dd[, 2] <- dd[, 2]^3
> > apply(dd, 2, BoxCoxLambda, simplify = TRUE)
> [1] 0.2 0.1
>
> Best,
> John
>
> > nor do I understand why
> > boxcox.lm apparently throws the error while boxcox.formula does not
> > (it also calls NextMethod()) So I would welcome clarification to clear
> > my clogged (cerebral) sinuses. :-)
> >
> >
> > Best,
> > Bert
> >
> >
> > On Sat, Jul 8, 2023 at 11:25 AM John Fox <jfox using mcmaster.ca> wrote:
> >>
> >> Dear Ron and Bert,
> >>
> >> First (and without considering why one would want to do this, e.g.,
> >> adding a start of 1 to the data), the following works for me:
> >>
> >> ------ snip ------
> >>
> >> > library(MASS)
> >>
> >> > BoxCoxLambda <- function(z){
> >> + b <- boxcox(z + 1 ~ 1,
> >> + lambda = seq(-5, 5, length.out = 101),
> >> + plotit = FALSE)
> >> + b$x[which.max(b$y)]
> >> + }
> >>
> >> > mrow <- 500
> >> > mcol <- 2
> >> > set.seed(12345)
> >> > dd <- matrix(rgamma(mrow*mcol, shape = 2, scale = 5), nrow = mrow, ncol =
> >> + mcol)
> >>
> >> > dd1 <- dd[, 1] # 1st column of dd
> >> > res <- boxcox(lm(dd1 + 1 ~ 1), lambda = seq(-5, 5, length.out = 101),
> >> plotit
> >> + = FALSE)
> >> > res$x[which.max(res$y)]
> >> [1] 0.2
> >>
> >> > apply(dd, 2, BoxCoxLambda, simplify = TRUE)
> >> [1] 0.2 0.2
> >>
> >> ------ snip ------
> >>
> >> One could also use the powerTransform() function in the car package,
> >> which in this context transforms towards *multi*normality:
> >>
> >> ------ snip ------
> >>
> >> > library(car)
> >> Loading required package: carData
> >>
> >> > powerTransform(dd + 1)
> >> Estimated transformation parameters
> >> Y1 Y2
> >> 0.1740200 0.2089925
> >>
> >> I hope this helps,
> >> John
> >>
> >> --
> >> John Fox, Professor Emeritus
> >> McMaster University
> >> Hamilton, Ontario, Canada
> >> web: https://www.john-fox.ca/
> >>
> >> On 2023-07-08 12:47 p.m., Bert Gunter wrote:
> >>> Caution: This email may have originated from outside the organization. Please exercise additional caution with any links and attachments.
> >>>
> >>>
> >>> No, I'm afraid I'm wrong. Something went wrong with my R session and gave
> >>> me incorrect answers. After restarting, I continued to get the same error
> >>> as you did with my supposed "fix." So just ignore what I said and sorry for
> >>> the noise.
> >>>
> >>> -- Bert
> >>>
> >>> On Sat, Jul 8, 2023 at 8:28 AM Bert Gunter <bgunter.4567 using gmail.com> wrote:
> >>>
> >>>> Try this for your function:
> >>>>
> >>>> BoxCoxLambda <- function(z){
> >>>> y <- z
> >>>> b <- boxcox(y + 1 ~ 1,lambda = seq(-5, 5, length.out = 61), plotit =
> >>>> FALSE)
> >>>> b$x[which.max(b$y)] # best lambda
> >>>> }
> >>>>
> >>>> ***I think*** (corrections and clarification strongly welcomed!) that `~`
> >>>> (the formula function) is looking for 'z' in the GlobalEnv, the caller of
> >>>> apply(), and not finding it. It finds 'y' here explicitly in the
> >>>> BoxCoxLambda environment.
> >>>>
> >>>> Cheers,
> >>>> Bert
> >>>>
> >>>>
> >>>>
> >>>> On Sat, Jul 8, 2023 at 4:28 AM Ron Crump via R-help <r-help using r-project.org>
> >>>> wrote:
> >>>>
> >>>>> Hi,
> >>>>>
> >>>>> Firstly, apologies as I have posted this on community.rstudio.com too.
> >>>>>
> >>>>> I want to optimise a Box-Cox transformation on columns of a matrix (ie, a
> >>>>> unique lambda for each column). So I wrote a function that includes the
> >>>>> call to MASS::boxcox in order that it can be applied to each column easily.
> >>>>> Except that I'm getting an error when calling the function. If I just
> >>>>> extract a column of the matrix and run the code not in the function, it
> >>>>> works. If I call the function either with an extracted column (ie dd1 in
> >>>>> the reprex below) or in a call to apply I get an error (see the reprex
> >>>>> below).
> >>>>>
> >>>>> I'm sure I'm doing something silly, but I can't see what it is. Any help
> >>>>> appreciated.
> >>>>>
> >>>>> library(MASS)
> >>>>>
> >>>>> # Find optimised Lambda for Boc-Cox transformation
> >>>>> BoxCoxLambda <- function(z){
> >>>>> b <- boxcox(lm(z+1 ~ 1), lambda = seq(-5, 5, length.out = 61), plotit
> >>>>> = FALSE)
> >>>>> b$x[which.max(b$y)] # best lambda
> >>>>> }
> >>>>>
> >>>>> mrow <- 500
> >>>>> mcol <- 2
> >>>>> set.seed(12345)
> >>>>> dd <- matrix(rgamma(mrow*mcol, shape = 2, scale = 5), nrow = mrow, ncol =
> >>>>> mcol)
> >>>>>
> >>>>> # Try it not using the BoxCoxLambda function:
> >>>>> dd1 <- dd[,1] # 1st column of dd
> >>>>> bb <- boxcox(lm(dd1+1 ~ 1), lambda = seq(-5, 5, length.out = 101), plotit
> >>>>> = FALSE)
> >>>>> print(paste0("1st column's lambda is ", bb$x[which.max(bb$y)]))
> >>>>> #> [1] "1st column's lambda is 0.2"
> >>>>>
> >>>>> # Calculate lambda for each column of dd
> >>>>> lambdas <- apply(dd, 2, BoxCoxLambda, simplify = TRUE)
> >>>>> #> Error in eval(predvars, data, env): object 'z' not found
> >>>>>
> >>>>> Created on 2023-07-08 with reprex v2.0.2
> >>>>>
> >>>>> Thanks for your time and help.
> >>>>>
> >>>>> Ron
> >>>>> [[alternative HTML version deleted]]
> >>>>>
> >>>>> ______________________________________________
> >>>>> R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
> >>>>> 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.
> >>>>>
> >>>>
> >>>
> >>> [[alternative HTML version deleted]]
> >>>
> >>> ______________________________________________
> >>> R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
> >>> 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.
> >>
More information about the R-help
mailing list