[R] Box-cox transformation

Ravi Varadhan ravi.varadhan at jhu.edu
Mon Jul 7 22:50:21 CEST 2014


Thank you.  It is very helpful.

Ravi

-----Original Message-----
From: Joshua Wiley [mailto:jwiley.psych at gmail.com] 
Sent: Monday, July 07, 2014 4:15 PM
To: Ravi Varadhan
Cc: r-help at r-project.org
Subject: Re: [R] Box-cox transformation

Dear Ravi,

In my previous example, I used the residuals, so:

sum [ (r_i / scaling)^2 ]

If you want to use the deviance from glm, that gives you:

sum [ r_i^2 ]

and since the scaling factor is just a constant for any given lambda, then the modification would be:

sum [ r_i^2 ] / ( scaling^2 )

and is given in the modified code below (posted back to R-help in case any else has this question).

Hope this helps,

Josh


##########################################
require(MASS)

myp <- function(y, lambda) (y^lambda-1)/lambda


lambda <- seq(-0.05, 0.45, len = 20)
N <- nrow(quine)
res <- matrix(numeric(0), nrow = length(lambda), 2, dimnames = list(NULL, c("Lambda", "LL")))

# scaling contant
C <- exp(mean(log(quine$Days+1)))

for(i in seq_along(lambda)) {
  SS <- deviance(glm(myp(Days + 1, lambda[i]) ~ Eth*Sex*Age*Lrn, data = quine))
  LL <- (- (N/2) * log(SS/((C^lambda[i])^2)))
  res[i, ] <- c(lambda[i], LL)
}

# box cox
boxcox(Days+1 ~ Eth*Sex*Age*Lrn, data = quine, lambda = lambda) # add our points on top to verify match points(res[, 1], res[,2], pch = 16)

##########################################

On Mon, Jul 7, 2014 at 11:57 PM, Ravi Varadhan <ravi.varadhan at jhu.edu> wrote:
> Dear Josh,
> Thank you very much.  I knew that the scaling had to be adjusted, but was not sure on how to do this.
>
> Can you please show me how to do this scaling with `glm'? In other words, how would I scale the deviance from glm?
>
> Thanks,
> Ravi
>
> -----Original Message-----
> From: Joshua Wiley [mailto:jwiley.psych at gmail.com]
> Sent: Sunday, July 06, 2014 11:34 PM
> To: Ravi Varadhan
> Cc: r-help at r-project.org
> Subject: Re: [R] Box-cox transformation
>
> Hi Ravi,
>
> Deviance is the SS in this case, but you need a normalizing constant adjusted by the lambda to put them on the same scale.  I modified your example below to simplify slightly and use the normalization (see the LL line).
>
> Cheers,
>
> Josh
>
> ######################################
>
> require(MASS)
>
> myp <- function(y, lambda) (y^lambda-1)/lambda
>
>
> lambda <- seq(-0.05, 0.45, len = 20)
> N <- nrow(quine)
> res <- matrix(numeric(0), nrow = length(lambda), 2, dimnames = 
> list(NULL, c("Lambda", "LL")))
>
> # scaling contant
> C <- exp(mean(log(quine$Days+1)))
>
> for(i in seq_along(lambda)) {
>   r <- resid(lm(myp(Days + 1, lambda[i]) ~ Eth*Sex*Age*Lrn, data = quine))
>   LL <- (- (N/2) * log(sum((r/(C^lambda[i]))^2)))
>   res[i, ] <- c(lambda[i], LL)
> }
>
> # box cox
> boxcox(Days+1 ~ Eth*Sex*Age*Lrn, data = quine, lambda = lambda) # add 
> our points on top to verify match points(res[, 1], res[,2], pch = 16)
>
> ######################################
>
>
>
> On Mon, Jul 7, 2014 at 11:33 AM, Ravi Varadhan <ravi.varadhan at jhu.edu> wrote:
>> Hi,
>>
>> I am trying to do Box-Cox transformation, but I am not sure how to do it correctly.  Here is an example showing what I am trying:
>>
>>
>>
>> # example from MASS
>>
>> require(MASS)
>> boxcox(Days+1 ~ Eth*Sex*Age*Lrn, data = quine,
>>        lambda = seq(-0.05, 0.45, len = 20))
>>
>> # Here is My attempt at getting the profile likelihood for the 
>> Box-Cox parameter lam <- seq(-0.05, 0.45, len = 20) dev <- rep(NA, 
>> length=20)
>>
>> for (i in 1:20) {
>> a <- lam[i]
>> ans <- glm(((Days+1)^a-1)/a ~ Eth*Sex*Age*Lrn, family=gaussian, data 
>> =
>> quine) dev[i] <- ans$deviance }
>>
>> plot(lam, dev, type="b", xlab="lambda", ylab="deviance")
>>
>> I am trying to create the profile likelihood for the Box-Cox parameter, but obviously I am not getting it right.  I am not sure that ans$deviance is the right thing to do.
>>
>> I would appreciate any guidance.
>>
>> Thanks & Best,
>> Ravi
>>
>>
>>
>>
>>         [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help at 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.
>
>
>
> --
> Joshua F. Wiley
> Ph.D. Student, UCLA Department of Psychology http://joshuawiley.com/ Senior Analyst, Elkhart Group Ltd.
> http://elkhartgroup.com
> Office: 260.673.5518



--
Joshua F. Wiley
Ph.D. Student, UCLA Department of Psychology http://joshuawiley.com/ Senior Analyst, Elkhart Group Ltd.
http://elkhartgroup.com
Office: 260.673.5518


More information about the R-help mailing list