FactoMineR::PCA() 中的 PCA 不返回与 base::prcomp() 相同的结果

MLE*_*LEN 4 r pca

进行一些 PCA 分析,并与FactoMineR函数的结果进行比较时,我没有得到相同的结果。一个例子PCAprcompbase

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)
Run Code Online (Sandbox Code Playgroud)

doc*_*ate 6

抱歉迟到了,该FactoMineR软件包使用相同的方法,svd()该方法应该与该prcomp()方法相似(但不相同),并且两者都列在 Q 模式下,Q 模式是进行 PCA 的首选方法,因为其数值精度。但请注意,我没有说相同,为什么?FactoMineR使用自己的 PCA 算法,计算组件数量,如下所示:

\n
ncp <- min(ncp, nrow(X) - 1, ncol(X))\n
Run Code Online (Sandbox Code Playgroud)\n

它清楚地告诉你为什么你得到的组件数量是 63 而不是prcomp()通常给出的 64。您的数据集是典型的基因组学数据,其中n行小于p基因列,上面的代码显然将采用列或行,以数量较少者为准。如果您遵循该svd()算法,它将返回 64 维而不是 63。

\n

要进一步探索源代码,请键入FactoMineR:::PCA.

\n

对于 Q 模式 ( svd,, ) 和 R 模式 ( , ) 之间的差异prcomp(),我建议访问此答案FactoMineR::PCA()eigen()princomp()

\n

旁注:因为prcomp()您希望在进行 PCA 之前传递center = T参数以便将数据居中。另一方面,缩放将使所有基因列的权重相等。

\n
pca_prcomp <- prcomp(df, center = T, scale. = T) # add center=T\n
Run Code Online (Sandbox Code Playgroud)\n

对于缩放,prcomp()使用while用作N除数。下面的代码将证明这一点(请参阅上面相同的链接答案):FactoMineR::PCA()N-1

\n
# this is the scaled data by scale()\ndf_scaled <- scale(df)\n\n# then you need to get the standardized data matrix from the output of the FactoMineR::PCR() function, which can be done easily as follows:\ndf_restored <- pca_facto$svd$U %*% diag(pca_facto$svd$vs) %*% t(pca_facto$svd$V)\n\n# the to make both FactoMineR::PCR() and scale() match up you need to do the correction\ndf_corrected <- df_restored * sqrt(63 / 64) # correct for sqrt(N-1/N)\n\nhead(df[, 1:5]) # glimpse the first five columns only!\nhead(df_scaled[, 1:5])\nhead(df_restored[, 1:5]) # glimpse the first five columns only!\nhead(df_corrected[, 1:5])\nround(head(df_scaled[, 1:5]), 3) == round(head(df_corrected[, 1:5]), 3) # TRUE\n\nR> head(df[, 1:5])\n       1      2      3      4      5\nV1 0.300  1.180  0.550  1.140 -0.265\nV2 0.680  1.290  0.170  0.380  0.465\nV3 0.940 -0.040 -0.170 -0.040 -0.605\nV4 0.280 -0.310  0.680 -0.810  0.625\nV5 0.485 -0.465  0.395  0.905  0.200\nV6 0.310 -0.030 -0.100 -0.460 -0.205\nR> head(df_scaled[, 1:5])\n       1        2      3      4      5\nV1 0.723  1.59461  1.315  1.345 -0.600\nV2 1.584  1.73979  0.438  0.649  0.905\nV3 2.173 -0.01609 -0.346  0.264 -1.301\nV4 0.678 -0.37256  1.615 -0.441  1.235\nV5 1.142 -0.57720  0.958  1.130  0.359\nV6 0.746 -0.00289 -0.185 -0.120 -0.476\nR> head(df_restored[, 1:5])\n      [,1]     [,2]   [,3]   [,4]   [,5]\n[1,] 0.729  1.60722  1.326  1.356 -0.605\n[2,] 1.596  1.75354  0.442  0.654  0.912\n[3,] 2.190 -0.01622 -0.349  0.266 -1.311\n[4,] 0.683 -0.37550  1.628 -0.444  1.244\n[5,] 1.151 -0.58176  0.965  1.139  0.361\n[6,] 0.752 -0.00291 -0.186 -0.121 -0.480\nR> head(df_corrected[, 1:5])\n      [,1]     [,2]   [,3]   [,4]   [,5]\n[1,] 0.723  1.59461  1.315  1.345 -0.600\n[2,] 1.584  1.73979  0.438  0.649  0.905\n[3,] 2.173 -0.01609 -0.346  0.264 -1.301\n[4,] 0.678 -0.37256  1.615 -0.441  1.235\n[5,] 1.142 -0.57720  0.958  1.130  0.359\n[6,] 0.746 -0.00289 -0.185 -0.120 -0.476\nR> round(head(df_scaled[, 1:5]), 3) == round(head(df_corrected[, 1:5]), 3)\n      1    2    3    4    5\nV1 TRUE TRUE TRUE TRUE TRUE\nV2 TRUE TRUE TRUE TRUE TRUE\nV3 TRUE TRUE TRUE TRUE TRUE\nV4 TRUE TRUE TRUE TRUE TRUE\nV5 TRUE TRUE TRUE TRUE TRUE\nV6 TRUE TRUE TRUE TRUE TRUE\n
Run Code Online (Sandbox Code Playgroud)\n

本书摘录

\n

还有FactoMineRFran\xc3\xa7ois Husson、S\xc3\xa9bastien 和 L\xc3\xaa J\xc3\xa9r\xc3\xb4me Pag\ 所著的名为“Exploratory Multivariate Analysis by Example using R”第二版的书。 xc3\xa8s。下面是本书第 55 页的摘录,其中讨论了与您的类似的基因组研究的数据集,其中n行 (43) 远少于p7407 列chicken.csv数据集,您可以在他们的网站以及数据集中查看更多信息本身可以从此链接下载。

\n

在此输入图像描述

\n