Without data, it's not possible to see what variables are being plotted. Try working from this
suppressPackageStartupMessages({
library(vegan)
})
## The recommended way of running NMDS (Minchin 1987)
##
data(dune)
# Global NMDS using monoMDS
sol <- metaMDS(dune)
#> Run 0 stress 0.1192678
#> Run 1 stress 0.1922241
#> Run 2 stress 0.1183186
#> ... New best solution
#> ... Procrustes: rmse 0.02027044 max resid 0.06496245
#> Run 3 stress 0.1192678
#> Run 4 stress 0.1922241
#> Run 5 stress 0.192224
#> Run 6 stress 0.2075713
#> Run 7 stress 0.1192679
#> Run 8 stress 0.1183186
#> ... Procrustes: rmse 3.656377e-06 max resid 7.680827e-06
#> ... Similar to previous best
#> Run 9 stress 0.1183186
#> ... Procrustes: rmse 1.743913e-06 max resid 5.607072e-06
#> ... Similar to previous best
#> Run 10 stress 0.1183186
#> ... Procrustes: rmse 5.908868e-06 max resid 1.872976e-05
#> ... Similar to previous best
#> Run 11 stress 0.1192678
#> Run 12 stress 0.1192678
#> Run 13 stress 0.1886532
#> Run 14 stress 0.1192678
#> Run 15 stress 0.1809577
#> Run 16 stress 0.1192679
#> Run 17 stress 0.1808911
#> Run 18 stress 0.1183186
#> ... Procrustes: rmse 8.668069e-06 max resid 2.714218e-05
#> ... Similar to previous best
#> Run 19 stress 0.1808911
#> Run 20 stress 0.1183186
#> ... Procrustes: rmse 2.455886e-06 max resid 8.038281e-06
#> ... Similar to previous best
#> *** Solution reached
sol
#>
#> Call:
#> metaMDS(comm = dune)
#>
#> global Multidimensional Scaling using monoMDS
#>
#> Data: dune
#> Distance: bray
#>
#> Dimensions: 2
#> Stress: 0.1183186
#> Stress type 1, weak ties
#> Two convergent solutions found after 20 tries
#> Scaling: centring, PC rotation, halfchange scaling
#> Species: expanded scores based on 'dune'
plot(sol, type="t")

## Start from previous best solution
sol <- metaMDS(dune, previous.best = sol)
#> Starting from 2-dimensional configuration
#> Run 0 stress 0.1183186
#> Run 1 stress 0.1192678
#> Run 2 stress 0.1183186
#> ... Procrustes: rmse 8.397218e-06 max resid 2.160064e-05
#> ... Similar to previous best
#> Run 3 stress 0.1192678
#> Run 4 stress 0.1192678
#> Run 5 stress 0.1183186
#> ... New best solution
#> ... Procrustes: rmse 2.69204e-06 max resid 8.154677e-06
#> ... Similar to previous best
#> Run 6 stress 0.1192679
#> Run 7 stress 0.1183186
#> ... Procrustes: rmse 6.199676e-06 max resid 2.021768e-05
#> ... Similar to previous best
#> Run 8 stress 0.1183186
#> ... Procrustes: rmse 3.444356e-06 max resid 1.167095e-05
#> ... Similar to previous best
#> Run 9 stress 0.1192678
#> Run 10 stress 0.1192678
#> Run 11 stress 0.1192678
#> Run 12 stress 0.1183186
#> ... Procrustes: rmse 3.667454e-06 max resid 1.242743e-05
#> ... Similar to previous best
#> Run 13 stress 0.1192678
#> Run 14 stress 0.1809578
#> Run 15 stress 0.1192678
#> Run 16 stress 0.2075713
#> Run 17 stress 0.1922241
#> Run 18 stress 0.1192678
#> Run 19 stress 0.1192679
#> Run 20 stress 0.1183186
#> ... Procrustes: rmse 1.418263e-06 max resid 2.796378e-06
#> ... Similar to previous best
#> *** Solution reached
## Use Arrhenius exponent 'z' as a binary dissimilarity measure
sol <- metaMDS(dune, distfun = betadiver, distance = "z")
#> Run 0 stress 0.1067169
#> Run 1 stress 0.1073148
#> Run 2 stress 0.1736237
#> Run 3 stress 0.1814422
#> Run 4 stress 0.107471
#> Run 5 stress 0.1067169
#> ... New best solution
#> ... Procrustes: rmse 2.836292e-06 max resid 7.464074e-06
#> ... Similar to previous best
#> Run 6 stress 0.2131677
#> Run 7 stress 0.1067169
#> ... Procrustes: rmse 2.666581e-06 max resid 6.15736e-06
#> ... Similar to previous best
#> Run 8 stress 0.1073148
#> Run 9 stress 0.1073148
#> Run 10 stress 0.107471
#> Run 11 stress 0.1067169
#> ... Procrustes: rmse 3.798976e-06 max resid 1.251695e-05
#> ... Similar to previous best
#> Run 12 stress 0.186847
#> Run 13 stress 0.1067169
#> ... Procrustes: rmse 8.847931e-06 max resid 2.444281e-05
#> ... Similar to previous best
#> Run 14 stress 0.1067169
#> ... Procrustes: rmse 4.341968e-06 max resid 1.047293e-05
#> ... Similar to previous best
#> Run 15 stress 0.107471
#> Run 16 stress 0.1067169
#> ... Procrustes: rmse 4.059188e-06 max resid 9.843856e-06
#> ... Similar to previous best
#> Run 17 stress 0.1824921
#> Run 18 stress 0.1073148
#> Run 19 stress 0.1073148
#> Run 20 stress 0.1644742
#> *** Solution reached
sol
#>
#> Call:
#> metaMDS(comm = dune, distance = "z", distfun = betadiver)
#>
#> global Multidimensional Scaling using monoMDS
#>
#> Data: dune
#> Distance: beta.z
#>
#> Dimensions: 2
#> Stress: 0.1067169
#> Stress type 1, weak ties
#> Two convergent solutions found after 20 tries
#> Scaling: centring, PC rotation, halfchange scaling
#> Species: expanded scores based on 'dune'
plot(sol$points, col=sol$species[,1],pch=16,cex=4)

Created on 2021-01-15 by the reprex package (v0.3.0.9001)