Plotting Indicator Species Analysis over time (indispecies) [closed]

Question

I want to plot the results of an indicator species analysis (package indispecies) based on the time the vegetation plots where sampled to see, if the indicator species change.

How can I get the indicator values out of the ISA for this treatment and the different years and visualize them?

I have this mock dataset with a species abundance matrix, a column ID for the plots, a column year for the year the plot was sampled and a column treatment for three treatments the plots are in. The year is only available for one of the treatments. I only want to plot the indicator species/time for the treatment IBR. The indicator species for the treatments SER and NER are fixed and do not depend on time.

  
## mock dataset
# IDs for the vegetation plots
ID <- paste0("PL", 1:50)

# some sppecies
species <- c(paste0("Species_", 1:100), paste0("Rare_Species_", 1:100))

# three treatments
treatments <- c(rep("IBR", 30), rep("SER", 10), rep("NER", 10))

# generate some species abundance data
generate_abundance <- function(id, treatment) {
  abundance <- runif(200, min = 0, max = 1)
  
  # Introduce some specific and rare and abundant species
  if (treatment == "IBR") {
    abundance[1:5] <- runif(5, min = 80, max = 100)  
    abundance[6:10] <- runif(5, min = 10, max = 20)  
    abundance[11:15] <- runif(5, min = 50, max = 60) 
    abundance[16:20] <- runif(5, min = 0, max = 1)   
    abundance[21:100] <- runif(80, min = 0, max = 1) 
  } else if (treatment == "SER") {
    abundance[11:15] <- runif(5, min = 80, max = 100) 
    abundance[16:20] <- runif(5, min = 10, max = 20)  
    abundance[1:5] <- runif(5, min = 50, max = 60)    
    abundance[6:10] <- runif(5, min = 0, max = 1)     
    abundance[21:100] <- runif(80, min = 0, max = 1)  
  } else if (treatment == "NER") {
    abundance[21:25] <- runif(5, min = 80, max = 100) 
    abundance[26:30] <- runif(5, min = 10, max = 20)  
    abundance[31:35] <- runif(5, min = 20, max = 50)  
    abundance[1:20] <- runif(20, min = 0, max = 1)    
    abundance[36:100] <- runif(65, min = 0, max = 1) 
  }
  
  # Add extremely rare species with mostly zero values
  for (i in 101:200) {
    if (runif(1) < 0.2) {  # 20% chance of occurrence
      abundance[i] <- runif(1, min = 0, max = 1)
    } else {
      abundance[i] <- 0
    }
  }
  
  total <- sum(abundance)
  
  # Abundance is between 100 and 150
  scaling_factor <- runif(1, min = 100, max = 150) / total
  scaled_abundance <- round(abundance * scaling_factor, 1)  
  
  return(scaled_abundance)
}

data <- t(sapply(1:50, function(x) generate_abundance(ID[x], treatments[x])))
df <- as.data.frame(data)

# Adding and ID and a treatment column
df <- cbind(ID = ID, df)
colnames(df)[-1] <- species
df <- cbind(Treatment = treatments, df)

# add year
years <- rep(NA, 50)
ibr_indices <- which(df$Treatment == "IBR")
years_pool <- sample(0:40, 10, replace = TRUE)
years_repeated <- unlist(lapply(years_pool, function(x) rep(x, sample(1:5, 1))))

years[ibr_indices] <- sample(years_repeated, 30)
df <- cbind(Year = years, df)

# ISA
species_data <- df[, 4:203]
treatment_data <- df$Treatment

result <- multipatt(species_data, treatment_data, func = "IndVal.g", control = how(nperm = 999))
summary(result)