使用 R 进行主成分分析。自动和手动结果不匹配

Principal component analysis using R. Automatic and manual results do not match

进行了两种不同的主成分分析方法,以使用下面的 Box1 的 R 代码分析以下数据 (ch082.dat)。
https://drive.google .com/file/d/1xykl6ln-bUnXIs-jIA3n5S3XgHjQbkWB/view?usp=共享

第一种方法使用旋转矩阵(参见 Box1 代码的“#rotated data”下的 'ans_mat')并且, 第二种方法使用 'pcomp' 函数(参见 Box1 代码的“#rotated data”下的 'rpca')。

然而,使用旋转矩阵的方法和使用'pcomp'函数的方法在答案上存在细微的差异。 使其匹配

我的问题

What should I do so that the result of the rotation matrix -based method matches the result of the'pcomp' function?

就我尝试使用各种数据(包括其他数据)而言,实际差异似乎仅限于比例偏移和镜像转换。

“ch082.dat”数据可以看到镜像反转(见图1); 看起来,在某些j中,“相关矩阵的第j个特征向量”和“prcomp函数输出值的第j列”的符号可能是相反的。如果特征值有一定程度的重叠,则差异可能比镜像反演更复杂。
图1

尽管对数据进行了集中化和标准化,但 Box2 的数据存在比例偏移(见图 2)。
图2

Box.1

#dataload
##Use the 'setwd' function to specify the directory containing 'ch082.dat'.
##For example, if you put this file directly under the C drive of your Windows PC, you can run the following command.

setwd("C:/") #Depending on where you put the file, you may need to change the path.
getwd()
w1<-read.table("ch082.dat",header = TRUE,row.names = 1,fileEncoding = "UTF-8")
w1

#Function for standardizing data
#Thanks to https://qiita.com/ohisama2/items/5922fac0c8a6c21fcbf8
standalize <- function(data)
{ for(i in length(data[1,]))
{
  x <- as.matrix(data[,i])
  y <- (x-mean(x)/sd(x))
  data[,i] <- y
}
  return(data)}


#Method using rotation matrix
z_=standalize(w1)
B_mat=cor(z_) #Compute correlation matrix
eigen_m <- eigen(B_mat)
sample_mat <- as.matrix(z_)
ans_mat=sample_mat

for(j in 1:length(sample_mat[1,])){
ans_mat[,j]=sample_mat%*%eigen_m$vectors[,j]
}

#Method using "rpca" function
rpca <- prcomp(w1,center=TRUE, scale=TRUE)

#eigen vectors
eigen_m$vectors
rpca

#rotated data
ans_mat
rpca$x

#Graph Plots
par(mfrow=c(1,2))
plot(
  ans_mat[,1], 
  ans_mat[,2], 
  main="Rotation using eigenvectors"
)
plot(rpca$x[,1], rpca$x[,2], 
     main="Principal component score")
par(mfrow=c(1,1))

#summary
summary(rpca)$importance

Box2.

sample_data <- data.frame(
  X = c(2,4, 6, 5,7, 8,10),
  Y = c(6,8,10,11,9,12,14)
)

X = c(2,4, 6, 5,7, 8,10)
Y = c(6,8,10,11,9,12,14)
plot(Y ~ X)

w1=sample_data 

参考
https://logics-of-blue.com/principal-components-analysis/ (用日文写的)

两组结果一致。首先,我们可以稍微简化您的代码。您不需要函数或 for 循环:

z_ <- scale(w1)
B_mat <- cor(z_)
eigen_m <- eigen(B_mat)
ans_mat <- z_ %*% eigen_m$vectors

现在 prcomp 版本

z_pca <- prcomp(z_)
z_pca$sdev^2    # Equals eigen_m$values
z_pca$rotation  # Equals eigen_m$vectors
z_pca$x         # Equals ans_mat

您的原始代码错误标记了 ans_mat 列。它们实际上是主成分分数。你可以用

解决这个问题
colnames(ans_mat) <- colnames(z_pca$x)

PC 负载(以及分数)在反射方面不是唯一定义的。换句话说,将一个分量中的所有载荷或分数乘以 -1 会翻转它们,但不会改变它们之间的关系。将 z_pca$x[, 1] 乘以 -1,绘图将匹配:

z_pca$x[, 1] <- z_pca$x[, 1] * -1
dev.new(width=10, height=6)
par(mfrow=c(1,2))
plot(ans_mat[,1], ans_mat[,2], main="Rotation using eigenvectors")
plot(z_pca$x[,1], z_pca$x[,2], main="Principal component score")