rm(list = ls())
library(tidyverse)
library(plyr)
library(ggplot2)
library(ggsci)
pca = prcomp(iris[,1:4]) # 进行PCA计算
var = pca[["sdev"]]^2 # 提取解释度
pc1 = round(var[1]/sum(var) * 100, 2)
pc2 = round(var[2]/sum(var) * 100, 2)
res.df = pca[["x"]] %>%
as.data.frame() %>%
mutate(Species = iris$Species)
# 置信椭圆
ggplot(res.df, aes(PC1, PC2, color = Species, shape = Species)) +
geom_point(size = 2)+
stat_ellipse(level = 0.95)+
scale_color_igv() +
scale_shape_manual(values = c(15,16,17)) +
theme_bw()+
labs(x = paste('PC1(',pc1,'%)',sep = ''),
y = paste('PC2(',pc2,'%)',sep = '')) +
theme(legend.position = c(0.86,0.85),
legend.title = element_blank(),
legend.background = element_blank(),
axis.text = element_text(color = 'black'),
axis.title = element_text(color = 'black'))
# 多边形
group_border = ddply(res.df, 'Species', function(df) df[chull(df[[1]], df[[2]]), ])
ggplot(res.df, aes(PC1, PC2, color = Species, fill = Species)) +
geom_polygon(data = group_border, alpha = 0.4, show.legend = F) +
geom_point(aes(shape = Species), size = 2)+
scale_color_igv() +
scale_shape_manual(values = c(15,16,17)) +
theme_bw()+
labs(x = paste('PC1(',pc1,'%)',sep = ''),
y = paste('PC2(',pc2,'%)',sep = '')) +
theme(legend.position = c(0.86,0.85),
legend.title = element_blank(),
legend.background = element_blank(),
axis.text = element_text(color = 'black'),
axis.title = element_text(color = 'black'))
