Представленные на графике векторы envfit не совпадают с баллами NMDS - PullRequest
1 голос
/ 11 декабря 2019

Я сделал график NMDS и построил свой envfit следующим образом

dataframe для mytable

sites=c("Site A","Site B","Site C","Site D","Site E","Site F","Site 
G","Site H","Site I","Site J","Site K","Site L","Site M","Site N","Site O","Site P","Site Q","Site R","Site S","Site T","Site U")
american.elm=c(41.91,10.11,2.62,5.31,7.51,9.72,17.44,9.06,19.83,30.81,62.6,21.29,20.7,28.68,27.69,34.89,35.65,3.87,12.68,1.58,2.97)
white.birch=c(7.07,15.89,26.77,15.61,14.59,6.33,2.23,11.66,21.49,20.15,7.61,23.29,0,0,0,0,0,0,0,56.09,42.34)
red.oak=c(0,0,0,0,0,0,0,0,0,0,0,6.02,0,0,0,0,0,0,0,0,0.05)
populus.grand=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0.11,0)
beech=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2.36,5.45)
sugar.maple=c(0.49,2.64,3.35,4.6,3.37,2,1.32,4.21,4.13,3.61,0.34,1.2,0,0,0,0,0,0,0,2.19,0.09)

mytable <- data.frame(sites,american.elm,red.oak,populus.grand,beech,sugar.maple)
mytable<-mytable[,2:ncol(mytable)]

Затем

library(vegan)

mytable.NMDS=metaMDS(mytable, distance = "jaccard", k = 4, trymax = 2000, autotransform=FALSE)

plot.mytable<-data.frame(mytable.NMDS$points)
plot.mytable
par(mar=c(3,3,2,5) ,mgp=c(1.8,0.75,0))
plot(plot.mytable$MDS1, plot.mytable$MDS2, pch=16, cex=1, col="black",
     xlab="NMDS1", ylab="NMDS2", cex.lab=1, cex.axis=1, main="", bty="L",
     mai=c(0,0,2,10), xlim=c(-1.5,1.3), ylim=c(-0.9,1))

fit <- envfit(mytable.NMDS, mytable, choices=c(1,2,3))
fit.plot = plot(fit, cex=1.3, col="red", xlim=c(-1.5,1.3), ylim=c(-1.2,1.2),
                xlab="NMDS1",ylab="NMDS2")

Это оценка NMDS деревьев

fit
# Table of the NMDS score of the trees
Trees=c("american.elm","red.oak","populus.grand","beech","sugar.maple")
Tree.NMDS1=c(-0.76538,-0.1533,0.36065,0.25411,0.49583)
Tree.NMDS2=c(-0.27961,0.06605,-0.51345,-0.79497,0.84299)
Tree.NMDS.scores=data.frame(Trees,Tree.NMDS1,Tree.NMDS2)
# Overlay the NMDS score on the plot
points(Tree.NMDS.scores$Tree.NMDS1,Tree.NMDS.scores$Tree.NMDS2,
       col="red", pch=16)

Я хотел бы знать, почему стрелки конца вектора не совпадают с оценками NMDS, полученными функцией envfit()?

1 Ответ

0 голосов
/ 12 декабря 2019

Значения, которые вы видите в таблице, представляют собой стандартизированные коэффициенты из линейной регрессии, используемые для проекции векторов в ординацию. Это направления для стрелок единичной длины. При построении графика мы масштабируем эти стрелки по квадратному корню из их корреляции. Как таковые стрелки с небольшими корреляциями представлены более короткими стрелками, чем стрелки с более сильными корреляциями. Вы можете получить эти масштабированные длины стрелок, используя метод scores():

> scores(fit, "vectors")
                    NMDS1       NMDS2      NMDS3
american.elm  -0.73129278 -0.26985224 -0.5479775
red.oak       -0.06624995  0.03042562  0.4270764
populus.grand  0.21774166 -0.31045377  0.4862402
beech          0.22772624 -0.70982231  0.4990966
sugar.maple    0.33541356  0.56604306 -0.1245767

Обратите внимание, однако, что это не фактические координаты стрелок на графике. Поскольку это просто направления, даже после того, как мы масштабировали длину отдельных стрелок по степени их корреляции с осями ординации, мы можем масштабировать все стрелки на одну и ту же величину, чтобы они лучше заполняли пространство графика.

Это все объясняется в ?envfit. Вот соответствующий раздел:

Вывод на печать непрерывных переменных (векторов) дает направляющие косинусы, которые являются координатами голов векторов единичной длины. В plot они масштабируются по их корреляции (квадратный корень столбца r2), так что слабые предикторы имеют более короткие стрелки, чем сильные предикторы. Вы можете увидеть масштабированные относительные длины, используя команду scores. Стрелки plot ted (и масштабированные) дополнительно подгоняются к текущему графику с использованием постоянного множителя: это сохранит относительную r2 масштабированную длину стрелок, но попытается заполнить текущий график. Вы можете увидеть множитель, используя ordiArrowMul(result_of_envfit), и установить его с аргументом arrow.mul.

Если мы следуем совету, который мы видим:

> ordiArrowMul(fit)
[1] 1.031244

Подразумевается, что мы умножаем масштабировали стрелки примерно на 3%.

> scrs <- scores(fit, "vectors", choices = 1:2)
> scrs * ordiArrowMul(fit)
                    NMDS1      NMDS2
american.elm  -0.73819522 -0.2723993
red.oak       -0.06687526  0.0307128
populus.grand  0.21979686 -0.3133840
beech          0.22987567 -0.7165221
sugar.maple    0.33857942  0.5713858

Собрав все это вместе с вашим кодом построения, вот как мы получаем стрелки, нарисованные plot.envfit:

plot(mytable.NMDS, display = "sites", type = "n")
points(mytable.NMDS, display = "sites", pch = 19, col = "black")
plot(fit, col = 'red')

## add the locations of arrow heads as blue points to see if the correspond
points(scrs * ordiArrowMul(fit), col = "blue")

И это то, что мы получаем

enter image description here


Обратите внимание, как я использовал существующие методы plot и функции экстрактора, такие как scores() для работы с объектами vegan и построения участков из составных частей. Делая это таким образом, вы избегаете i) необходимости печатать информацию, которая уже доступна для вас, и ii) становиться укушенной, когда мы меняем внутреннее представление веганских объектов или когда значения, хранимые внутри, являются действительно рабочими данными, которые требуют последующей трансформации / обработкичтобы получить действительные или интерпретируемые значения. По возможности избегайте использования $, чтобы рыться в предметах.

...