3

我通常使用该函数执行主成分分析,并使用(或仅使用提取)prcomp以一种奇特的方式绘制结果。ggbiplotggplot2pca.obj$x

像这样:

#install_github("vqv/ggbiplot")
library(ggbiplot)
data(iris)
pca.obj <- prcomp(iris[,1:4], center=TRUE, scale.=TRUE)
P <- ggbiplot(pca.obj,
         obs.scale = 1, 
         var.scale=1,
         ellipse=T,
         circle=F,
         varname.size=3,
         var.axes=T,
         groups=iris$Species, #no need for coloring, I'm making the points invisible
         alpha=0) #invisible points, I add them below
P$layers <- c(geom_point(aes(color=iris$Species), cex=5), P$layers) #add geom_point in a layer underneath (only way I have to change the size of the points in ggbiplot)
png(filename="test.png", height=600, width=600)
print(#or ggsave()
    P
)
dev.off()

测试1

但是,现在我面临具有一定数量 NA 的数据,并且我正在使用pcaMethodspca包中的包装函数,应用该方法(一种能够处理少量缺失值的迭代方法)。nipals

pca返回 class 的对象pcaRes,并ggbiplot返回以下错误:

#introduce NAs
iris$Sepal.Length[sample(1:150, 5)] <- NA
iris$Sepal.Width[sample(1:150, 5)] <- NA
iris$Petal.Length[sample(1:150, 5)] <- NA
iris$Petal.Width[sample(1:150, 5)] <- NA
#pca.obj2 <- prcomp(iris[,1:4], center=TRUE, scale.=TRUE) #cannot use prcomp with NAs
#source("https://bioconductor.org/biocLite.R")
#biocLite("pcaMethods")
library(pcaMethods)
pca.obj2 <- pca(iris[,1:4], method="nipals", nPcs=3, center=TRUE, scale.=TRUE)
class(pca.obj2)
ggbiplot(pca.obj2)

ggbiplot(pca.obj2) 中的错误:需要类 prcomp、princomp、PCA 或 lda 的对象

我的问题是:

如何应用ggbiplotpcaRes对象?

如何将此对象转换为prcomp对象?

我可以使用另一个函数而不是ggbiplot接受一个pcaRes对象来获得相同类型的绘图吗?

我应该用变量的平均值替换 NA 值并prcomp像往常一样应用函数吗?

非常感谢!

4

1 回答 1

1

首先,找到一个可以处理 NA 的 PCA 包是件好事。

由于ggbiplot不会接受pcaRes对象,我们可以使用通过获取的数据pcaRes并将其潜入原始prcomp对象中。

显然,您的真实数据已经包含这些NA值,因此我们将从该数据集开始并将它们换成一些虚拟值,以允许我们运行第一个prcomp pca.

iris_na<-iris

iris_na$Sepal.Length[sample(1:150, 5)] <- NA
iris_na$Sepal.Width[sample(1:150, 5)] <- NA
iris_na$Petal.Length[sample(1:150, 5)] <- NA
iris_na$Petal.Width[sample(1:150, 5)] <- NA

iris_dummy<-iris_na

iris_dummy[is.na(iris_dummy)]<-7777 #swap out your NAs with a dummy number so prcomp will run

pca然后我们像你一样运行第一个:

pca.obj <- prcomp(iris_dummy[,1:4], center=TRUE, scale.=TRUE)

该对象有 5 个分量,x(分数)、rotation(载荷)、sdev(标准差)centerscale。虽然我怀疑只有分数和载荷被 使用ggbiplot,但为了确定,我们将它们全部交换掉。

查看分数组件向我们展示了在函数pca.obj$x中计算了四个主组件。prcomp

head(pca.obj$x)

#           PC1        PC2         PC3         PC4
#[1,] -2.656740  0.3176722  0.03763067 -0.04122948
#[2,] -2.688275 -0.1821744  0.19912795  0.07297624
#[3,] -2.862673 -0.1447518 -0.02134749 -0.02462359
#[4,] -2.718294 -0.3189371 -0.03318459 -0.11675762
#[5,] -2.700864  0.3274887 -0.07503096 -0.11347939
#[6,] -2.252918  0.7436711 -0.14611455 -0.08218007

因此,当我们使用 运行下一个 pca 时pcaRes,我们确保指定使用nPcs参数计算 4 个主成分。这里我们使用的是真实数据,其中包含NAs.

pca.obj2 <- pca(iris_na[,1:4], method="nipals", nPcs=4, center=TRUE, scale.=TRUE)

然后只需将pcaRes值换出prcomp值并将其传递给ggbiplot

pca.obj$x<-pca.obj2@scores 

pca.obj$rotation<-pca.obj2@loadings 

pca.obj$sdev<-pca.obj2@sDev

pca.obj$center<-pca.obj2@center

pca.obj$scale<-pca.obj2@scale

P2 <- ggbiplot(pca.obj,
              obs.scale = 1, 
              var.scale=1,
              ellipse=T,
              circle=F,
              varname.size=3,
              var.axes=T,
              groups=iris$Species, 
              alpha=0) 
P2$layers <- c(geom_point(aes(color=iris$Species), cex=5), P2$layers)

在此处输入图像描述

于 2018-04-12T04:53:58.977 回答