如何处理和变换变量？（三）

幂变换也可以稳定分布。下面的代码模拟了四个不同组的服从泊松分布的响应，其平均值分别为 5、10、15 和 20。

set.seed(523)
x <- rep(c(5,10,15,20), each=500)
y <- rpois(2000, lambda=x)

对于泊松分布来说，方差等于均值，因此均值20的组的方差要比均值5的组的方差大得多。

比较原始泊松随机变量与平方根变换变量的箱形图。

par(mfrow=c(1,2))
boxplot(y ~ x, main="Spread increases with level", xlab="level",
        ylab="Y ~ Poisson(level)")
boxplot(sqrt(y) ~ x, main="Spread is stabilised", xlab="level",
        ylab=expression(sqrt(Y))) 

 

从第一组方框图中观察各组响应的分布是如何扩大的。第二组方框图显示的是经平方根转换后的响应分布。各组内的分布明显更加均匀。换句话说，我们可以观察到，经过变换后的分布确实更集中了。

我们可以通过尝试不同的幂变换和检查图来找到一个能稳定分布的幂变换。Car软件包中的函数 spreadLevelPlot() 提供了一种更系统的方法。 我们添加了一个 0.5 的起始值，这样零计数就不会在对数变换中产生问题。函数 spreadLevelPlot() 可以生成表格和曲线图。在这里我们介绍一些常用的术语，“lower hinge”是各组样本的下四分位数。“upper hinge”是上四分位数。它们定义了Box图中方框的界限。“hinge spread”是四分位数之间的范围（IQR，四分位数范围 (IQR) 是第三个四分位数 (Q3) 和第一个四分位数 (Q1) 之间的差：IQR=Q3−Q1）。

spreadLevelPlot() 函数使用两个坐标轴的对数变换，绘制出hinge spread与中位数的对比图，然后对所绘制的点进行直线拟合。如果该直线的斜率为 b，则建议的幂变换（打印在表格底部）为 p = 1 - b。 

library(car)
par(mfrow=c(1,1))
spreadLevelPlot(y~x, start=0.5)

 

 在这种情况下，建议的幂变换是0.1125654。

p=0.1125654
boxplot(y^p ~ x, main="Suggested", xlab="level", ylab=expression(Y^p))

 

建议的幂变换和平方根变换相比，展现的是数据的一种极度的压缩，这个情况下的稳定程度已经非常强了。

（未完待续）