筆記內容：
- 用ggplot2繪制scatterplot with confidence ellipses
- envfit的R實現(xiàn)及釋意
- 用ggpubr繪制scatterplot with confidence ellipses
- ggsci的調色板彩蛋！！QAQ
在微生物β-diversity分析中常用距離矩陣(unifrac)做PcoA聚類分析，以觀察不同組間物種構成的差異。結合針對距離矩陣的MANOVA統(tǒng)計檢驗結果，總結出組間物種構成差異的顯著性。通過將各點PcoA結果坐標繪制在二維平面圖上，并加上confidence ellipses進行可視化。
##使用內置數(shù)據(jù)集iris演示
library(ggplot2)
library(ade4)
data(iris)
iris <- data.frame(iris)
iris_sub <- iris[,1:4]
group <- as.factor(iris$Species)
mds <- cmdscale(dist(iris_sub), k = 2, eig = TRUE)
mds_point <- data.frame(mds$points) # 得到各樣本的坐標
colnames(mds_point) <- c('X1','X2')
eig <- mds$eig
用ggplot2繪制scatterplot with confidence ellipses
color <- c(brewer.pal(3,"Set1"))
ggplot(mds_point, aes(x = X1, y = X2, color = group)) +
geom_point(aes(color = group), size = 4, alpha = 0.6) +
stat_ellipse(aes(x = X1, y = X2, fill = group), geom = "polygon", alpha = 1/2, levles = 0.95) +
# geom用于設置填充形狀，alpha設置透明度。不設置則為實心填充，遮蓋橢圓中的點, levels設置confidence ellipses的置信區(qū)間, 在0-1范圍內。levels越小橢圓面積越小，涵蓋的點越集中。
# 不需要填充的時候去掉fill及對fill的補充參數(shù)geom,alpha等即可
scale_fill_manual(values= color) +
scale_color_manual(values = color)
# 顏色可以自己設置，或者直接用scale_color_brewer()
置信橢圓的算法復雜，背后有很多繁雜的數(shù)學原理，在這里不深究，只是了解置信橢圓在本例中用于高維數(shù)據(jù)中，像95%置信區(qū)間一樣，在一組數(shù)據(jù)中隨機抽取一個樣本，其落在這個區(qū)間（域）內的概率為95%。在圖中展現(xiàn)了數(shù)據(jù)集中的范圍。再結合MANOVA統(tǒng)計檢驗結果，能更直觀的看出有差異的組別。
envfit的R實現(xiàn)及釋意
library(vegan)
fit <- envfit(mds, iris_sub,permutations = 999)
fit_val <- scores(fit, display = c("vectors"))
fit_val <- fit_val*vegan::ordiArrowMul(fit_val, fill = 1.5)
## fit的結果output:
## P值表示顯著性，Dim1, Dim2的坐標用于接下來把向量添加到PcoA二維圖中。
> fit$vectors
Dim1 Dim2 r2 Pr(>r)
Sepal.Length 0.48219 0.87607 0.9579 0.001 ***
Sepal.Width -0.11499 0.99337 0.8400 0.001 ***
Petal.Length 0.98013 -0.19836 0.9981 0.001 ***
Petal.Width 0.97852 -0.20615 0.9366 0.001 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Permutation: free
Number of permutations: 999
##ggplot添加向量
ggplot(mds_point, aes(x = X1, y = X2, color = group)) +
geom_point(size = 4, alpha = .6) +
stat_ellipse(aes(x = X1, y = X2, fill = group), geom = "polygon", alpha = 1/2) +
scale_fill_manual(values=color) +
scale_color_manual(values = color) +
geom_segment(data=data.frame(fit_val),
aes(x=0,y=0,xend=Dim1, yend=Dim2),
arrow=arrow(length=unit(0.2,"cm")),
color='black',alpha=1) +
geom_label_repel(data=data.frame(fit_val), aes(Dim1, Dim2, label=rownames(fit_val)),
color='black',alpha=1,
segment.color = 'grey35',
point.padding = unit(0.1,"lines")) +
labs(x = paste("PCoA 1 (", format(100*eig[1]/sum(eig), digits = 4), "%)",sep = ""),
y = paste("PCoA 2 (", format(100*eig[2]/sum(eig), digits = 4), "%)",sep = ""))
每個向量代表一個數(shù)據(jù)類型為連續(xù)的環(huán)境變量，在這里即4中不同的width和length。
每個向量的方向環(huán)境變量發(fā)生變化最快的方向；向量的長度表示環(huán)境變量與分組的關系。其與PcoA繪制在一起的意義為結合關鍵環(huán)境變量解讀，通過向量的方向，長短及顯著性，了解哪些環(huán)境變量的變化會趨向構成組間差異，且這樣的差異有意義。起到篩選關鍵環(huán)境變量的作用。在下圖中，可以解讀為Petal.Length和Petal.Width值越大，樣本越趨近versicolor和virginica這兩個物種類別。
用ggpubr繪制scatterplot with confidence ellipses
library(ggpubr)
mds_point <- cbind(mds_point, iris$Species)
colnames(mds_point)[3] <- "group"
ggscatter(mds_point, x= "X1", y = "X2",
color = "group", palette = "Set1", # 任何存在的palette都可以，不僅僅是brewer.pal中的，可以調用ggsci中的各種sci-fi主題顏色包！
ellipse = TRUE, # 設置是否需要confidence ellipses
mean.point = TRUE, star.plot = TRUE, # 設置confidence ellipses中心是否與所有點連線
ellipse.level = 0.95, # 設置confidence level可以調整橢圓的大小
ggtheme = theme_minimal()) +
labs(x = paste("PCoA 1 (", format(100*eig[1]/sum(eig), digits = 4), "%)",sep = ""),
y = paste("PCoA 2 (", format(100*eig[2]/sum(eig), digits = 4), "%)",sep = "")) +
ggtitle("PcoA")
ggpubr中的ggscatter(), ggboxplot()...等等函數(shù)，在繪制較復雜設置較多的圖時比較方便，參數(shù)一目了然。如果用ggplot2可能會稍繁瑣一些。要根據(jù)實際情況選擇使用哪種。
彩蛋?。。。。。。。。AQ
ggsci包就是各種sci-fi主題的調色盤?。。。∮腥俗隽宋乙恢毕胱龅氖拢。。。。?！
除了辛普森一家，星際迷航，也有各種大牌期刊里圖表的配色，與ggplot2,ggpubr無縫兼容，直接使用scale_color_rickandmorty()及scale_fill_rickandmorty()
參考鏈接：
http://ggplot2.tidyverse.org/reference/stat_ellipse.html
https://cran.r-project.org/web/packages/ggsci/vignettes/ggsci.html
vegan包中對ordination的文檔:
https://cran.r-project.org/web/packages/vegan/vignettes/intro-vegan.pdf
