导读

非度量多维尺度分析(nonmetric multidimensional scaling, NMDS)，是一种简介的梯度分析方法，也是基于距离或者相异性矩阵。与其它主要用于最大化变异和一致性的方法不一样，NMDS是一种排序方法。这对于距离缺失的数据有优势，只要先办法确定对象之间的位置关系，便可以进行NMDS分析。NMDS的计算过程会以代码的形式贴在下方，供大家参考

数据格式

NMDS与PCoA一样，NMDS可以基于任何类型距离、相异性矩阵对象(样方)进行排序。当然也可以是原始数据矩阵。这里我用的是weighted unifrac距离矩阵数据

图1 weighted 距离矩阵

分析代码

NMDS排序分析可以通过生态学分析R包vegan中的metaMDS()函数实现。因为输入metaNMDS()的数据可以是原始数据矩阵，也可以是距离矩阵，这里拿上面列举的数据做示范。

运行NMDS分析

rm(list = ls())

path

setwd(path)

library(vegan) ## 加载包

weight_dm

sep = "\t",check.names = F) ## 加载weighted距离矩阵数据

meta_info

#group

set.seed(1234) ## 设置随机种子，以便结果可以重复

weight_nmds

weight_nmds$stress # 压力值，一般小于0.1比较好，但是也要根据所选择的的主成分数目来看

图2 strees的scree图，stress随着主成分数的增加而减少

图片来源GUSTA ME[2]

作图

tiff("NMDS.tiff",width=1000,height=1000)

p

cex=1.5, font=1, pch=5) ##NMDS作图###

legend('topleft', legend=levels(meta_info$group), col=c("purple","green","blue","red"),pch=5,cex=1.5,box.lty =1)###NMDS添加Legend###

with(weight_nmds,ordiellipse(weight_nmds, meta_info$group,display = "sites", kind = "se", conf = 0.95, lwd=2,cex=0.8,lty=1,col=c("purple","green","blue","red"),font=2,label = FALSE))

#dev.off()

#env

#head(env)

#ef1

#p

#ef1

dev.off()

图3 常规方法作图

site.scores

site.scores

head(site.scores)

library("FactoMineR")

library("factoextra")

p

geom_point(aes(shape=group,color=group,fill=group),size=3)

p

p

panel.grid.minor = element_blank(),

legend.position="top",

axis.title.x = element_text(size = 16),

axis.text.x = element_text(angle=0,color="black",vjust = 0.95,hjust = 0.95,size=12),

axis.title.y = element_text(size=16),

axis.text.y = element_text(size=14,color="black"),

strip.text.x = element_text(size=18),

legend.text = element_text(size=14),

legend.title = element_text(size=16))

p

p

ggsave("NMDS.tiff", height=8, width=8, units="in")

图4 ggplot2作图

参考

摘要：beta多样性分析网络上已经有很多，大部分出图都是基于R的基础包，很难看，而且很多图形不够完善，没有做方差分析，更没有将其直接输到图上，本文主要做NMDS分析并做一张完善的高质量图片，提取stress值，推荐适合NMDS结果的差异分析并通过命令展示在图形上，最后加上置信区间椭圆。

NMDS分析，网络上已近有很多相关教程分享其原理，与其他排序(PCA、PCoA、CCA、RDA)

方法的不同之处，简单来讲NMDS也是一种使用物种组成数据的排序称作非限制性排序；NMDS基于距离算法，优于PCA、PCoA、CCA、RDA的地方在于当样本或者物种数量过多的时候使用NMDS会更加准确；

基于R语言，开始作图：安装package，两个重要的包

install.packages("vegan")#微生态常用包，一定安装上

install.packages("ggplot2")

调用package

library("vegan")

#设定工作路径

setwd("E:/马兄作图")

#读取数据，一份otu.table文件和一份分组信息文件

design

otu=read.table("16s_otu_tablejidai.txt",row.names=1, header=T,sep="\t")

#如果数据调整为列名是OTU，行名是样本名

otu=t(otu)

#查看一下

head(otu)

#标准化，当然method有很多，可以跟换?decostand查看其它method

vare.hel

#计算矩阵

vare.dis

#使用NMDS的方法

vare.mds

到此分析过程就算完成了，使用基础包plot()作图当然已不是我们最佳的选择了，图形很基础，我就不展示了，下面基于ggplot做一张完善的图#首先提取前两轴坐标

point = scores(vare.mds)

#将分组文件和得分文件合并

index = merge(design, point,by="row.names",all=F)

计算Stress值

Stress值是反映模型合适程度的指标，NMDS会多次打乱数据计算Stress值，知道找到最合适的模型，也就是最低的Stress值；理想状况下，Stress值为0，一般Stress值低于0.1较为合理(本数据这个值偏高一些)

