做了一些主成分分析,并与FactoMineR函数、PCA和base的prcomp结果进行了比较,得到了不同的结果。一个例子
library(ISLR)
library(FactoMineR)
data("NCI60")
df <- NCI60$data
pca_prcomp <- prcomp(df, scale. = T)
pca_facto <- FactoMineR::PCA(df, scale.unit = T, graph = F, ncp = 65)
# One column is missing
dim(pca_prcomp$x)
dim(pca_facto$ind$coord)
# Values are similiare - but not the same
head(pca_prcomp$x[, 1:2])
head(pca_facto$ind$coord[, 1:2])
# Using scale function - does not return same values
pca_facto_scale <- PCA(scale(df), scale.unit = F, graph = F, ncp = 65)
head(pca_facto$ind$coord[, 1:2], 3)
head(pca_facto_scale$ind$coord[, 1:2], 3)发布于 2020-10-30 11:46:45
遗憾的是,FactoMineR包使用了与prcomp()方法类似(但不完全相同)的svd()方法,这两种方法都列在Q模式下,这是进行主成分分析的首选方法。但请注意,我没有说相同,为什么?FactoMineR使用自己的算法对PCA计算组件数,如下所示:
ncp <- min(ncp, nrow(X) - 1, ncol(X))这清楚地说明了为什么您有63个组件,而不是prcomp()通常提供的64个组件。您的数据集是典型的基因组数据,其中n行小于p列的基因,而上面的代码将明显地接受列或行,无论哪个列或行有较少的数目。如果您遵循svd()算法,它将返回64维,而不是63维。
要进一步探索源代码,请输入FactoMineR:::PCA。
对于Q模式(svd,prcomp(),FactoMineR::PCA())和R模式(eigen(),princomp())之间的差异,我建议访问这个回答。
侧请注意:对于prcomp(),您希望传递center = T参数,以便在进行prcomp()之前对数据进行中心化。另一方面,缩放会给你所有的基因列同等的重量。
pca_prcomp <- prcomp(df, center = T, scale. = T) # add center=T对于缩放,prcomp()使用N作为除数,而FactoMineR::PCA()使用N-1。下面的代码将证明这一点(参见上面相同的链接答案):
# this is the scaled data by scale()
df_scaled <- scale(df)
# then you need to get the standardized data matrix from the output of the FactoMineR::PCR() function, which can be done easily as follows:
df_restored <- pca_facto$svd$U %*% diag(pca_facto$svd$vs) %*% t(pca_facto$svd$V)
# the to make both FactoMineR::PCR() and scale() match up you need to do the correction
df_corrected <- df_restored * sqrt(63 / 64) # correct for sqrt(N-1/N)
head(df[, 1:5]) # glimpse the first five columns only!
head(df_scaled[, 1:5])
head(df_restored[, 1:5]) # glimpse the first five columns only!
head(df_corrected[, 1:5])
round(head(df_scaled[, 1:5]), 3) == round(head(df_corrected[, 1:5]), 3) # TRUE
R> head(df[, 1:5])
1 2 3 4 5
V1 0.300 1.180 0.550 1.140 -0.265
V2 0.680 1.290 0.170 0.380 0.465
V3 0.940 -0.040 -0.170 -0.040 -0.605
V4 0.280 -0.310 0.680 -0.810 0.625
V5 0.485 -0.465 0.395 0.905 0.200
V6 0.310 -0.030 -0.100 -0.460 -0.205
R> head(df_scaled[, 1:5])
1 2 3 4 5
V1 0.723 1.59461 1.315 1.345 -0.600
V2 1.584 1.73979 0.438 0.649 0.905
V3 2.173 -0.01609 -0.346 0.264 -1.301
V4 0.678 -0.37256 1.615 -0.441 1.235
V5 1.142 -0.57720 0.958 1.130 0.359
V6 0.746 -0.00289 -0.185 -0.120 -0.476
R> head(df_restored[, 1:5])
[,1] [,2] [,3] [,4] [,5]
[1,] 0.729 1.60722 1.326 1.356 -0.605
[2,] 1.596 1.75354 0.442 0.654 0.912
[3,] 2.190 -0.01622 -0.349 0.266 -1.311
[4,] 0.683 -0.37550 1.628 -0.444 1.244
[5,] 1.151 -0.58176 0.965 1.139 0.361
[6,] 0.752 -0.00291 -0.186 -0.121 -0.480
R> head(df_corrected[, 1:5])
[,1] [,2] [,3] [,4] [,5]
[1,] 0.723 1.59461 1.315 1.345 -0.600
[2,] 1.584 1.73979 0.438 0.649 0.905
[3,] 2.173 -0.01609 -0.346 0.264 -1.301
[4,] 0.678 -0.37256 1.615 -0.441 1.235
[5,] 1.142 -0.57720 0.958 1.130 0.359
[6,] 0.746 -0.00289 -0.185 -0.120 -0.476
R> round(head(df_scaled[, 1:5]), 3) == round(head(df_corrected[, 1:5]), 3)
1 2 3 4 5
V1 TRUE TRUE TRUE TRUE TRUE
V2 TRUE TRUE TRUE TRUE TRUE
V3 TRUE TRUE TRUE TRUE TRUE
V4 TRUE TRUE TRUE TRUE TRUE
V5 TRUE TRUE TRUE TRUE TRUE
V6 TRUE TRUE TRUE TRUE TRUE图书摘录
此外,还有弗朗索瓦·胡森(Fran Ois Husson )、塞巴斯蒂安( Sébastien )和莱热梅·帕热斯( LêJér me Pagès )的第二版“探索性多变量分析”(探索性多变量分析--译注)。以下是该书第55页的摘录,其中讨论了一项与你的基因组研究相似的数据集,其中n行(43)远远低于p 7407列chicken.csv数据集,你可以在他们的网站上看到更多的信息,以及数据集本身可以从这个链接下载。

发布于 2018-01-17 10:21:38
特征值与用于执行主成分分析的奇异值方法之间的差别很大(一些细节请参见这个伟大的答案 )。
来自?prcomp
计算是通过数据矩阵的奇异值分解来完成的,而不是在协方差矩阵上使用“特征”。这通常是数值精度的首选方法。
来自?常设仲裁院:
返回一个列表,包括: eig:包含所有特征值、方差百分比和累积方差百分比的矩阵。
https://stackoverflow.com/questions/48298385
复制相似问题