重采样未产生主成分分析的预期结果

Question

我正在尝试使用以下代码来使用带替换的重采样（如 bootstrap）生成主成分分析的置信区间。我正在使用 iris 数据集的前 4 列：

prcomp 函数产生以下输出：

> mydf = iris[1:4]
> print(prcomp(mydf))
Standard deviations:
[1] 2.0562689 0.4926162 0.2796596 0.1543862

Rotation:
                     PC1         PC2         PC3        PC4
Sepal.Length  0.36138659 -0.65658877  0.58202985  0.3154872
Sepal.Width  -0.08452251 -0.73016143 -0.59791083 -0.3197231
Petal.Length  0.85667061  0.17337266 -0.07623608 -0.4798390
Petal.Width   0.35828920  0.07548102 -0.54583143  0.7536574

使用带替换的重采样：

> times = 1000
> ll = list()
> for(i in 1:times) {
+ tempdf =  mydf[sample(nrow(mydf), replace = TRUE), ]
+ ll[[length(ll)+1]] = prcomp(tempdf)$rotation
+ }
> 
> dd = data.frame(apply(simplify2array(ll), 1:2, mean))
> print(dd)
                      PC1          PC2          PC3          PC4
Sepal.Length  0.005574165 -0.039480258  0.044537991  0.007778055
Sepal.Width  -0.002587333 -0.040273812 -0.050793200 -0.005473271
Petal.Length  0.015681233  0.010952361 -0.005769051 -0.011351172
Petal.Width   0.006513656  0.008296928 -0.041805210  0.019109323

确定较低的置信区间：

> ddlower = data.frame(apply(simplify2array(ll), 1:2, quantile, probs=0.025))
> print(ddlower)
                    PC1        PC2        PC3        PC4
Sepal.Length -0.3859257 -0.7274809 -0.6560139 -0.3807826
Sepal.Width  -0.1127749 -0.7907801 -0.6818251 -0.3941001
Petal.Length -0.8633386 -0.2058064 -0.1333520 -0.4919584
Petal.Width  -0.3702979 -0.1328146 -0.6203322 -0.8088710

确定置信区间上限：

> ddupper = data.frame(apply(simplify2array(ll), 1:2, quantile, probs=0.975))
> print(ddupper)
                   PC1       PC2       PC3       PC4
Sepal.Length 0.3860431 0.7250412 0.6632126 0.3831889
Sepal.Width  0.1111863 0.7993649 0.6758156 0.3987939
Petal.Length 0.8638549 0.2106540 0.1318556 0.4915670
Petal.Width  0.3721362 0.1510708 0.6246988 0.8083421

我发现加载值非常不同。此外，所有变量和成分的置信区间都在 0 左右。我还检查了其他（大型）数据集，结果非常相似。从这些置信区间 none 的载荷与 0 明显不同。代码中显然存在一些错误，但我似乎找不到它。感谢您的帮助。

Answer 1

鉴于未定义特征向量的符号（您可以翻转配置并获得相同的结果），在 signed[= 上形成置信区间没有意义44=] 加载值。

而是计算加载的绝对值的置信区间，而不是有符号值。

想一想当 Sepal.Length 的特征向量从 ~ -0.3 翻转到 ~ +0.3 时，你的区间会发生什么变化？从绝对尺寸的角度考虑时，两种情况下的负载相似。但是，当您查看实际的带符号值时，加载平均为 0 是合乎逻辑的，因为您平均了很多 ~-0.3s 和 ~0.3s。

为了形象化您最初尝试失败的原因，运行：

set.seed(1) mydf <- iris[1:4] times <- 1000 ll <- vector(mode = "list", length = times) for (i in seq_len(times)) { tempdf <- mydf[sample(nrow(mydf), replace = TRUE), ] ll[[i]] <- prcomp(tempdf)$rotation }

这实际上是您的代码，已根据我的感受进行了修改。现在在 PC1 上提取 Sepal.Length 的负载并绘制值的直方图：

hist(sapply(ll, `[`, 1, 1))

产生

而是计算加载的 绝对值 的置信区间，而不是有符号值。

例如

set.seed(1) mydf <- iris[1:4] times <- 1000 ll <- vector(mode = "list", length = times) for (i in seq_len(times)) { tempdf <- mydf[sample(nrow(mydf), replace = TRUE), ] ll[[i]] <- abs(prcomp(tempdf)$rotation) ## NOTE: abs(...) }

这给出：

> data.frame(apply(simplify2array(ll), 1:2, quantile, probs = 0.025)) PC1 PC2 PC3 PC4 Sepal.Length 0.33066830 0.578558222 0.45955051 0.2252653 Sepal.Width 0.05211013 0.623424084 0.49591685 0.2351746 Petal.Length 0.84823899 0.133137927 0.01226608 0.4607265 Petal.Width 0.34284824 0.007403214 0.44932031 0.6780493 > data.frame(apply(simplify2array(ll), 1:2, quantile, probs = 0.975)) PC1 PC2 PC3 PC4 Sepal.Length 0.3891499 0.7443276 0.6690553 0.3898237 Sepal.Width 0.1186205 0.7988607 0.7010495 0.4083784 Petal.Length 0.8653324 0.2153410 0.1450756 0.4933340 Petal.Width 0.3742441 0.1645692 0.6350899 0.8154254

重采样未产生主成分分析的预期结果

Resampling not producing expected result of principal component analysis

r

statistics-bootstrap