orders <- read.csv("condensed_order.csv", header = TRUE)

鏈接: https://pan.baidu.com/s/1dQIwcr_gE7v6xFJVKpJ9Wg?pwd=xkeh 提取碼: xkeh 復(fù)制這段內(nèi)容后打開百度網(wǎng)盤手機(jī)App，操作更方便哦

head(orders)

#?First load the vegan package

library(vegan)

nmds_results <- metaMDS(comm = orders[ , 4:11],? # Define the community data?

? ? ? ? ? ? ? ? ? ? ? ? distance = "bray",? ? ? ?# Specify a bray-curtis distance

? ? ? ? ? ? ? ? ? ? ? ? try = 100)? ? ? ? ? ? ? ?# Number of iterations?

nmds_results

##?

## Call:

## metaMDS(comm = orders[, 4:11], distance = "bray", try = 100)?

##?

## global Multidimensional Scaling using monoMDS

##?

## Data:? ? ?wisconsin(sqrt(orders[, 4:11]))?

## Distance: bray?

##?

## Dimensions: 2?

## Stress:? ? ?0.1756999?

## Stress type 1, weak ties

## Two convergent solutions found after 100 tries

## Scaling: centring, PC rotation, halfchange scaling?

## Species: expanded scores based on 'wisconsin(sqrt(orders[, 4:11]))'

#第二步， jump into plot

library(ggplot2)

library(viridis)

#First create a data frame of the scores from the individual sites.

# This data frame will contain x and y values for where sites are located.

data_scores=as.data.frame(scores(nmds_results, display = "sites"))

# Now add the extra aquaticSiteType column?

data_scores <- cbind(data_scores, orders[, 14])

colnames(data_scores)[3] <- "aquaticSiteType"??

# Next, we can add the scores for species data

species_scores <- as.data.frame(scores(nmds_results, "species"))

# Add a column equivalent to the row name to create species labels

species_scores$species <- rownames(species_scores)

# Now we can build the plot!

ggplot() +

? geom_text(data = species_scores, aes(x = NMDS1, y = NMDS2, label = species),

? ? ? ? ? ? alpha = 0.5, size = 8) +?

? geom_point(data = data_scores, aes(x = NMDS1, y = NMDS2,?

? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?color = aquaticSiteType), size = 3) +

? scale_color_manual(values = inferno(15)[c(3, 8, 11)],

? ? ? ? ? ? ? ? ? ? ?name = "Aquatic System Type") +

? annotate(geom = "label", x = -1, y = 1.25, size = 6,

? ? ? ? ? ?label = paste("Stress: ", round(nmds_results$stress, digits = 3))) +

? theme_minimal() +

? theme(legend.position = "top",

? ? ? ? text = element_text(size = 16))

# 調(diào)整

ggplot() +

? geom_text(data = species_scores, aes(x = NMDS1, y = NMDS2, label = species),

? ? ? ? ? ? alpha = 0.6, size = 6) + #alpha字體顏色深度

? geom_point(data = data_scores, aes(x = NMDS1, y = NMDS2,?

? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?color = aquaticSiteType), size = 3) +

? scale_color_manual(values = inferno(15)[c(3, 8, 11)],

? ? ? ? ? ? ? ? ? ? ?name = "Aquatic System Type") +

? annotate(geom = "label", x = -1, y = 1.25, size = 6,

? ? ? ? ? ?label = paste("Stress: ", round(nmds_results$stress, digits = 3))) +

? theme_minimal() +

? theme(legend.position = "top",

? ? ? ? text = element_text(size = 16))+

? scale_colour_manual(values = c("orange", "steelblue"))+

? xlim(-2,1)+ ylim(-1.5,1.5)

Code and data source: https://cougrstats.wordpress.com/2019/12/11/non-metric-multidimensional-scaling-nmds-in-r/