先上圖

image.png
注:填充顏色由AI處理。
(樣本和分組文件如下)(這里只是展示數(shù)據(jù)形式并不是原數(shù)據(jù)(分組名與代碼中不一樣))

image.png

image.png
vegan分析數(shù)據(jù)
行名樣本,列名物種,內(nèi)容豐度。
#加載R包
library(yyplot)
library(ggbiplot)
library(vegan)
library(ggplot2)
#導(dǎo)入樣本和分組文件
sample <- read.table("sample.txt",sep = "\t",header = T)
group <- read.table("group1.txt",sep = "\t",header = T,row.names = 1)
group <- group[match(rownames(sample),rownames(group)),] #匹配行名,很重要
#nmds分析
nmds1 <- metaMDS(sample, distance = 'bray', k = 2)
summary(nmds1)
#提取數(shù)據(jù)
nmds1.stress <- nmds1$stress
nmds1.point <- data.frame(nmds1$point)
nmds1.species <- data.frame(nmds1$species)
sample_site <- nmds1.point[1:2]
sample_site$names <- rownames(sample_site)
colnames(sample_site)[1:2] <- c('NMDS1', 'NMDS2')
#合并分組數(shù)據(jù)
sample_site <- cbind(sample_site,group)
繪圖
#NMDS圖繪制
nmds_plot <- ggplot() +
geom_point(data = sample_site, aes(NMDS1, NMDS2,color = site , shape = depth), size = 5, alpha = 0.8) + #可在這里修改點(diǎn)的透明度、大小
#scale_shape_manual(values = c(17, 16)) + #可在這里修改點(diǎn)的形狀
#scale_color_manual(values = c('red', 'blue')) + #可在這里修改點(diǎn)的顏色
geom_ord_ellipse(aes(sample_site$NMDS1,sample_site$NMDS2,group= sample_site$site2), ##添加0.8置信橢圓
ellipse_pro = 0.8,linetype=2,size=0.7,color='firebrick')+
geom_ord_ellipse(aes(sample_site$NMDS1,sample_site$NMDS2,color= sample_site$site2,group= sample_site$site2),
ellipse_pro = 0.9,linetype=3,size=1)+ ##添加0.9置信橢圓
theme(panel.grid = element_blank(), panel.background = element_rect(color = 'black', fill = 'transparent')) + #去掉背景
theme(legend.key = element_rect(fill = 'transparent'), legend.title = element_blank()) + #去掉圖例標(biāo)題及標(biāo)簽背景
labs(x = 'NMDS axis1', y = 'NMDS axis2', title = paste('Stress =', round(nmds1$stress, 4))) +
theme(plot.title = element_text(hjust = 0.5))+ #標(biāo)題居中
theme(panel.background = element_blank(),axis.line = element_line(color = "black"))#去上右邊框
nmds_plot
注意?。eom_ord_ellipse()函數(shù)中的參數(shù)需要data$列名 而不是直接用列名
geom_ord_ellipse()函數(shù)參數(shù)中的fill好像只能填充一種顏色。默認(rèn)無填充。
ggplot2自帶的stat_ellipse() 同樣0.8參數(shù)情況下會(huì)繪制不同大小的圈。
用法:
p+stat_ellipse(aes(fill=group),geom="polygon",level=0.8,alpha=0.1,show.lengend=F)
#或者
p+stat_ellipse(aes(color=group),level=0.8,linetype=2,show.lengend=F)
以下是NMDS結(jié)果的粗略評估
#stress記錄了NMDS排序分析的應(yīng)力函數(shù)值,points記錄了各樣本的排序坐標(biāo),
#species記錄了各物種(OTU)的排序坐標(biāo)。
#而在NMDS排序分析中,盡可能選擇較低的應(yīng)力函數(shù)值。一般情況下,應(yīng)力函數(shù)值的值不要大于0.2。
write.csv(nmds1.point, 'nmds.sample.csv')
#簡單繪圖評估
nmds_plot <- nmds1
nmds_plot$species <- {nmds_plot$species}[1:10, ]
stressplot(nmds_plot, main = 'Shepard 圖')
gof <- goodness(nmds_plot)
plot(nmds_plot,type = 't', main = '擬合度')
points(nmds_plot, display = 'sites', cex = gof * 200, col = 'red')
關(guān)鍵詞:nmds,ggplot2,yyplot,gglayer,ggbiplot,github,geom_ord_ellipse
參考網(wǎng)站:
https://guangchuangyu.github.io/cn/2018/01/geom-ord-ellipse/
https://github.com/GuangchuangYu/gglayer#geom_ord_ellipse
https://mp.weixin.qq.com/s/4vw6chPvrz_tLmLXHQ7ThQ
https://mp.weixin.qq.com/s/_XH6u2hPv_JJi_MbGbeoHA
https://blog.csdn.net/qq_38854576/article/details/83024468/