如何在 R 中使用 Box-Cox 幂变换

Posted 2023-02-16

【中文标题】如何在 R 中使用 Box-Cox 幂变换

【英文标题】：how to use the Box-Cox power transformation in R

【发布时间】：2016-03-04 03:25:19

【问题描述】：

我需要将一些数据转换为“正常形状”，并且我了解到 Box-Cox 可以识别用于转换数据的指数。

根据我的理解

car::boxCoxVariable(y)

用于线性模型中的响应变量，并且

MASS::boxcox(object)

用于公式或拟合模型对象。所以，因为我的数据是数据框的变量，所以我发现我可以使用的唯一函数是：

car::powerTransform(dataframe$variable, family="bcPower")

正确吗？还是我错过了什么？

第二个问题是我拿到后该怎么办

Estimated transformation parameters
dataframe$variable
0.6394806

我应该简单地将变量乘以这个值吗？ 我这样做了：

aaa = 0.6394806
dataframe$variable2 = (dataframe$variable)*aaa

然后我运行了 shapiro-wilks 检验的正态性，但我的数据似乎也不服从正态分布：

shapiro.test(dataframe$variable2)
data:  dataframe$variable2
W = 0.97508, p-value < 2.2e-16

【问题讨论】：

我发现这是一个很好的文档，其中包含清晰的 R 代码和示例（以及其他转换）手册中的章节 Transforming Data Summary and Analysis of Extension Program Evaluation in R

@Valentin 书中提到的很好的解释。非常感谢！

【参考方案1】：

Box 和 Cox (1964) 提出了一系列旨在减少线性模型中错误的非正态性的变换。事实证明，这样做通常也会减少非线性。

这里是对原始工作和之后所有工作的一个很好的总结：http://www.ime.usp.br/~abe/lista/pdfm9cJKUmFZp.pdf

但是，您会注意到，控制 lambda 幂变换选择的对数似然函数取决于基础模型的残差平方和（SO 上没有 LaTeX - 请参阅参考资料），因此没有变换无需模型即可应用。

一个典型的应用如下：

library(MASS)

# generate some data
set.seed(1)
n <- 100
x <- runif(n, 1, 5)
y <- x^3 + rnorm(n)

# run a linear model
m <- lm(y ~ x)

# run the box-cox transformation
bc <- boxcox(y ~ x)

(lambda <- bc$x[which.max(bc$y)])
[1] 0.4242424

powerTransform <- function(y, lambda1, lambda2 = NULL, method = "boxcox") 

  boxcoxTrans <- function(x, lam1, lam2 = NULL) 

    # if we set lambda2 to zero, it becomes the one parameter transformation
    lam2 <- ifelse(is.null(lam2), 0, lam2)

    if (lam1 == 0L) 
      log(y + lam2)
     else 
      (((y + lam2)^lam1) - 1) / lam1
    
  

  switch(method
         , boxcox = boxcoxTrans(y, lambda1, lambda2)
         , tukey = y^lambda1
  )



# re-run with transformation
mnew <- lm(powerTransform(y, lambda) ~ x)

# QQ-plot
op <- par(pty = "s", mfrow = c(1, 2))
qqnorm(m$residuals); qqline(m$residuals)
qqnorm(mnew$residuals); qqline(mnew$residuals)
par(op)

正如您所看到的，这不是灵丹妙药——只有一些数据可以被有效地转换（通常小于 -2 或大于 2 的 lambda 表示您不应该使用该方法）。与任何统计方法一样，在实施前请谨慎使用。

要使用两参数 Box-Cox 转换，请使用 geoR 包查找 lambda：

library("geoR")
bc2 <- boxcoxfit(x, y, lambda2 = TRUE)

lambda1 <- bc2$lambda[1]
lambda2 <- bc2$lambda[2]

编辑：@Yui-Shiuan 指出的 Tukey 和 Box-Cox 实现的合并已修复。

【讨论】：

您可能会指出，通过使用模型lm(y ~ 1)（在本例中为bc <- boxcox(variable ~ 1, data=dataframe)），使用MASS 可以获得相同的答案。 powerTransform() 给出了“正确”的 lambda，但数据中的某些情况使得仅使用 Box-Cox 无法强制正常。

很好的答案！请问您为什么在答案的开头强调“错误”？这是因为要转换响应变量（例如 y in y ~ x_1 + x_2）而不是协变量（x_1 或 x_2），还是可以同时转换协变量？

【参考方案2】：

根据论文中的Box-cox变换公式Box,George E. P.；考克斯，D.R. (1964)。 “转换分析”，我认为 mlegge 的帖子可能需要稍微编辑一下。转换后的 y 应该是 (y^(lambda)-1)/lambda 而不是 y^(lambda)。 （其实y^(lambda)被称为Tukey变换，是另一种不同的变换公式。） 所以，代码应该是：

(trans <- bc$x[which.max(bc$y)])
[1] 0.4242424
# re-run with transformation
mnew <- lm(((y^trans-1)/trans) ~ x) # Instead of mnew <- lm(y^trans ~ x) 

更多信息

在 R 中通过 boxcox() 正确实现 Box-Cox 变换公式：https://www.r-bloggers.com/on-box-cox-transform-in-regression-models/

Box-Cox 变换和 Tukey 变换之间的一个很好的比较。 http://onlinestatbook.com/2/transformations/box-cox.html

还可以在 Wikipedia 上找到 Box-Cox 变换公式： en.wikipedia.org/wiki/Power_transform#Box.E2.80.93Cox_transformation

如果我误解了，请纠正我。

【讨论】：

感谢您指出这一点（提供出色的文档！）。我已经更新了我的答案以尝试解决这个问题

【参考方案3】：

目前可以使用包geoR 将BoxCox 转换应用于数据，而无需任何底层模型。具体来说，您可以使用函数 boxcoxfit() 来寻找最佳参数，然后使用函数 BCtransform() 预测转换后的变量。

【参考方案4】：

如果我只想传输响应变量 y 而不是指定 x 的线性模型，例如我想传输/标准化数据列表，我可以为 x 取 1，然后对象变成线性模型：

library(MASS)
y = rf(500,30,30)
hist(y,breaks = 12)
result = boxcox(y~1, lambda = seq(-5,5,0.5))
mylambda = result$x[which.max(result$y)]
mylambda
y2 = (y^mylambda-1)/mylambda
hist(y2)

【讨论】：

有 4 个不同的变量，所有 4 个直方图都显示非正态分布，这个解决方案帮助我将它们分别带到正态分布。