vare.mds

结果输出中寻找下面stress结果：

#显着性检验；anosim本质是基于排名的算法更加适合NMDSanosim.result

summary(anosim.result)

结果输出，得到下面这两个值：

#tiff输出图形，适合大部分出版刊物，入门级别分辩率300,18*14的长宽；

tiff(file="beta_bray_NMDS.tif", res = 300, compression ="none", width=180,height=140,units= "mm")

#开始出图，将上面得到的三个指标在图中更换stress，R，p，不多说，代码如下：

library("ggplot2")

p = ggplot(index, aes(x=NMDS1, y=NMDS2, color=SampleType)) +

geom_point(alpha=.7, size=2) +

labs(x=paste("NMDS1"),

y=paste("NMDS1"),

title="")#置信区间当然要加上，有三种方式，线条类型也可以更改

p+stat_ellipse(type = "t", linetype = 2)+

annotate("text",x=-1.07,y=-1,parse=TRUE,size=4,label=""R:"*0.7031",family="serif",fontface="italic",colour="darkred")+

annotate("text",x=-1.1,y=-0.9,parse=TRUE,size=4,label=""p:"*0.001",family="serif",fontface="italic",colour="darkred")+

annotate("text",x=-1,y=-0.78,parse=TRUE,size=4,label=""Stress:"*0.1208",family="serif",fontface="italic",colour="darkred")dev.off()

成图展示：

当然这还达不到对于拥有轻微强迫症的我对于图行要求：

请看下一篇文章

用R语言做非度量多维尺度分析(NMDS)

2018-11-24

原创 2017-03-10 OmicShare.com 基迪奥生物 基迪奥生物

微信号 gene-denovo

功能介绍 广州基迪奥生物官方公众平台，小圆每天分享各种生信软件使用技巧， 搜罗好玩的科研生活，定期还有掉节操的生物界八卦分享，让科研变得有意思咯~

NMDS(Non-metric Multidimensional scaling)简单来讲就是降维，比如你有20个样本，首先要计算样本与样本间的距离。在微生物分析过程中，这个距离的选择有很多种。

比如有：manhattan、euclidean、canberra、bray、kulczynski、jaccard、gower、altGower、morisita、horn、mountford、raup、binomial、chao、cao、mahalanobis。默认是bray。常用的其实还有jaccard。

非度量多维尺度法是一种将多维空间的研究对象(样本或变量)简化到低维空间进行定位、分析和归类，同时又保留对象间原始关系的数据分析方法。适用于无法获得研究对象间精确的相似性或相异性数据，仅能得到他们之间等级关系数据的情形。其基本特征是将对象间的相似性或相异性数据看成点间距离的单调函数，在保持原始数据次序关系的基础上，用新的相同次序的数据列替换原始数据进行度量型多维尺度分析。

换句话说，当资料不适合直接进行变量型多维尺度分析时，对其进行变量变换，再采用变量型多维尺度分析，对原始资料而言，就称之为非度量型多维尺度分析。其特点是根据样品中包含的物种信息，以点的形式反映在多维空间上，而对不同样品间的差异程度，则是通过点与点间的距离体现的，最终获得样品的空间定位点图。

用R语言实现非度量型多维尺度分析的过程如下：

library(vegan)

data(varespec)

data(varechem)

## The variables are automatically scaled

rankindex(varechem, varespec)

euc    man      gow     bra      kul

0.2396330 0.2735087 0.2288358 0.2837910 0.2839834

#结果显示kul最优，但是其实bra指数也挺好。

> vare.dis

[1] 0.1117280 0.9539592

> require(MASS)

#载入需要的程辑包：MASS

> vare.mds0

initial  value 18.026495

iter   5 value 10.095483

final  value 10.020469

converged

> ordiplot(vare.mds0, type = "t")

Warning message:

In ordiplot(vare.mds0, type = "t") : Species scores not available

> stressplot(vare.mds0, vare.dis)

vare.mds

vare.mds

plot(vare.mds, type = "t")

参考文献：

【1】A microbita signature associated with experimental food allergy promotes allergic senitization and anaphylaxis. The Journal of Allergy and Clinical Immunology.Volume 131, Issue 1 , Pages 201-212, January 2013.

阅读

精选留言

该文章作者已设置需关注才可以留言

该文章作者已设置需关注才可以留言

加载中

以上留言由公众号筛选后显示

了解留言功能详情

微信扫一扫

关注该公众号

免责声明：本文仅代表文章作者的个人观点，与本站无关。其原创性、真实性以及文中陈述文字和内容未经本站证实，对本文以及其中全部或者部分内容文字的真实性、完整性和原创性本站不作任何保证或承诺，请读者仅作参考，并自行核实相关内容。

http://www.pinlue.com/style/images/nopic.gif