R 的 prcomp()

写在前面：与 PCA 有关的一些名词解释，参考：https://zhuanlan.zhihu.com/p/414190170 。

R 的 prcomp() 进行 PCA 降维时，是通过对 (中心化和可能标准化的) 数据矩阵进行奇异值分解计算的，而不是使用协方差矩阵的特征值（eigen）。
与 princomp() 不同，prcomp() 的方差是用除数 N-1 计算的。
prcomp() 的结果是一个 prcomp 类对象，是一个列表，包括以下组分：

• sdev：主成分的标准差（协方差矩阵/相关性矩阵的特征值的平方根，虽然计算是用数据矩阵的奇异值完成的）。
• rotation：variable loadings（变量载荷） 矩阵，列包含 eigenvectors（特征向量），对应to的 princomp 返回的 loadings。
• x：如果 retx=TRUE，返回 rotated data，即（中心化(和标准化)的）数据乘以 rotation 矩阵。因此，cov(x) 是对角矩阵 diag(sdev^2)。
• center，scale：使用的中心化和缩放，或者 FALSE。

注意：rotation matrix 的列的符号是任意的，可能在不同的 PCA 程序，甚至不同的 R builds 之间都有所不同。因此，最好将 PCA 模型保存下来，以便在不同的环境下使用。
实例：

# 
file <- './input.data.txt'

data_frame <- read.table(file, header = T, sep = "\t", row.names = 1, check.names=F)
data_tr <- data_frame[1:35, ]
data_te <- data_frame[36:60, ]

# 基于训练集进行 pca 降维
pca_model <- prcomp(data_tr, center=TRUE, scale.=TRUE)
# 使用训练集的参数，对测试集进行降维
pred_tr <- predict(pca_model, data_te)

• 训练集和测试集的维度

print(dim(data_tr))
print(dim(data_te))

image.png

• 主成分的标准差，可以用于计算我们常说的每个 PC 的贡献度（每个主成分能解释多少原始数据中的变异，即方差。参考）

print(pca_model$sdev)

pca_model$sdev

计算每个特征的贡献度：

print(pca_model$sdev^2)

pca_model$sdev^2

# sum(pca_model$sdev^2)=39，所有 PC 可解释所有变量
# pca_model$sdev^2/sum(pca_model$sdev^2)，每个 PC 可解释变量占比
ratio <- pca_model$sdev^2/sum(pca_model$sdev^2)
# 累积贡献率
cumsum(ratio)

• 变量载荷矩阵，每一列代表的是特征向量。

print(pca_model$rotation)

rotation 矩阵

• 降维后的矩阵

print(pca_model$x)
# 等价于：
a = scale(data_tr)
a %*% pca_model$rotation

降维后的矩阵

• 每个特征的平均值

print(pca_model$center)

平均值

• 每个特征的标准差，样本标准差，无偏估计

print(pca_model$scale)

标准差