r - biplot princomp 和 biplot 中的颜色

Question

如何使用 princomp 和 biplot 为 biplot 中的 iris 数据集的种类着色不同的颜色。最好的

data(iris)


fit <- princomp(iris[,c(1:3)], cor=TRUE)
biplot(fit)

Answer 1

不幸的是，它不支持它。您必须编写自己的 biplot 函数并添加可能为每个样本赋予不同的颜色，源代码非常简单。

https://github.com/SurajGupta/r-source/blob/master/src/library/stats/R/biplot.R

或者使用更现代的功能，例如 autoplot

autoplot( fit, data=iris, colour="Species", loadings=TRUE )

如此处所示：

https://cran.r-project.org/web/packages/ggfortify/vignettes/plot_pca.html

Answer 2

该函数用于渲染双标图，stats:::biplot.princomp并且stats:::biplot.default不允许不同点使用多种颜色或不同颜色。最简单的解决方案是使用一个包，例如另一个答案中提到的 ggfortify ：

library(ggfortify)
autoplot( fit, data=iris, colour="Species", loadings=TRUE )

或者factoextra：

library(factoextra)
fviz_pca_biplot(fit, col.ind = iris$Species)

最后一个选项，是你重写 biplot 函数，如下所示，col1作为数据点的颜色向量，col2作为载荷的颜色：

biplot_col = function (x, y, var.axes = TRUE,col1,col2, cex = 0.8, 
    xlabs = NULL, ylabs = NULL, expand = 1, xlim = NULL, ylim = NULL, 
    arrow.len = 0.1, main = NULL, sub = NULL, xlab = NULL, ylab = NULL, 
    ...) 
{
    n <- nrow(x)
    p <- nrow(y)
    xlabs <- as.character(1L:n)
    dimnames(x) <- list(xlabs, dimnames(x)[[2L]])
    ylabs <- dimnames(y)[[1L]]
    ylabs <- as.character(ylabs)
    dimnames(y) <- list(ylabs, dimnames(y)[[2L]])
    
    unsigned.range <- function(x) c(-abs(min(x, na.rm = TRUE)), 
        abs(max(x, na.rm = TRUE)))
    rangx1 <- unsigned.range(x[, 1L])
    rangx2 <- unsigned.range(x[, 2L])
    rangy1 <- unsigned.range(y[, 1L])
    rangy2 <- unsigned.range(y[, 2L])
    if (missing(xlim) && missing(ylim)) 
        xlim <- ylim <- rangx1 <- rangx2 <- range(rangx1, rangx2)
    else if (missing(xlim)) 
        xlim <- rangx1
    else if (missing(ylim)) 
        ylim <- rangx2
    ratio <- max(rangy1/rangx1, rangy2/rangx2)/expand
    on.exit(par(op))
    op <- par(pty = "s")
    if (!is.null(main)) 
        op <- c(op, par(mar = par("mar") + c(0, 0, 1, 0)))
    plot(x, type = "n", xlim = xlim, ylim = ylim, col = col1, 
        xlab = xlab, ylab = ylab, sub = sub, main = main, ...)
    text(x, xlabs, cex = cex[1L], col = col1, ...)
    par(new = TRUE)
    dev.hold()
    on.exit(dev.flush(), add = TRUE)
    plot(y, axes = FALSE, type = "n", xlim = xlim * ratio, ylim = ylim * 
        ratio, xlab = "", ylab = "", col = col1, ...)
    axis(3, col = col2, ...)
    axis(4, col = col2, ...)
    box(col = col1)
    text(y, labels = ylabs, cex = cex[2L], col = col2, ...)
    if (var.axes) 
        arrows(0, 0, y[, 1L] * 0.8, y[, 2L] * 0.8, col = col2, 
            length = arrow.len)
    invisible()
}

然后像这样绘制：

lam <- fit$sdev[1:2]
lam <- lam * sqrt(fit$n.obs)
scores <- fit$scores
species2col = c("#c15050","#d97642","#d49d42")
names(species2col) = unique(iris$Species)
col1 = species2col[as.character(iris$Species)]
col2 = "#693c72"

par(mar=rep(2.2,4))
biplot_col(t(t(scores[,1:2])/lam), t(t(fit$loadings[,1:2]) * lam),
           col1 = col1, col2 = col2)