導讀
Bray-Curtis NMDS降維分析菌群結構。
一倘感、輸入數(shù)據(jù)
1.1 菌屬豐度矩陣
df = data.frame(abs(round(matrix(rnorm(729, 100, 50), 27, 27))))
rownames(df) = paste(rep(paste(rep(c("A", "B", "C"), each=3), rep(c("I","II","III")), sep=""), each=3), 1:3, sep="")
colnames(df) = paste("genus", 1:27, sep="")
1.2 樣品分組
df_group = paste(rep(c("A", "B", "C"), each=9), rep(c("I","II","III"), each=3), sep="")
data.frame(rownames(df), df_group)
二放坏、BC NMDS函數(shù)
library(vegan)
library(ggplot2)
library(ape)
# nmds
nmds_f = function(input, group)
{
nmds = metaMDS(input, distance="bray", k=2)
nmds_point = data.frame(nmds$points)
nmds_point$Category = group
nmds_point$Group = substring(nmds_point[, 3], 1, 1)
nmds_point$Mark = substring(nmds_point[, 3], 2)
# substring(x, start, end) no end is all
write.csv(nmds_point, file="bray_curtis_nmds.csv")
result = ggplot(nmds_point, aes(x=MDS1, y=MDS2, color = Group)) +
geom_point(aes(shape = Mark), size = 4) +
stat_ellipse(level = 0.95, show.legend = F, aes(fill = Category)) +
theme_classic() + # 經典，只有xy軸
labs(title = paste('Bray-Curtis NMDS Stress =', round(nmds$stress, 4))) +
scale_color_manual(values=c("A" = "#EE6363", "B" = "#1E90FF", "C" = "#00C5CD")) +
scale_shape_manual(values=c(1, 2, 0)) +
guides(color = guide_legend(order=1), shape = guide_legend(order=2))
# guides(color=guide_legend(title = NULL)) +
# theme(legend.title=element_blank()) +
# geom_hline(yintercept = 0, size=0.1, linetype="dashed") +
# geom_vline(xintercept = 0, size=0.1, linetype="dashed") +
# geom_text(aes(label = rownames(nmds_point))) +
# theme(panel.grid=element_blank(), panel.background=element_rect(color='black', fill='transparent')) +
# theme_bw() + # 黑框老玛，無背景
# theme(legend.title = element_text(size = 13, face = "bold"), plot.title = element_text(hjust = 0.5))
ggsave(result, filename="bray_curtis_nmds.pdf")
}
三淤年、分析及結果
nmds_f(df, df_group)
ellipse圈圈沒出來，數(shù)據(jù)是模擬的蜡豹，沒事
3.1 結果文件
3.2 結果表
3.3 結果圖
四麸粮、分析舉例
