# 引入包
install_if_missing <- function(pkg , BiocM=FALSE ) {
  if (!require(pkg, character.only = TRUE)) {
    if( BiocM ){
      BiocManager::install(pkg ,update=FALSE )
    }else{
      install.packages(pkg, dependencies = TRUE)
    }
    suppressPackageStartupMessages( library(pkg, character.only = TRUE) )
  }
}
required_packages <- c("ggpubr")
invisible(lapply(required_packages, install_if_missing))

sangon <- c('#CA0000','#0087FF','#00BA1D','#CF00FF','#00DBE2','#FFAF00','#0017F4','#006012','#6F1A1A','#B5CF00','#009E73','#F0E442','#D55E00','#CC79A7','#FF8A8A','#AA6400','#50008A','#00FF58','#E175FF','#FFCC99','#009999','#CC0066','#C0C0C0','#666600','#505050','#56B4E9','#0072B2','#33FF33','#99004c','#CCFF99','#660066','#9370DB','#D8BFD8','#BC8F8F','#2F4F4F','#FF6347','#877878','#CD5C5C','#FF0000','#00FF00','#000080','#0D7339','#3D90D9','#D98236','#DF514F','#8CCDA3','#F5E2D7','#FACBBE','#EA4885','#E69F00')

name <- 'genus' 
data <- read.table('genus_abundance.xls', header=T, sep='\t', check.names=F, row.names=1, quote="")
sample_list <- read.table('Sample_list.txt', sep='\t', header=F, check.names=F, colClasses='character', na.strings='')

data <- data[, match(sample_list[, 2], colnames(data)), drop=F]
tdata <- t(data)

pca <- prcomp(tdata, scale=FALSE)
imp <- summary(pca)$importance

write.table(pca$x, paste0(name, '_pca_axis.xls'), sep='\t', col.names=NA, quote=F)
write.table(pca$rotation, paste0(name, '_pca_rotation.xls'), sep='\t', col.names=NA, quote=F)
write.table(imp, paste0(name, '_pca_importance.xls'), sep='\t', col.names=NA, quote=F)

xlab <- paste0('PC1 (', round(imp[2,1]*100,2), '%)')
ylab <- paste0('PC2 (', round(imp[2,2]*100,2), '%)')

si <- as.data.frame(pca$x)
legend2 <- 'right' 
color <- 'Samples'
si$Samples <- rownames(si) 
si$Samples <- factor(si$Samples, levels=unique(si$Samples))
ff <- as.numeric(si$Samples)
legend3 <- unique(si$Samples)
si$Groups <- sample_list[match(rownames(si), sample_list[,2]), 1]
si$Groups <- factor(si$Groups, levels=unique(si$Groups))
ff <- as.numeric(si$Groups)
legend3 <- unique(si$Groups)
color <- 'Groups'
n <- nlevels(factor(si$Groups))
mean.point <- 'T' 
ellipse <- 'T'
ellipse.border.remove <- 'F'
label <- 'Samples'
pointch <- rep(15:18, 15)
if (n >= 50) {
  sangon <- rep('#CA0000', n)
  pointch <- rep(15, n)
}



p <- ggscatter(data.frame(si), x='PC1', y='PC2', color=color, shape=color, sangon=sangon, ellipse=ellipse, ellipse.level=0.8, ellipse.type='confidence',  ellipse.alpha=0.1, mean.point=FALSE, star.plot=FALSE, ellipse.border.remove=ellipse.border.remove, repel=TRUE, label=label, legend=legend2, xlab=xlab, ylab=ylab, show.legend.text = FALSE) + scale_shape_manual(values=pointch[1:n])
p <- p +  border() + geom_vline(xintercept=0, linetype=2, colour='grey') + geom_hline(yintercept=0, linetype=2, colour='grey')

p