【问题标题】:How does R calculate the PCA ellipses?R如何计算PCA椭圆?
【发布时间】:2021-11-05 22:55:48
【问题描述】:

R 如何知道将置信椭圆放置在 PCA 图的位置? 我有一个使用 iris 数据集的最小代码:

library(factoextra)
a<-data.matrix(iris[-5])
b<-prcomp(a, scale. = TRUE, center = TRUE)
fviz_pca_ind(b,
             col.ind = iris$Species,
             addEllipses = TRUE)

我知道我可以使用b$x 找到绘图坐标。我也知道我可以用b$center 找到集群中心。如何从数据中重新导出椭圆?

【问题讨论】:

    标签: r pca


    【解决方案1】:

    如果您在谈论如何,it eventually 会调用 ggplot2::stat_ellipse

    如果你想要坐标,就像其他ggplot对象一样,你可以用ggplot_build提取数据

    library(factoextra)
    a<-data.matrix(iris[-5])
    b<-prcomp(a, scale. = TRUE, center = TRUE)
    p <- fviz_pca_ind(b,
                 col.ind = iris$Species,
                 addEllipses = TRUE)
    
    ell <- ggplot2::ggplot_build(p)$data[[2]]
    
    head(ell)
    #    colour    fill         x           y group PANEL size linetype alpha
    # 1 #F8766D #F8766D -1.697756 -0.06395559     1     1  0.5        1   0.1
    # 2 #F8766D #F8766D -1.701694  0.22197334     1     1  0.5        1   0.1
    # 3 #F8766D #F8766D -1.713449  0.50017215     1     1  0.5        1   0.1
    # 4 #F8766D #F8766D -1.732842  0.76642364     1     1  0.5        1   0.1
    # 5 #F8766D #F8766D -1.759579  1.01669171     1     1  0.5        1   0.1
    # 6 #F8766D #F8766D -1.793255  1.24718254     1     1  0.5        1   0.1
    
    p + geom_point(aes(x, y, color = factor(group)), data = ell, size = 4)
    

    【讨论】:

    • 哈哈!很高兴看到我不是唯一一个将函数追溯到其源头的人!不过花了一段时间。
    【解决方案2】:

    如果你一直追踪代码,你会发现省略号只是geom_polygonsstat = "ellipse"创建的,即它们是由ggplot中的stat_ellipse计算的。

    我们可以通过仅使用基数 R 和 ggplot 重新创建绘图来展示这一点。以下是一个完全可重现的示例:

    library(ggplot2)
    
    b <- prcomp(iris[-5], scale. = TRUE, center = TRUE)
    df <- as.data.frame(predict(b)[,1:2])
    df$Species <- iris$Species
    
    
    ggplot(df, aes(PC1, PC2, color = Species)) + 
      geom_point() +
      theme_bw() +
      geom_polygon(stat = "ellipse", aes(fill = Species), alpha = 0.3)
    

    最终,stat_ellipse 从与cars::dataEllipse 相同的方法获取数据,所以如果你想要椭圆的原始坐标,你可以这样做:

    e <- car::dataEllipse(df$PC1, df$PC2, df$Species)
    
    

    并像这样获得第 95 个百分位法线数据椭圆坐标:

    e$setosa$`0.95`
    #>               x           y
    #>  [1,] -2.167825  2.06328716
    #>  [2,] -2.104642  2.04546589
    #>  [3,] -2.043166  1.99227221
    #>  [4,] -1.984331  1.90451250
    #>  [5,] -1.929028  1.78351710
    #>  [6,] -1.878095  1.63112017
    #>  [7,] -1.832305  1.44963190
    #>  [8,] -1.792351  1.24180347
    #>  [9,] -1.758839  1.01078534
    #> [10,] -1.732278  0.76007952
    #> [11,] -1.713069  0.49348644
    #> [12,] -1.701504  0.21504739
    #> [13,] -1.697759 -0.07101678
    #> [14,] -1.701889 -0.36036963
    #> [15,] -1.713833 -0.64862486
    #> [16,] -1.733410 -0.93141283
    #> [17,] -1.760322 -1.20444675
    #> [18,] -1.794162 -1.46358770
    #> [19,] -1.834417 -1.70490738
    #> [20,] -1.880476 -1.92474763
    #> [21,] -1.931641 -2.11977588
    #> [22,] -1.987137 -2.28703571
    #> [23,] -2.046123 -2.42399164
    #> [24,] -2.107703 -2.52856754
    #> [25,] -2.170946 -2.59917816
    #> [26,] -2.234892 -2.63475311
    #> [27,] -2.298571 -2.63475311
    #> [28,] -2.361018 -2.59917816
    #> [29,] -2.421288 -2.52856754
    #> [30,] -2.478465 -2.42399164
    #> [31,] -2.531684 -2.28703571
    #> [32,] -2.580138 -2.11977588
    #> [33,] -2.623091 -1.92474763
    #> [34,] -2.659894 -1.70490738
    #> [35,] -2.689988 -1.46358770
    #> [36,] -2.712917 -1.20444675
    #> [37,] -2.728333 -0.93141283
    #> [38,] -2.736002 -0.64862486
    #> [39,] -2.735809 -0.36036963
    #> [40,] -2.727757 -0.07101678
    #> [41,] -2.711966  0.21504739
    #> [42,] -2.688678  0.49348644
    #> [43,] -2.658244  0.76007952
    #> [44,] -2.621126  1.01078534
    #> [45,] -2.577888  1.24180347
    #> [46,] -2.529183  1.44963190
    #> [47,] -2.475751  1.63112017
    #> [48,] -2.418401  1.78351710
    #> [49,] -2.358004  1.90451250
    #> [50,] -2.295473  1.99227221
    #> [51,] -2.231758  2.04546589
    #> [52,] -2.167825  2.06328716
    

    reprex package (v2.0.0) 于 2021-11-05 创建

    【讨论】:

      猜你喜欢
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 2012-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      相关资源
      最近更新 更多