使用 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?
就我尝试使用各种数据(包括其他数据)而言,实际差异似乎仅限于比例偏移和镜像转换。
- 基于旋转矩阵的方法的结果显示在左侧面板中。
- 基于 pcomp 函数的方法的结果显示在右侧面板中。
“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")
进行了两种不同的主成分分析方法,以使用下面的 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?
就我尝试使用各种数据(包括其他数据)而言,实际差异似乎仅限于比例偏移和镜像转换。
- 基于旋转矩阵的方法的结果显示在左侧面板中。
- 基于 pcomp 函数的方法的结果显示在右侧面板中。
“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")