5  Co-occurance simulations

# Load required packages
# Data manipulation
library(here)
library(readxl)
library(tidyverse)
library(DT)
# Modeling
library(mrIML)
library(tidymodels)
library(future)
library(finetune)
library(flashlight)
# Plotting
library(igraph)
library(ggnetwork)
library(cowplot)
library(patchwork)

library(seqtime)
library(Hmsc)

# Load custom functions
source(here("R-functions", "select_variables.r"))
source(here("R-functions", "plot_network.r"))

# Set random seed
set.seed(109)

# Set up parallel processing
#n_cores <- parallel::detectCores()
#plan("multisession", workers = n_cores - 2)
plan("sequential")

# Load simulated data
source(
  here("sim-study", "symetric-sims.r")
)

Generate a co-ocurance matrix, M.

network_size <- 20
k_average <- 4

simulated_network <- generateM_specific_type(
  nn = network_size,
  k_ave = k_average,
  type.network = "random",
  type.interact = "random",
  interact.str.max = 0.4,
  mix.compt.ratio = 0.5
)

Warning: `erdos.renyi.game()` was deprecated in igraph 0.8.0.
ℹ Please use `sample_gnm()` instead.

Warning: `get.adjacency()` was deprecated in igraph 2.0.0.
ℹ Please use `as_adjacency_matrix()` instead.

Warning: `get.edgelist()` was deprecated in igraph 2.0.0.
ℹ Please use `as_edgelist()` instead.

true_network <- simulated_network[[1]]
M <- simulated_network[[2]]

plot population dynamics

# initial abundance
y <- rpois(network_size, lambda = 100)
# growth rates
r <- runif(network_size)
res <- glv(network_size, M, r, y)
tsplot(10 * res[, 20:1000], time.given = T)

Generate dataset on species abundance using the GLV model and convert to presence-absence data.

data <- generateDataSet(
  900,
  M,
  count = network_size * 10000,
  mode = 4
) %>%
  t() %>%
  as.data.frame() %>%
  mutate(
    across(everything(), ~ ifelse(. < 1, yes = 0, no = 1))
  )

Prepare data and fit mrIML model.

Y <- filterRareCommon(data, lower = 0.01, higher = 0.99)
X1 <- Y

model_rf <- rand_forest(
  trees = 100,
  mode = "classification",
  mtry = tune(),
  min_n = tune()
) %>%
  set_engine("randomForest")

yhats_rf_sim <- mrIMLpredicts(
  Y = Y,
  X = NULL,
  X1 = X1,
  Model = model_rf,
  prop = 0.7,
  k = 5,
  racing = TRUE
)

i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"

→ A | warning: No control observations were detected in `truth` with control level '1'.

There were issues with some computations   A: x1

There were issues with some computations   A: x10

i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"

Perform bootstrapping.

bs_sim <- mrBootstrap(yhats_rf_sim)

Extract co-occurance network.

assoc_net <- mrCoOccurNet(bs_sim)

Plot the network alongside the truth.

mrIML_mat <- assoc_net %>%
  filter(mean_strength > 0.1) %>%
  mutate(
    mean_strength_dir = ifelse(
      direction == "negative",
      yes = -mean_strength,
      no = mean_strength
    ),
    across(contains("taxa"), ~ sub("^sp", "", .))
  ) %>%
  graph_from_data_frame(
    directed = FALSE,
    vertices = sub("^sp", "", names(Y))
  ) %>%
  get.adjacency(attr = "mean_strength_dir", sparse = FALSE)

#plots
g_mrIML <- graph_from_adjacency_matrix(
  mrIML_mat,
  mode     = "undirected",
  weighted = T,
  diag     = F
)

edge_colors_mrIML <- ifelse(E(g_mrIML)$weight < 0, yes = "red", no = "blue")
edge_widths_mrIML <- abs(E(g_mrIML)$weight) * 10 # Adjust the scaling factor as needed
matching_indices <- as.numeric(rownames(mrIML_mat))

true_mat <- M[matching_indices, matching_indices]
rownames(true_mat) <- rownames(mrIML_mat)
colnames(true_mat) <- rownames(mrIML_mat)
g_true <- graph_from_adjacency_matrix(
  true_mat,
  mode = 'undirected',
  weighted = T,
  diag = F
)

Warning: The `adjmatrix` argument of `graph_from_adjacency_matrix()` must be symmetric
with mode = "undirected" as of igraph 1.6.0.
ℹ Use mode = "max" to achieve the original behavior.

edge_colors_true <- ifelse(E(g_true)$weight < 0, yes = "red", no = "blue")
edge_widths_true <- abs(E(g_true)$weight) * 5 # Adjust the scaling factor as needed

layout1 <- layout_nicely(g_true) # This produces a warning but not a big deal

Warning: Non-positive edge weight found, ignoring all weights during graph
layout.

Fit HMSC model and extract residual species associations.

studyDesign <- data.frame(sample = factor(1:nrow(Y)))
rL <- HmscRandomLevel(units = studyDesign$sample)
m_hmsc <- Hmsc(
  Y           = as.matrix(Y),
  XFormula    = ~1,
  XData       = data.frame(intercept = rep(1, nrow(Y))),
  distr       = "probit",
  studyDesign = studyDesign,
  ranLevels   = list("sample" = rL)
)
m_hmsc <- sampleMcmc(
  m_hmsc,
  thin      = 50,
  samples   = 1000,
  transient = 5000,
  nChains   = 4,
  verbose   = 0
)

mpost <- convertToCodaObject(m_hmsc)
psrf <- gelman.diag(mpost$Omega[[1]], multivariate = FALSE)$psrf
psrf

                               Point est. Upper C.I.
Omega1[sp1 (S1), sp1 (S1)]       1.002579   1.005567
Omega1[sp2 (S2), sp1 (S1)]       1.003338   1.011642
Omega1[sp3 (S3), sp1 (S1)]       1.007633   1.025017
Omega1[sp4 (S4), sp1 (S1)]       1.111253   1.322265
Omega1[sp7 (S5), sp1 (S1)]       1.022965   1.026343
Omega1[sp8 (S6), sp1 (S1)]       1.000758   1.004006
Omega1[sp9 (S7), sp1 (S1)]       1.072644   1.098432
Omega1[sp10 (S8), sp1 (S1)]      1.000802   1.004249
Omega1[sp11 (S9), sp1 (S1)]      1.000078   1.001790
Omega1[sp12 (S10), sp1 (S1)]     1.001323   1.003427
Omega1[sp14 (S11), sp1 (S1)]     1.004115   1.011201
Omega1[sp16 (S12), sp1 (S1)]     1.002281   1.006504
Omega1[sp18 (S13), sp1 (S1)]     1.019635   1.055160
Omega1[sp19 (S14), sp1 (S1)]     1.001909   1.003836
Omega1[sp1 (S1), sp2 (S2)]       1.003338   1.011642
Omega1[sp2 (S2), sp2 (S2)]       1.016959   1.021997
Omega1[sp3 (S3), sp2 (S2)]       1.000540   1.000757
Omega1[sp4 (S4), sp2 (S2)]       1.042006   1.096878
Omega1[sp7 (S5), sp2 (S2)]       1.101174   1.125807
Omega1[sp8 (S6), sp2 (S2)]       1.001505   1.004524
Omega1[sp9 (S7), sp2 (S2)]       1.127384   1.172535
Omega1[sp10 (S8), sp2 (S2)]      1.003408   1.008399
Omega1[sp11 (S9), sp2 (S2)]      1.000579   1.001177
Omega1[sp12 (S10), sp2 (S2)]     1.003208   1.011223
Omega1[sp14 (S11), sp2 (S2)]     1.003227   1.004472
Omega1[sp16 (S12), sp2 (S2)]     1.007575   1.020445
Omega1[sp18 (S13), sp2 (S2)]     1.002368   1.006750
Omega1[sp19 (S14), sp2 (S2)]     1.001879   1.004671
Omega1[sp1 (S1), sp3 (S3)]       1.007633   1.025017
Omega1[sp2 (S2), sp3 (S3)]       1.000540   1.000757
Omega1[sp3 (S3), sp3 (S3)]       1.004420   1.012804
Omega1[sp4 (S4), sp3 (S3)]       1.019943   1.024236
Omega1[sp7 (S5), sp3 (S3)]       1.132974   1.260815
Omega1[sp8 (S6), sp3 (S3)]       1.004705   1.014012
Omega1[sp9 (S7), sp3 (S3)]       1.162444   1.363484
Omega1[sp10 (S8), sp3 (S3)]      1.003650   1.007105
Omega1[sp11 (S9), sp3 (S3)]      1.012262   1.028198
Omega1[sp12 (S10), sp3 (S3)]     1.006543   1.021140
Omega1[sp14 (S11), sp3 (S3)]     1.005164   1.014684
Omega1[sp16 (S12), sp3 (S3)]     1.004465   1.005216
Omega1[sp18 (S13), sp3 (S3)]     1.014271   1.041018
Omega1[sp19 (S14), sp3 (S3)]     1.010009   1.023095
Omega1[sp1 (S1), sp4 (S4)]       1.111253   1.322265
Omega1[sp2 (S2), sp4 (S4)]       1.042006   1.096878
Omega1[sp3 (S3), sp4 (S4)]       1.019943   1.024236
Omega1[sp4 (S4), sp4 (S4)]       1.137633   1.402078
Omega1[sp7 (S5), sp4 (S4)]       1.041194   1.044866
Omega1[sp8 (S6), sp4 (S4)]       1.026640   1.049034
Omega1[sp9 (S7), sp4 (S4)]       1.053854   1.057850
Omega1[sp10 (S8), sp4 (S4)]      1.031952   1.070159
Omega1[sp11 (S9), sp4 (S4)]      1.026858   1.054462
Omega1[sp12 (S10), sp4 (S4)]     1.090760   1.257927
Omega1[sp14 (S11), sp4 (S4)]     1.024250   1.046065
Omega1[sp16 (S12), sp4 (S4)]     1.040775   1.086127
Omega1[sp18 (S13), sp4 (S4)]     1.111860   1.323821
Omega1[sp19 (S14), sp4 (S4)]     1.021086   1.034652
Omega1[sp1 (S1), sp7 (S5)]       1.022965   1.026343
Omega1[sp2 (S2), sp7 (S5)]       1.101174   1.125807
Omega1[sp3 (S3), sp7 (S5)]       1.132974   1.260815
Omega1[sp4 (S4), sp7 (S5)]       1.041194   1.044866
Omega1[sp7 (S5), sp7 (S5)]       1.966827   5.189187
Omega1[sp8 (S6), sp7 (S5)]       1.190070   1.576896
Omega1[sp9 (S7), sp7 (S5)]       3.064993   7.452516
Omega1[sp10 (S8), sp7 (S5)]      1.121440   1.195891
Omega1[sp11 (S9), sp7 (S5)]      1.160646   1.414712
Omega1[sp12 (S10), sp7 (S5)]     1.025766   1.026427
Omega1[sp14 (S11), sp7 (S5)]     1.198901   1.378766
Omega1[sp16 (S12), sp7 (S5)]     1.160138   1.405056
Omega1[sp18 (S13), sp7 (S5)]     1.006218   1.008522
Omega1[sp19 (S14), sp7 (S5)]     1.248532   1.870695
Omega1[sp1 (S1), sp8 (S6)]       1.000758   1.004006
Omega1[sp2 (S2), sp8 (S6)]       1.001505   1.004524
Omega1[sp3 (S3), sp8 (S6)]       1.004705   1.014012
Omega1[sp4 (S4), sp8 (S6)]       1.026640   1.049034
Omega1[sp7 (S5), sp8 (S6)]       1.190070   1.576896
Omega1[sp8 (S6), sp8 (S6)]       1.021919   1.038877
Omega1[sp9 (S7), sp8 (S6)]       1.203800   1.647567
Omega1[sp10 (S8), sp8 (S6)]      1.004540   1.011818
Omega1[sp11 (S9), sp8 (S6)]      1.005415   1.014445
Omega1[sp12 (S10), sp8 (S6)]     1.000903   1.004430
Omega1[sp14 (S11), sp8 (S6)]     1.003423   1.010058
Omega1[sp16 (S12), sp8 (S6)]     1.000155   1.000243
Omega1[sp18 (S13), sp8 (S6)]     1.002682   1.008011
Omega1[sp19 (S14), sp8 (S6)]     1.006978   1.018010
Omega1[sp1 (S1), sp9 (S7)]       1.072644   1.098432
Omega1[sp2 (S2), sp9 (S7)]       1.127384   1.172535
Omega1[sp3 (S3), sp9 (S7)]       1.162444   1.363484
Omega1[sp4 (S4), sp9 (S7)]       1.053854   1.057850
Omega1[sp7 (S5), sp9 (S7)]       3.064993   7.452516
Omega1[sp8 (S6), sp9 (S7)]       1.203800   1.647567
Omega1[sp9 (S7), sp9 (S7)]       1.774716   4.065773
Omega1[sp10 (S8), sp9 (S7)]      1.164105   1.321146
Omega1[sp11 (S9), sp9 (S7)]      1.185213   1.520832
Omega1[sp12 (S10), sp9 (S7)]     1.068747   1.070710
Omega1[sp14 (S11), sp9 (S7)]     1.181455   1.336070
Omega1[sp16 (S12), sp9 (S7)]     1.189615   1.527176
Omega1[sp18 (S13), sp9 (S7)]     1.033633   1.036806
Omega1[sp19 (S14), sp9 (S7)]     1.257982   1.918974
Omega1[sp1 (S1), sp10 (S8)]      1.000802   1.004249
Omega1[sp2 (S2), sp10 (S8)]      1.003408   1.008399
Omega1[sp3 (S3), sp10 (S8)]      1.003650   1.007105
Omega1[sp4 (S4), sp10 (S8)]      1.031952   1.070159
Omega1[sp7 (S5), sp10 (S8)]      1.121440   1.195891
Omega1[sp8 (S6), sp10 (S8)]      1.004540   1.011818
Omega1[sp9 (S7), sp10 (S8)]      1.164105   1.321146
Omega1[sp10 (S8), sp10 (S8)]     1.005202   1.007473
Omega1[sp11 (S9), sp10 (S8)]     1.001362   1.004189
Omega1[sp12 (S10), sp10 (S8)]    1.000376   1.001105
Omega1[sp14 (S11), sp10 (S8)]    1.006741   1.011182
Omega1[sp16 (S12), sp10 (S8)]    1.001469   1.005031
Omega1[sp18 (S13), sp10 (S8)]    1.001620   1.003381
Omega1[sp19 (S14), sp10 (S8)]    1.002564   1.004368
Omega1[sp1 (S1), sp11 (S9)]      1.000078   1.001790
Omega1[sp2 (S2), sp11 (S9)]      1.000579   1.001177
Omega1[sp3 (S3), sp11 (S9)]      1.012262   1.028198
Omega1[sp4 (S4), sp11 (S9)]      1.026858   1.054462
Omega1[sp7 (S5), sp11 (S9)]      1.160646   1.414712
Omega1[sp8 (S6), sp11 (S9)]      1.005415   1.014445
Omega1[sp9 (S7), sp11 (S9)]      1.185213   1.520832
Omega1[sp10 (S8), sp11 (S9)]     1.001362   1.004189
Omega1[sp11 (S9), sp11 (S9)]     1.000765   1.001519
Omega1[sp12 (S10), sp11 (S9)]    1.001025   1.004123
Omega1[sp14 (S11), sp11 (S9)]    1.001387   1.004330
Omega1[sp16 (S12), sp11 (S9)]    1.002140   1.004723
Omega1[sp18 (S13), sp11 (S9)]    1.005241   1.013241
Omega1[sp19 (S14), sp11 (S9)]    1.004879   1.011304
Omega1[sp1 (S1), sp12 (S10)]     1.001323   1.003427
Omega1[sp2 (S2), sp12 (S10)]     1.003208   1.011223
Omega1[sp3 (S3), sp12 (S10)]     1.006543   1.021140
Omega1[sp4 (S4), sp12 (S10)]     1.090760   1.257927
Omega1[sp7 (S5), sp12 (S10)]     1.025766   1.026427
Omega1[sp8 (S6), sp12 (S10)]     1.000903   1.004430
Omega1[sp9 (S7), sp12 (S10)]     1.068747   1.070710
Omega1[sp10 (S8), sp12 (S10)]    1.000376   1.001105
Omega1[sp11 (S9), sp12 (S10)]    1.001025   1.004123
Omega1[sp12 (S10), sp12 (S10)]   1.003634   1.006544
Omega1[sp14 (S11), sp12 (S10)]   1.004044   1.012693
Omega1[sp16 (S12), sp12 (S10)]   1.002004   1.005730
Omega1[sp18 (S13), sp12 (S10)]   1.019006   1.044453
Omega1[sp19 (S14), sp12 (S10)]   1.000648   1.001188
Omega1[sp1 (S1), sp14 (S11)]     1.004115   1.011201
Omega1[sp2 (S2), sp14 (S11)]     1.003227   1.004472
Omega1[sp3 (S3), sp14 (S11)]     1.005164   1.014684
Omega1[sp4 (S4), sp14 (S11)]     1.024250   1.046065
Omega1[sp7 (S5), sp14 (S11)]     1.198901   1.378766
Omega1[sp8 (S6), sp14 (S11)]     1.003423   1.010058
Omega1[sp9 (S7), sp14 (S11)]     1.181455   1.336070
Omega1[sp10 (S8), sp14 (S11)]    1.006741   1.011182
Omega1[sp11 (S9), sp14 (S11)]    1.001387   1.004330
Omega1[sp12 (S10), sp14 (S11)]   1.004044   1.012693
Omega1[sp14 (S11), sp14 (S11)]   1.035342   1.037011
Omega1[sp16 (S12), sp14 (S11)]   1.000216   1.001936
Omega1[sp18 (S13), sp14 (S11)]   1.009780   1.021970
Omega1[sp19 (S14), sp14 (S11)]   1.002029   1.005521
Omega1[sp1 (S1), sp16 (S12)]     1.002281   1.006504
Omega1[sp2 (S2), sp16 (S12)]     1.007575   1.020445
Omega1[sp3 (S3), sp16 (S12)]     1.004465   1.005216
Omega1[sp4 (S4), sp16 (S12)]     1.040775   1.086127
Omega1[sp7 (S5), sp16 (S12)]     1.160138   1.405056
Omega1[sp8 (S6), sp16 (S12)]     1.000155   1.000243
Omega1[sp9 (S7), sp16 (S12)]     1.189615   1.527176
Omega1[sp10 (S8), sp16 (S12)]    1.001469   1.005031
Omega1[sp11 (S9), sp16 (S12)]    1.002140   1.004723
Omega1[sp12 (S10), sp16 (S12)]   1.002004   1.005730
Omega1[sp14 (S11), sp16 (S12)]   1.000216   1.001936
Omega1[sp16 (S12), sp16 (S12)]   1.025021   1.031176
Omega1[sp18 (S13), sp16 (S12)]   1.004005   1.006651
Omega1[sp19 (S14), sp16 (S12)]   1.001040   1.002031
Omega1[sp1 (S1), sp18 (S13)]     1.019635   1.055160
Omega1[sp2 (S2), sp18 (S13)]     1.002368   1.006750
Omega1[sp3 (S3), sp18 (S13)]     1.014271   1.041018
Omega1[sp4 (S4), sp18 (S13)]     1.111860   1.323821
Omega1[sp7 (S5), sp18 (S13)]     1.006218   1.008522
Omega1[sp8 (S6), sp18 (S13)]     1.002682   1.008011
Omega1[sp9 (S7), sp18 (S13)]     1.033633   1.036806
Omega1[sp10 (S8), sp18 (S13)]    1.001620   1.003381
Omega1[sp11 (S9), sp18 (S13)]    1.005241   1.013241
Omega1[sp12 (S10), sp18 (S13)]   1.019006   1.044453
Omega1[sp14 (S11), sp18 (S13)]   1.009780   1.021970
Omega1[sp16 (S12), sp18 (S13)]   1.004005   1.006651
Omega1[sp18 (S13), sp18 (S13)]   1.052033   1.113861
Omega1[sp19 (S14), sp18 (S13)]   1.002157   1.003967
Omega1[sp1 (S1), sp19 (S14)]     1.001909   1.003836
Omega1[sp2 (S2), sp19 (S14)]     1.001879   1.004671
Omega1[sp3 (S3), sp19 (S14)]     1.010009   1.023095
Omega1[sp4 (S4), sp19 (S14)]     1.021086   1.034652
Omega1[sp7 (S5), sp19 (S14)]     1.248532   1.870695
Omega1[sp8 (S6), sp19 (S14)]     1.006978   1.018010
Omega1[sp9 (S7), sp19 (S14)]     1.257982   1.918974
Omega1[sp10 (S8), sp19 (S14)]    1.002564   1.004368
Omega1[sp11 (S9), sp19 (S14)]    1.004879   1.011304
Omega1[sp12 (S10), sp19 (S14)]   1.000648   1.001188
Omega1[sp14 (S11), sp19 (S14)]   1.002029   1.005521
Omega1[sp16 (S12), sp19 (S14)]   1.001040   1.002031
Omega1[sp18 (S13), sp19 (S14)]   1.002157   1.003967
Omega1[sp19 (S14), sp19 (S14)]   1.010227   1.018360

if (any(psrf[, "Point est."] > 1.1)) {
  cat("::: {.callout-warning}\n")
  cat("**HMSC convergence concern:** One or more Omega parameters have Rhat > 1.1,")
  cat(" indicating incomplete chain mixing. Interpret these results cautiously.\n")
  cat(":::\n")
}

Warning

HMSC convergence concern: One or more Omega parameters have Rhat > 1.1, indicating incomplete chain mixing. Interpret these results cautiously.

assoc <- computeAssociations(m_hmsc)[[1]]
OmegaCor <- assoc$mean
OmegaCor[assoc$support > 0.05 & assoc$support < 0.95] <- 0
hmsc_species_names <- paste0("sp", rownames(mrIML_mat))
hmsc_mat <- OmegaCor[hmsc_species_names, hmsc_species_names]
rownames(hmsc_mat) <- rownames(mrIML_mat)
colnames(hmsc_mat) <- rownames(mrIML_mat)

g_hmsc <- graph_from_adjacency_matrix(
  hmsc_mat,
  mode     = "undirected",
  weighted = TRUE,
  diag     = FALSE
)
edge_colors_hmsc <- ifelse(E(g_hmsc)$weight < 0, yes = "red", no = "blue")
edge_widths_hmsc <- abs(E(g_hmsc)$weight) * 10

par_old <- par(mfrow = c(1, 3))

plot(
  g_mrIML,
  layout = layout1,
  edge.color = edge_colors_mrIML,
  edge.width = edge_widths_mrIML,
  isolates = TRUE,
  main = "mrIML co-occurrence"
)

plot(
  g_hmsc,
  layout = layout1,
  edge.color = edge_colors_hmsc,
  edge.width = edge_widths_hmsc,
  isolates = TRUE,
  main = "HMSC co-occurrence"
)

plot(
  g_true,
  layout = layout1,
  edge.color = edge_colors_true,
  isolates = TRUE,
  edge.width = edge_widths_true,
  main = "Simulated (truth)"
)

par(par_old)

Compare the learned networks with the truth.

dist_mat_true <- 1 - true_mat
dist_mat_mrIML <- 1 - mrIML_mat
mantel_mrIML <- vegan::mantel(
  as.dist(dist_mat_true),
  as.dist(dist_mat_mrIML),
  method = "spearman"
)
mantel_mrIML

Mantel statistic based on Spearman's rank correlation rho 

Call:
vegan::mantel(xdis = as.dist(dist_mat_true), ydis = as.dist(dist_mat_mrIML),      method = "spearman") 

Mantel statistic r: 0.2781 
      Significance: 0.003 

Upper quantiles of permutations (null model):
  90%   95% 97.5%   99% 
0.124 0.170 0.183 0.213 
Permutation: free
Number of permutations: 999

dist_mat_hmsc <- 1 - hmsc_mat
mantel_hmsc <- vegan::mantel(
  as.dist(dist_mat_true),
  as.dist(dist_mat_hmsc),
  method = "spearman"
)
mantel_hmsc

Mantel statistic based on Spearman's rank correlation rho 

Call:
vegan::mantel(xdis = as.dist(dist_mat_true), ydis = as.dist(dist_mat_hmsc),      method = "spearman") 

Mantel statistic r: 0.1158 
      Significance: 0.089 

Upper quantiles of permutations (null model):
  90%   95% 97.5%   99% 
0.105 0.179 0.188 0.243 
Permutation: free
Number of permutations: 999

5.1 Mutualistic

simulated_mutual_network <- generateM_specific_type(
  nn = network_size,
  k_ave = k_average,
  type.network = "random",
  type.interact = "mutual",
  interact.str.max = 0.3,
  mix.compt.ratio = 0.5
)
true_mutual_network <- simulated_mutual_network[[1]]
M_mutual <- simulated_mutual_network[[2]]

data_mutual <- generateDataSet(
  900,
  M_mutual,
  count = network_size * 10000,
  mode = 4
) %>%
  t() %>%
  as.data.frame() %>%
  mutate(
    across(everything(), ~ ifelse(. < 1, yes = 0, no = 1))
  )

# Prepare data
Y_mutual <- filterRareCommon(data_mutual, lower = 0.01, higher = 0.99)
X1_mutual <- Y_mutual
# Fit mrIML obj
yhats_rf_sim_mutual <- mrIMLpredicts(
  Y = Y_mutual,
  X = NULL,
  X1 = X1_mutual,
  Model = model_rf,
  prop = 0.7,
  k = 5,
  racing = FALSE
)

i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"

→ A | warning: No event observations were detected in `truth` with event level '0'.

There were issues with some computations   A: x10
There were issues with some computations   A: x10

i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"

→ A | warning: No event observations were detected in `truth` with event level '0'.

There were issues with some computations   A: x10
There were issues with some computations   A: x10

i Creating pre-processing data to finalize 1 unknown parameter: "mtry"

# Bootstrapping
bs_sim_mutual <- mrBootstrap(yhats_rf_sim_mutual)
# Extract network
assoc_net_mutual <- mrCoOccurNet(bs_sim_mutual)

mrIML_mat_mutual <- assoc_net_mutual %>%
  filter(mean_strength > 0.1) %>%
  mutate(
    mean_strength_dir = ifelse(
      direction == "negative",
      yes = -mean_strength,
      no = mean_strength
    ),
    across(contains("taxa"), ~ sub("^sp", "", .))
  ) %>%
  graph_from_data_frame(
    directed = FALSE,
    vertices = sub("^sp", "", names(Y_mutual))
  ) %>%
  get.adjacency(attr = "mean_strength_dir", sparse = FALSE)

#plots
g_mrIML_mutual <- graph_from_adjacency_matrix(
  mrIML_mat_mutual,
  mode = 'undirected',
  weighted = T,
  diag = F
)

edge_colors_mrIML_mutual <- ifelse(
  E(g_mrIML_mutual)$weight < 0,
  yes = "red",
  no = "blue"
)
edge_widths_mrIML_mutual <- abs(E(g_mrIML_mutual)$weight) * 10 # Adjust the scaling factor as needed
matching_indices_mutual <- as.numeric(rownames(mrIML_mat_mutual))

true_mat_mutual <- M_mutual[matching_indices_mutual, matching_indices_mutual]
rownames(true_mat_mutual) <- rownames(mrIML_mat_mutual)
colnames(true_mat_mutual) <- rownames(mrIML_mat_mutual)
g_true_mutual <- graph_from_adjacency_matrix(
  true_mat_mutual,
  mode = 'undirected',
  weighted = T,
  diag = F
)
edge_colors_true_mutual <- ifelse(
  E(g_true_mutual)$weight < 0,
  yes = "red",
  no = "blue"
)
edge_widths_true_mutual <- abs(E(g_true_mutual)$weight) * 5 # Adjust the scaling factor as needed

layout1_mutual <- layout_nicely(g_true_mutual)

studyDesign_mutual <- data.frame(sample = factor(1:nrow(Y_mutual)))
rL_mutual <- HmscRandomLevel(units = studyDesign_mutual$sample)
m_hmsc_mutual <- Hmsc(
  Y           = as.matrix(Y_mutual),
  XFormula    = ~1,
  XData       = data.frame(intercept = rep(1, nrow(Y_mutual))),
  distr       = "probit",
  studyDesign = studyDesign_mutual,
  ranLevels   = list("sample" = rL_mutual)
)
m_hmsc_mutual <- sampleMcmc(
  m_hmsc_mutual,
  thin      = 50,
  samples   = 1000,
  transient = 5000,
  nChains   = 4,
  verbose   = 0
)

mpost_mutual <- convertToCodaObject(m_hmsc_mutual)
psrf_mutual <- gelman.diag(mpost_mutual$Omega[[1]], multivariate = FALSE)$psrf
psrf_mutual

                               Point est. Upper C.I.
Omega1[sp1 (S1), sp1 (S1)]      1.0011363   1.001994
Omega1[sp3 (S2), sp1 (S1)]      1.0181461   1.041691
Omega1[sp4 (S3), sp1 (S1)]      1.0009263   1.003813
Omega1[sp5 (S4), sp1 (S1)]      1.0005589   1.002397
Omega1[sp6 (S5), sp1 (S1)]      1.0044104   1.008312
Omega1[sp8 (S6), sp1 (S1)]      1.0010156   1.003846
Omega1[sp10 (S7), sp1 (S1)]     1.0000184   1.000972
Omega1[sp11 (S8), sp1 (S1)]     1.0012724   1.002793
Omega1[sp12 (S9), sp1 (S1)]     1.0032118   1.011042
Omega1[sp13 (S10), sp1 (S1)]    1.0041198   1.011264
Omega1[sp19 (S11), sp1 (S1)]    1.0010673   1.004851
Omega1[sp1 (S1), sp3 (S2)]      1.0181461   1.041691
Omega1[sp3 (S2), sp3 (S2)]      1.0749143   1.136986
Omega1[sp4 (S3), sp3 (S2)]      1.0035913   1.008305
Omega1[sp5 (S4), sp3 (S2)]      1.0140460   1.032510
Omega1[sp6 (S5), sp3 (S2)]      1.0098730   1.021577
Omega1[sp8 (S6), sp3 (S2)]      1.0120561   1.026569
Omega1[sp10 (S7), sp3 (S2)]     1.0100769   1.025636
Omega1[sp11 (S8), sp3 (S2)]     1.0055805   1.010506
Omega1[sp12 (S9), sp3 (S2)]     1.0064728   1.013953
Omega1[sp13 (S10), sp3 (S2)]    1.0155625   1.038022
Omega1[sp19 (S11), sp3 (S2)]    1.0055355   1.012416
Omega1[sp1 (S1), sp4 (S3)]      1.0009263   1.003813
Omega1[sp3 (S2), sp4 (S3)]      1.0035913   1.008305
Omega1[sp4 (S3), sp4 (S3)]      1.0018322   1.005321
Omega1[sp5 (S4), sp4 (S3)]      1.0024260   1.006268
Omega1[sp6 (S5), sp4 (S3)]      1.0024981   1.005516
Omega1[sp8 (S6), sp4 (S3)]      1.0008856   1.002416
Omega1[sp10 (S7), sp4 (S3)]     1.0006326   1.002709
Omega1[sp11 (S8), sp4 (S3)]     1.0026823   1.004594
Omega1[sp12 (S9), sp4 (S3)]     1.0016982   1.003932
Omega1[sp13 (S10), sp4 (S3)]    1.0025411   1.006440
Omega1[sp19 (S11), sp4 (S3)]    1.0007479   1.003199
Omega1[sp1 (S1), sp5 (S4)]      1.0005589   1.002397
Omega1[sp3 (S2), sp5 (S4)]      1.0140460   1.032510
Omega1[sp4 (S3), sp5 (S4)]      1.0024260   1.006268
Omega1[sp5 (S4), sp5 (S4)]      1.0002730   1.000965
Omega1[sp6 (S5), sp5 (S4)]      1.0027950   1.004593
Omega1[sp8 (S6), sp5 (S4)]      1.0002256   1.001564
Omega1[sp10 (S7), sp5 (S4)]     1.0003927   1.000733
Omega1[sp11 (S8), sp5 (S4)]     1.0017665   1.004178
Omega1[sp12 (S9), sp5 (S4)]     1.0057859   1.012058
Omega1[sp13 (S10), sp5 (S4)]    1.0029100   1.010169
Omega1[sp19 (S11), sp5 (S4)]    1.0010933   1.004760
Omega1[sp1 (S1), sp6 (S5)]      1.0044104   1.008312
Omega1[sp3 (S2), sp6 (S5)]      1.0098730   1.021577
Omega1[sp4 (S3), sp6 (S5)]      1.0024981   1.005516
Omega1[sp5 (S4), sp6 (S5)]      1.0027950   1.004593
Omega1[sp6 (S5), sp6 (S5)]      1.0221378   1.035102
Omega1[sp8 (S6), sp6 (S5)]      1.0040114   1.008724
Omega1[sp10 (S7), sp6 (S5)]     1.0042488   1.008228
Omega1[sp11 (S8), sp6 (S5)]     1.0009598   1.004519
Omega1[sp12 (S9), sp6 (S5)]     1.0054569   1.014189
Omega1[sp13 (S10), sp6 (S5)]    1.0034441   1.007928
Omega1[sp19 (S11), sp6 (S5)]    1.0029505   1.007068
Omega1[sp1 (S1), sp8 (S6)]      1.0010156   1.003846
Omega1[sp3 (S2), sp8 (S6)]      1.0120561   1.026569
Omega1[sp4 (S3), sp8 (S6)]      1.0008856   1.002416
Omega1[sp5 (S4), sp8 (S6)]      1.0002256   1.001564
Omega1[sp6 (S5), sp8 (S6)]      1.0040114   1.008724
Omega1[sp8 (S6), sp8 (S6)]      1.0003373   1.001392
Omega1[sp10 (S7), sp8 (S6)]     0.9994481   1.000111
Omega1[sp11 (S8), sp8 (S6)]     1.0007958   1.002746
Omega1[sp12 (S9), sp8 (S6)]     1.0052133   1.014016
Omega1[sp13 (S10), sp8 (S6)]    1.0023884   1.008486
Omega1[sp19 (S11), sp8 (S6)]    1.0017878   1.006507
Omega1[sp1 (S1), sp10 (S7)]     1.0000184   1.000972
Omega1[sp3 (S2), sp10 (S7)]     1.0100769   1.025636
Omega1[sp4 (S3), sp10 (S7)]     1.0006326   1.002709
Omega1[sp5 (S4), sp10 (S7)]     1.0003927   1.000733
Omega1[sp6 (S5), sp10 (S7)]     1.0042488   1.008228
Omega1[sp8 (S6), sp10 (S7)]     0.9994481   1.000111
Omega1[sp10 (S7), sp10 (S7)]    1.0004874   1.002000
Omega1[sp11 (S8), sp10 (S7)]    1.0012993   1.001673
Omega1[sp12 (S9), sp10 (S7)]    1.0010482   1.003524
Omega1[sp13 (S10), sp10 (S7)]   1.0019386   1.006017
Omega1[sp19 (S11), sp10 (S7)]   1.0008398   1.003451
Omega1[sp1 (S1), sp11 (S8)]     1.0012724   1.002793
Omega1[sp3 (S2), sp11 (S8)]     1.0055805   1.010506
Omega1[sp4 (S3), sp11 (S8)]     1.0026823   1.004594
Omega1[sp5 (S4), sp11 (S8)]     1.0017665   1.004178
Omega1[sp6 (S5), sp11 (S8)]     1.0009598   1.004519
Omega1[sp8 (S6), sp11 (S8)]     1.0007958   1.002746
Omega1[sp10 (S7), sp11 (S8)]    1.0012993   1.001673
Omega1[sp11 (S8), sp11 (S8)]    1.0018240   1.003689
Omega1[sp12 (S9), sp11 (S8)]    1.0039650   1.009737
Omega1[sp13 (S10), sp11 (S8)]   1.0006810   1.003274
Omega1[sp19 (S11), sp11 (S8)]   1.0036261   1.007786
Omega1[sp1 (S1), sp12 (S9)]     1.0032118   1.011042
Omega1[sp3 (S2), sp12 (S9)]     1.0064728   1.013953
Omega1[sp4 (S3), sp12 (S9)]     1.0016982   1.003932
Omega1[sp5 (S4), sp12 (S9)]     1.0057859   1.012058
Omega1[sp6 (S5), sp12 (S9)]     1.0054569   1.014189
Omega1[sp8 (S6), sp12 (S9)]     1.0052133   1.014016
Omega1[sp10 (S7), sp12 (S9)]    1.0010482   1.003524
Omega1[sp11 (S8), sp12 (S9)]    1.0039650   1.009737
Omega1[sp12 (S9), sp12 (S9)]    1.0080927   1.018717
Omega1[sp13 (S10), sp12 (S9)]   1.0032941   1.011700
Omega1[sp19 (S11), sp12 (S9)]   1.0041497   1.012408
Omega1[sp1 (S1), sp13 (S10)]    1.0041198   1.011264
Omega1[sp3 (S2), sp13 (S10)]    1.0155625   1.038022
Omega1[sp4 (S3), sp13 (S10)]    1.0025411   1.006440
Omega1[sp5 (S4), sp13 (S10)]    1.0029100   1.010169
Omega1[sp6 (S5), sp13 (S10)]    1.0034441   1.007928
Omega1[sp8 (S6), sp13 (S10)]    1.0023884   1.008486
Omega1[sp10 (S7), sp13 (S10)]   1.0019386   1.006017
Omega1[sp11 (S8), sp13 (S10)]   1.0006810   1.003274
Omega1[sp12 (S9), sp13 (S10)]   1.0032941   1.011700
Omega1[sp13 (S10), sp13 (S10)]  1.0043786   1.010953
Omega1[sp19 (S11), sp13 (S10)]  1.0049489   1.015112
Omega1[sp1 (S1), sp19 (S11)]    1.0010673   1.004851
Omega1[sp3 (S2), sp19 (S11)]    1.0055355   1.012416
Omega1[sp4 (S3), sp19 (S11)]    1.0007479   1.003199
Omega1[sp5 (S4), sp19 (S11)]    1.0010933   1.004760
Omega1[sp6 (S5), sp19 (S11)]    1.0029505   1.007068
Omega1[sp8 (S6), sp19 (S11)]    1.0017878   1.006507
Omega1[sp10 (S7), sp19 (S11)]   1.0008398   1.003451
Omega1[sp11 (S8), sp19 (S11)]   1.0036261   1.007786
Omega1[sp12 (S9), sp19 (S11)]   1.0041497   1.012408
Omega1[sp13 (S10), sp19 (S11)]  1.0049489   1.015112
Omega1[sp19 (S11), sp19 (S11)]  1.0046386   1.010554

if (any(psrf_mutual[, "Point est."] > 1.1)) {
  cat("::: {.callout-warning}\n")
  cat("**HMSC convergence concern:** One or more Omega parameters have Rhat > 1.1,")
  cat(" indicating incomplete chain mixing. Interpret these results cautiously.\n")
  cat(":::\n")
}

assoc_mutual <- computeAssociations(m_hmsc_mutual)[[1]]
OmegaCor_mutual <- assoc_mutual$mean
OmegaCor_mutual[assoc_mutual$support > 0.05 & assoc_mutual$support < 0.95] <- 0
hmsc_species_names_mutual <- paste0("sp", rownames(mrIML_mat_mutual))
hmsc_mat_mutual <- OmegaCor_mutual[hmsc_species_names_mutual, hmsc_species_names_mutual]
rownames(hmsc_mat_mutual) <- rownames(mrIML_mat_mutual)
colnames(hmsc_mat_mutual) <- rownames(mrIML_mat_mutual)

g_hmsc_mutual <- graph_from_adjacency_matrix(
  hmsc_mat_mutual,
  mode     = "undirected",
  weighted = TRUE,
  diag     = FALSE
)
edge_colors_hmsc_mutual <- ifelse(E(g_hmsc_mutual)$weight < 0, yes = "red", no = "blue")
edge_widths_hmsc_mutual <- abs(E(g_hmsc_mutual)$weight) * 10

par_old <- par(mfrow = c(1, 3))

plot(
  g_mrIML_mutual,
  layout = layout1_mutual,
  edge.color = edge_colors_mrIML_mutual,
  edge.width = edge_widths_mrIML_mutual,
  isolates = TRUE,
  main = "mrIML co-occurrence"
)

plot(
  g_hmsc_mutual,
  layout = layout1_mutual,
  edge.color = edge_colors_hmsc_mutual,
  edge.width = edge_widths_hmsc_mutual,
  isolates = TRUE,
  main = "HMSC co-occurrence"
)

plot(
  g_true_mutual,
  layout = layout1_mutual,
  edge.color = edge_colors_true_mutual,
  isolates = TRUE,
  edge.width = edge_widths_true_mutual,
  main = "Simulated (truth)"
)

par(par_old)

dist_mat_true_mutual <- 1 - true_mat_mutual
dist_mat_mrIML_mutual <- 1 - mrIML_mat_mutual
mantel_mrIML_mutual <- vegan::mantel(
  as.dist(dist_mat_true_mutual),
  as.dist(dist_mat_mrIML_mutual),
  method = "spearman"
)
mantel_mrIML_mutual

Mantel statistic based on Spearman's rank correlation rho 

Call:
vegan::mantel(xdis = as.dist(dist_mat_true_mutual), ydis = as.dist(dist_mat_mrIML_mutual),      method = "spearman") 

Mantel statistic r: 0.08057 
      Significance: 0.281 

Upper quantiles of permutations (null model):
  90%   95% 97.5%   99% 
0.194 0.242 0.297 0.357 
Permutation: free
Number of permutations: 999

dist_mat_hmsc_mutual <- 1 - hmsc_mat_mutual
mantel_hmsc_mutual <- vegan::mantel(
  as.dist(dist_mat_true_mutual),
  as.dist(dist_mat_hmsc_mutual),
  method = "spearman"
)
mantel_hmsc_mutual

Mantel statistic based on Spearman's rank correlation rho 

Call:
vegan::mantel(xdis = as.dist(dist_mat_true_mutual), ydis = as.dist(dist_mat_hmsc_mutual),      method = "spearman") 

Mantel statistic r: 0.3157 
      Significance: 0.023 

Upper quantiles of permutations (null model):
  90%   95% 97.5%   99% 
0.202 0.261 0.303 0.363 
Permutation: free
Number of permutations: 999

5.2 Competitive

simulated_compt_network <- generateM_specific_type(
  nn = network_size,
  k_ave = k_average,
  type.network = "random",
  type.interact = "compt",
  interact.str.max = 0.5,
  mix.compt.ratio = 0.5
)
true_compt_network <- simulated_compt_network[[1]]
M_compt <- simulated_compt_network[[2]]
diag(M_compt) <- -0.5 # Reduce the level of self regulation

data_compt <- generateDataSet(
  900,
  M_compt,
  count = network_size * 10000,
  mode = 4
) %>%
  t() %>%
  as.data.frame() %>%
  mutate(
    across(everything(), ~ ifelse(. < 1, yes = 0, no = 1))
  )

# Prepare data
Y_compt <- filterRareCommon(data_compt, lower = 0.01, higher = 0.99)
X1_compt <- Y_compt
# Fit mrIML obj
yhats_rf_sim_compt <- mrIMLpredicts(
  Y = Y_compt,
  X = NULL,
  X1 = X1_compt,
  Model = model_rf,
  prop = 0.7,
  k = 5,
  racing = FALSE
)

i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"

# Bootstrapping
bs_sim_compt <- mrBootstrap(yhats_rf_sim_compt)
# Extract network
assoc_net_compt <- mrCoOccurNet(bs_sim_compt)

mrIML_mat_compt <- assoc_net_compt %>%
  filter(mean_strength > 0.1) %>%
  mutate(
    mean_strength_dir = ifelse(
      direction == "negative",
      yes = -mean_strength,
      no = mean_strength
    ),
    across(contains("taxa"), ~ sub("^sp", "", .))
  ) %>%
  graph_from_data_frame(
    directed = FALSE,
    vertices = sub("^sp", "", names(Y_compt))
  ) %>%
  get.adjacency(attr = "mean_strength_dir", sparse = FALSE)

#plots
g_mrIML_compt <- graph_from_adjacency_matrix(
  mrIML_mat_compt,
  mode = 'undirected',
  weighted = T,
  diag = F
)

edge_colors_mrIML_compt <- ifelse(
  E(g_mrIML_compt)$weight < 0,
  yes = "red",
  no = "blue"
)
edge_widths_mrIML_compt <- abs(E(g_mrIML_compt)$weight) * 10 # Adjust the scaling factor as needed
matching_indices_compt <- as.numeric(rownames(mrIML_mat_compt))

true_mat_compt <- M_compt[matching_indices_compt, matching_indices_compt]
rownames(true_mat_compt) <- rownames(mrIML_mat_compt)
colnames(true_mat_compt) <- rownames(mrIML_mat_compt)
g_true_compt <- graph_from_adjacency_matrix(
  true_mat_compt,
  mode = 'undirected',
  weighted = T,
  diag = F
)
edge_colors_true_compt <- ifelse(
  E(g_true_compt)$weight < 0,
  yes = "red",
  no = "blue"
)
edge_widths_true_compt <- abs(E(g_true_compt)$weight) * 5 # Adjust the scaling factor as needed

layout1_compt <- layout_nicely(g_true_compt)

Warning: Non-positive edge weight found, ignoring all weights during graph
layout.

studyDesign_compt <- data.frame(sample = factor(1:nrow(Y_compt)))
rL_compt <- HmscRandomLevel(units = studyDesign_compt$sample)
m_hmsc_compt <- Hmsc(
  Y           = as.matrix(Y_compt),
  XFormula    = ~1,
  XData       = data.frame(intercept = rep(1, nrow(Y_compt))),
  distr       = "probit",
  studyDesign = studyDesign_compt,
  ranLevels   = list("sample" = rL_compt)
)
m_hmsc_compt <- sampleMcmc(
  m_hmsc_compt,
  thin      = 50,
  samples   = 1000,
  transient = 5000,
  nChains   = 4,
  verbose   = 0
)

mpost_compt <- convertToCodaObject(m_hmsc_compt)
psrf_compt <- gelman.diag(mpost_compt$Omega[[1]], multivariate = FALSE)$psrf
psrf_compt

                               Point est. Upper C.I.
Omega1[sp1 (S1), sp1 (S1)]       1.270422   1.762505
Omega1[sp2 (S2), sp1 (S1)]       2.256754   3.953556
Omega1[sp3 (S3), sp1 (S1)]       2.356771   4.164917
Omega1[sp4 (S4), sp1 (S1)]       1.086477   1.172172
Omega1[sp5 (S5), sp1 (S1)]       2.079887   3.663132
Omega1[sp6 (S6), sp1 (S1)]       1.451965   2.195345
Omega1[sp7 (S7), sp1 (S1)]       1.012400   1.037824
Omega1[sp8 (S8), sp1 (S1)]       2.302071   4.125137
Omega1[sp9 (S9), sp1 (S1)]       1.409530   2.065599
Omega1[sp10 (S10), sp1 (S1)]     1.275245   1.295672
Omega1[sp11 (S11), sp1 (S1)]     1.415438   3.391775
Omega1[sp12 (S12), sp1 (S1)]     1.548997   2.287311
Omega1[sp13 (S13), sp1 (S1)]     1.076130   1.204723
Omega1[sp14 (S14), sp1 (S1)]     1.638921   2.510918
Omega1[sp15 (S15), sp1 (S1)]     1.284190   1.753065
Omega1[sp16 (S16), sp1 (S1)]     1.043077   1.121326
Omega1[sp17 (S17), sp1 (S1)]     1.094172   1.248602
Omega1[sp18 (S18), sp1 (S1)]     2.113410   3.441301
Omega1[sp19 (S19), sp1 (S1)]     2.064257   3.295152
Omega1[sp20 (S20), sp1 (S1)]     2.044315   3.589435
Omega1[sp1 (S1), sp2 (S2)]       2.256754   3.953556
Omega1[sp2 (S2), sp2 (S2)]       2.483049   4.433643
Omega1[sp3 (S3), sp2 (S2)]       3.782437   7.077664
Omega1[sp4 (S4), sp2 (S2)]       1.251716   1.714372
Omega1[sp5 (S5), sp2 (S2)]       3.763737   7.411932
Omega1[sp6 (S6), sp2 (S2)]       1.893039   3.175918
Omega1[sp7 (S7), sp2 (S2)]       2.617260   4.347640
Omega1[sp8 (S8), sp2 (S2)]       3.124237   5.896341
Omega1[sp9 (S9), sp2 (S2)]       1.307867   1.857076
Omega1[sp10 (S10), sp2 (S2)]     1.390661   2.816312
Omega1[sp11 (S11), sp2 (S2)]     1.146678   1.365786
Omega1[sp12 (S12), sp2 (S2)]     1.037929   1.076708
Omega1[sp13 (S13), sp2 (S2)]     1.666691   2.478975
Omega1[sp14 (S14), sp2 (S2)]     1.262411   1.720752
Omega1[sp15 (S15), sp2 (S2)]     1.036195   1.107134
Omega1[sp16 (S16), sp2 (S2)]     2.057829   3.322602
Omega1[sp17 (S17), sp2 (S2)]     1.755395   2.718741
Omega1[sp18 (S18), sp2 (S2)]     1.936293   3.157408
Omega1[sp19 (S19), sp2 (S2)]     1.598978   2.457425
Omega1[sp20 (S20), sp2 (S2)]     3.597732   6.921324
Omega1[sp1 (S1), sp3 (S3)]       2.356771   4.164917
Omega1[sp2 (S2), sp3 (S3)]       3.782437   7.077664
Omega1[sp3 (S3), sp3 (S3)]       2.601721   4.737706
Omega1[sp4 (S4), sp3 (S3)]       1.093398   1.269600
Omega1[sp5 (S5), sp3 (S3)]       3.166908   6.205153
Omega1[sp6 (S6), sp3 (S3)]       1.906478   3.234905
Omega1[sp7 (S7), sp3 (S3)]       2.112127   3.355252
Omega1[sp8 (S8), sp3 (S3)]       2.598851   4.836893
Omega1[sp9 (S9), sp3 (S3)]       1.353945   1.984707
Omega1[sp10 (S10), sp3 (S3)]     1.248999   1.250016
Omega1[sp11 (S11), sp3 (S3)]     1.212425   1.596020
Omega1[sp12 (S12), sp3 (S3)]     1.084074   1.216852
Omega1[sp13 (S13), sp3 (S3)]     1.310456   1.755893
Omega1[sp14 (S14), sp3 (S3)]     1.359958   1.981669
Omega1[sp15 (S15), sp3 (S3)]     1.028616   1.072121
Omega1[sp16 (S16), sp3 (S3)]     1.467756   2.108053
Omega1[sp17 (S17), sp3 (S3)]     1.299711   1.728927
Omega1[sp18 (S18), sp3 (S3)]     2.254668   3.829278
Omega1[sp19 (S19), sp3 (S3)]     1.843128   2.993526
Omega1[sp20 (S20), sp3 (S3)]     4.484489   8.623411
Omega1[sp1 (S1), sp4 (S4)]       1.086477   1.172172
Omega1[sp2 (S2), sp4 (S4)]       1.251716   1.714372
Omega1[sp3 (S3), sp4 (S4)]       1.093398   1.269600
Omega1[sp4 (S4), sp4 (S4)]       1.503676   3.878657
Omega1[sp5 (S5), sp4 (S4)]       1.050261   1.070709
Omega1[sp6 (S6), sp4 (S4)]       1.124692   1.310379
Omega1[sp7 (S7), sp4 (S4)]       1.701717   3.508329
Omega1[sp8 (S8), sp4 (S4)]       1.097803   1.288252
Omega1[sp9 (S9), sp4 (S4)]       1.039390   1.073012
Omega1[sp10 (S10), sp4 (S4)]     1.831792   3.576396
Omega1[sp11 (S11), sp4 (S4)]     1.273935   1.769959
Omega1[sp12 (S12), sp4 (S4)]     1.478085   2.726861
Omega1[sp13 (S13), sp4 (S4)]     1.481443   2.617217
Omega1[sp14 (S14), sp4 (S4)]     1.064958   1.173735
Omega1[sp15 (S15), sp4 (S4)]     1.194046   1.529358
Omega1[sp16 (S16), sp4 (S4)]     1.618682   3.222337
Omega1[sp17 (S17), sp4 (S4)]     1.529137   2.742858
Omega1[sp18 (S18), sp4 (S4)]     1.082042   1.217280
Omega1[sp19 (S19), sp4 (S4)]     1.172287   1.498882
Omega1[sp20 (S20), sp4 (S4)]     1.136485   1.403853
Omega1[sp1 (S1), sp5 (S5)]       2.079887   3.663132
Omega1[sp2 (S2), sp5 (S5)]       3.763737   7.411932
Omega1[sp3 (S3), sp5 (S5)]       3.166908   6.205153
Omega1[sp4 (S4), sp5 (S5)]       1.050261   1.070709
Omega1[sp5 (S5), sp5 (S5)]       2.161167   4.334407
Omega1[sp6 (S6), sp5 (S5)]       1.831710   3.133652
Omega1[sp7 (S7), sp5 (S5)]       1.732751   2.886528
Omega1[sp8 (S8), sp5 (S5)]       2.813188   5.366644
Omega1[sp9 (S9), sp5 (S5)]       1.426270   2.197006
Omega1[sp10 (S10), sp5 (S5)]     1.047692   1.140589
Omega1[sp11 (S11), sp5 (S5)]     1.037431   1.050601
Omega1[sp12 (S12), sp5 (S5)]     2.337627   4.145845
Omega1[sp13 (S13), sp5 (S5)]     1.051015   1.064944
Omega1[sp14 (S14), sp5 (S5)]     1.531714   2.459877
Omega1[sp15 (S15), sp5 (S5)]     1.347792   1.995293
Omega1[sp16 (S16), sp5 (S5)]     1.036694   1.038282
Omega1[sp17 (S17), sp5 (S5)]     1.060420   1.138031
Omega1[sp18 (S18), sp5 (S5)]     2.635934   4.860072
Omega1[sp19 (S19), sp5 (S5)]     2.177819   3.909196
Omega1[sp20 (S20), sp5 (S5)]     3.326483   6.505529
Omega1[sp1 (S1), sp6 (S6)]       1.451965   2.195345
Omega1[sp2 (S2), sp6 (S6)]       1.893039   3.175918
Omega1[sp3 (S3), sp6 (S6)]       1.906478   3.234905
Omega1[sp4 (S4), sp6 (S6)]       1.124692   1.310379
Omega1[sp5 (S5), sp6 (S6)]       1.831710   3.133652
Omega1[sp6 (S6), sp6 (S6)]       1.208992   1.619918
Omega1[sp7 (S7), sp6 (S6)]       1.101758   1.276177
Omega1[sp8 (S8), sp6 (S6)]       1.823828   3.125586
Omega1[sp9 (S9), sp6 (S6)]       1.301789   1.834714
Omega1[sp10 (S10), sp6 (S6)]     1.586783   6.133118
Omega1[sp11 (S11), sp6 (S6)]     1.409938   3.818275
Omega1[sp12 (S12), sp6 (S6)]     1.217145   1.565752
Omega1[sp13 (S13), sp6 (S6)]     1.016142   1.016768
Omega1[sp14 (S14), sp6 (S6)]     1.319697   1.818658
Omega1[sp15 (S15), sp6 (S6)]     1.528251   2.843589
Omega1[sp16 (S16), sp6 (S6)]     1.006542   1.007695
Omega1[sp17 (S17), sp6 (S6)]     1.019085   1.021026
Omega1[sp18 (S18), sp6 (S6)]     1.703217   2.600789
Omega1[sp19 (S19), sp6 (S6)]     1.685327   2.550552
Omega1[sp20 (S20), sp6 (S6)]     1.874537   3.205330
Omega1[sp1 (S1), sp7 (S7)]       1.012400   1.037824
Omega1[sp2 (S2), sp7 (S7)]       2.617260   4.347640
Omega1[sp3 (S3), sp7 (S7)]       2.112127   3.355252
Omega1[sp4 (S4), sp7 (S7)]       1.701717   3.508329
Omega1[sp5 (S5), sp7 (S7)]       1.732751   2.886528
Omega1[sp6 (S6), sp7 (S7)]       1.101758   1.276177
Omega1[sp7 (S7), sp7 (S7)]       1.322426   1.795485
Omega1[sp8 (S8), sp7 (S7)]       2.103686   3.589119
Omega1[sp9 (S9), sp7 (S7)]       1.236605   1.596334
Omega1[sp10 (S10), sp7 (S7)]     2.203257  10.512896
Omega1[sp11 (S11), sp7 (S7)]     3.475752  16.963435
Omega1[sp12 (S12), sp7 (S7)]     1.271909   1.663453
Omega1[sp13 (S13), sp7 (S7)]     1.000450   1.000853
Omega1[sp14 (S14), sp7 (S7)]     1.685951   2.585173
Omega1[sp15 (S15), sp7 (S7)]     1.343045   1.911503
Omega1[sp16 (S16), sp7 (S7)]     1.225268   1.568044
Omega1[sp17 (S17), sp7 (S7)]     1.037715   1.110113
Omega1[sp18 (S18), sp7 (S7)]     1.968825   3.072294
Omega1[sp19 (S19), sp7 (S7)]     2.276410   3.666242
Omega1[sp20 (S20), sp7 (S7)]     2.076216   3.445866
Omega1[sp1 (S1), sp8 (S8)]       2.302071   4.125137
Omega1[sp2 (S2), sp8 (S8)]       3.124237   5.896341
Omega1[sp3 (S3), sp8 (S8)]       2.598851   4.836893
Omega1[sp4 (S4), sp8 (S8)]       1.097803   1.288252
Omega1[sp5 (S5), sp8 (S8)]       2.813188   5.366644
Omega1[sp6 (S6), sp8 (S8)]       1.823828   3.125586
Omega1[sp7 (S7), sp8 (S8)]       2.103686   3.589119
Omega1[sp8 (S8), sp8 (S8)]       1.790947   3.116747
Omega1[sp9 (S9), sp8 (S8)]       1.391020   2.109206
Omega1[sp10 (S10), sp8 (S8)]     1.105432   1.227558
Omega1[sp11 (S11), sp8 (S8)]     1.042637   1.126218
Omega1[sp12 (S12), sp8 (S8)]     1.663671   2.577455
Omega1[sp13 (S13), sp8 (S8)]     1.108450   1.322537
Omega1[sp14 (S14), sp8 (S8)]     1.463845   2.272318
Omega1[sp15 (S15), sp8 (S8)]     1.207462   1.598785
Omega1[sp16 (S16), sp8 (S8)]     1.197246   1.525229
Omega1[sp17 (S17), sp8 (S8)]     1.047342   1.128882
Omega1[sp18 (S18), sp8 (S8)]     2.293735   4.094142
Omega1[sp19 (S19), sp8 (S8)]     2.115426   3.705113
Omega1[sp20 (S20), sp8 (S8)]     4.726713   9.203475
Omega1[sp1 (S1), sp9 (S9)]       1.409530   2.065599
Omega1[sp2 (S2), sp9 (S9)]       1.307867   1.857076
Omega1[sp3 (S3), sp9 (S9)]       1.353945   1.984707
Omega1[sp4 (S4), sp9 (S9)]       1.039390   1.073012
Omega1[sp5 (S5), sp9 (S9)]       1.426270   2.197006
Omega1[sp6 (S6), sp9 (S9)]       1.301789   1.834714
Omega1[sp7 (S7), sp9 (S9)]       1.236605   1.596334
Omega1[sp8 (S8), sp9 (S9)]       1.391020   2.109206
Omega1[sp9 (S9), sp9 (S9)]       1.066910   1.188546
Omega1[sp10 (S10), sp9 (S9)]     1.359199   3.285121
Omega1[sp11 (S11), sp9 (S9)]     1.443967   4.430039
Omega1[sp12 (S12), sp9 (S9)]     1.106646   1.172803
Omega1[sp13 (S13), sp9 (S9)]     1.230998   1.641228
Omega1[sp14 (S14), sp9 (S9)]     1.039665   1.118837
Omega1[sp15 (S15), sp9 (S9)]     1.214715   1.458027
Omega1[sp16 (S16), sp9 (S9)]     1.157480   1.400326
Omega1[sp17 (S17), sp9 (S9)]     1.178468   1.465912
Omega1[sp18 (S18), sp9 (S9)]     1.113840   1.316262
Omega1[sp19 (S19), sp9 (S9)]     1.075776   1.213683
Omega1[sp20 (S20), sp9 (S9)]     1.398873   2.126414
Omega1[sp1 (S1), sp10 (S10)]     1.275245   1.295672
Omega1[sp2 (S2), sp10 (S10)]     1.390661   2.816312
Omega1[sp3 (S3), sp10 (S10)]     1.248999   1.250016
Omega1[sp4 (S4), sp10 (S10)]     1.831792   3.576396
Omega1[sp5 (S5), sp10 (S10)]     1.047692   1.140589
Omega1[sp6 (S6), sp10 (S10)]     1.586783   6.133118
Omega1[sp7 (S7), sp10 (S10)]     2.203257  10.512896
Omega1[sp8 (S8), sp10 (S10)]     1.105432   1.227558
Omega1[sp9 (S9), sp10 (S10)]     1.359199   3.285121
Omega1[sp10 (S10), sp10 (S10)]   4.332217  26.045394
Omega1[sp11 (S11), sp10 (S10)]   5.924001  36.296716
Omega1[sp12 (S12), sp10 (S10)]   2.095047   9.013630
Omega1[sp13 (S13), sp10 (S10)]   1.349508   3.070748
Omega1[sp14 (S14), sp10 (S10)]   1.974568   9.687445
Omega1[sp15 (S15), sp10 (S10)]   6.249137  38.287672
Omega1[sp16 (S16), sp10 (S10)]   2.202570   8.011191
Omega1[sp17 (S17), sp10 (S10)]   1.268834   1.317881
Omega1[sp18 (S18), sp10 (S10)]   1.462826   4.627038
Omega1[sp19 (S19), sp10 (S10)]   1.586783   5.706480
Omega1[sp20 (S20), sp10 (S10)]   1.192098   1.398793
Omega1[sp1 (S1), sp11 (S11)]     1.415438   3.391775
Omega1[sp2 (S2), sp11 (S11)]     1.146678   1.365786
Omega1[sp3 (S3), sp11 (S11)]     1.212425   1.596020
Omega1[sp4 (S4), sp11 (S11)]     1.273935   1.769959
Omega1[sp5 (S5), sp11 (S11)]     1.037431   1.050601
Omega1[sp6 (S6), sp11 (S11)]     1.409938   3.818275
Omega1[sp7 (S7), sp11 (S11)]     3.475752  16.963435
Omega1[sp8 (S8), sp11 (S11)]     1.042637   1.126218
Omega1[sp9 (S9), sp11 (S11)]     1.443967   4.430039
Omega1[sp10 (S10), sp11 (S11)]   5.924001  36.296716
Omega1[sp11 (S11), sp11 (S11)]   3.827910  22.742716
Omega1[sp12 (S12), sp11 (S11)]   2.053221   7.622777
Omega1[sp13 (S13), sp11 (S11)]   2.463439  12.218637
Omega1[sp14 (S14), sp11 (S11)]   2.270148  11.565373
Omega1[sp15 (S15), sp11 (S11)]   7.937467  48.667801
Omega1[sp16 (S16), sp11 (S11)]   2.082228   6.695805
Omega1[sp17 (S17), sp11 (S11)]   2.564397  11.608945
Omega1[sp18 (S18), sp11 (S11)]   1.355140   2.891993
Omega1[sp19 (S19), sp11 (S11)]   1.418321   3.353283
Omega1[sp20 (S20), sp11 (S11)]   1.016420   1.022556
Omega1[sp1 (S1), sp12 (S12)]     1.548997   2.287311
Omega1[sp2 (S2), sp12 (S12)]     1.037929   1.076708
Omega1[sp3 (S3), sp12 (S12)]     1.084074   1.216852
Omega1[sp4 (S4), sp12 (S12)]     1.478085   2.726861
Omega1[sp5 (S5), sp12 (S12)]     2.337627   4.145845
Omega1[sp6 (S6), sp12 (S12)]     1.217145   1.565752
Omega1[sp7 (S7), sp12 (S12)]     1.271909   1.663453
Omega1[sp8 (S8), sp12 (S12)]     1.663671   2.577455
Omega1[sp9 (S9), sp12 (S12)]     1.106646   1.172803
Omega1[sp10 (S10), sp12 (S12)]   2.095047   9.013630
Omega1[sp11 (S11), sp12 (S12)]   2.053221   7.622777
Omega1[sp12 (S12), sp12 (S12)]   1.825208   3.773916
Omega1[sp13 (S13), sp12 (S12)]   1.877132   3.266352
Omega1[sp14 (S14), sp12 (S12)]   1.388664   2.250585
Omega1[sp15 (S15), sp12 (S12)]   2.043197   3.805658
Omega1[sp16 (S16), sp12 (S12)]   1.420419   2.005956
Omega1[sp17 (S17), sp12 (S12)]   1.755644   3.033060
Omega1[sp18 (S18), sp12 (S12)]   1.653418   2.704892
Omega1[sp19 (S19), sp12 (S12)]   1.859729   3.346476
Omega1[sp20 (S20), sp12 (S12)]   1.583294   2.333122
Omega1[sp1 (S1), sp13 (S13)]     1.076130   1.204723
Omega1[sp2 (S2), sp13 (S13)]     1.666691   2.478975
Omega1[sp3 (S3), sp13 (S13)]     1.310456   1.755893
Omega1[sp4 (S4), sp13 (S13)]     1.481443   2.617217
Omega1[sp5 (S5), sp13 (S13)]     1.051015   1.064944
Omega1[sp6 (S6), sp13 (S13)]     1.016142   1.016768
Omega1[sp7 (S7), sp13 (S13)]     1.000450   1.000853
Omega1[sp8 (S8), sp13 (S13)]     1.108450   1.322537
Omega1[sp9 (S9), sp13 (S13)]     1.230998   1.641228
Omega1[sp10 (S10), sp13 (S13)]   1.349508   3.070748
Omega1[sp11 (S11), sp13 (S13)]   2.463439  12.218637
Omega1[sp12 (S12), sp13 (S13)]   1.877132   3.266352
Omega1[sp13 (S13), sp13 (S13)]   1.230577   1.581025
Omega1[sp14 (S14), sp13 (S13)]   1.691725   2.979961
Omega1[sp15 (S15), sp13 (S13)]   1.192938   1.516749
Omega1[sp16 (S16), sp13 (S13)]   1.093998   1.258677
Omega1[sp17 (S17), sp13 (S13)]   1.244172   1.608766
Omega1[sp18 (S18), sp13 (S13)]   1.945498   3.129905
Omega1[sp19 (S19), sp13 (S13)]   2.222872   3.873230
Omega1[sp20 (S20), sp13 (S13)]   1.115087   1.337931
Omega1[sp1 (S1), sp14 (S14)]     1.638921   2.510918
Omega1[sp2 (S2), sp14 (S14)]     1.262411   1.720752
Omega1[sp3 (S3), sp14 (S14)]     1.359958   1.981669
Omega1[sp4 (S4), sp14 (S14)]     1.064958   1.173735
Omega1[sp5 (S5), sp14 (S14)]     1.531714   2.459877
Omega1[sp6 (S6), sp14 (S14)]     1.319697   1.818658
Omega1[sp7 (S7), sp14 (S14)]     1.685951   2.585173
Omega1[sp8 (S8), sp14 (S14)]     1.463845   2.272318
Omega1[sp9 (S9), sp14 (S14)]     1.039665   1.118837
Omega1[sp10 (S10), sp14 (S14)]   1.974568   9.687445
Omega1[sp11 (S11), sp14 (S14)]   2.270148  11.565373
Omega1[sp12 (S12), sp14 (S14)]   1.388664   2.250585
Omega1[sp13 (S13), sp14 (S14)]   1.691725   2.979961
Omega1[sp14 (S14), sp14 (S14)]   1.028897   1.067099
Omega1[sp15 (S15), sp14 (S14)]   1.515274   3.858793
Omega1[sp16 (S16), sp14 (S14)]   1.513683   2.312818
Omega1[sp17 (S17), sp14 (S14)]   1.535215   2.524891
Omega1[sp18 (S18), sp14 (S14)]   1.016173   1.049527
Omega1[sp19 (S19), sp14 (S14)]   1.003182   1.007500
Omega1[sp20 (S20), sp14 (S14)]   1.478544   2.298121
Omega1[sp1 (S1), sp15 (S15)]     1.284190   1.753065
Omega1[sp2 (S2), sp15 (S15)]     1.036195   1.107134
Omega1[sp3 (S3), sp15 (S15)]     1.028616   1.072121
Omega1[sp4 (S4), sp15 (S15)]     1.194046   1.529358
Omega1[sp5 (S5), sp15 (S15)]     1.347792   1.995293
Omega1[sp6 (S6), sp15 (S15)]     1.528251   2.843589
Omega1[sp7 (S7), sp15 (S15)]     1.343045   1.911503
Omega1[sp8 (S8), sp15 (S15)]     1.207462   1.598785
Omega1[sp9 (S9), sp15 (S15)]     1.214715   1.458027
Omega1[sp10 (S10), sp15 (S15)]   6.249137  38.287672
Omega1[sp11 (S11), sp15 (S15)]   7.937467  48.667801
Omega1[sp12 (S12), sp15 (S15)]   2.043197   3.805658
Omega1[sp13 (S13), sp15 (S15)]   1.192938   1.516749
Omega1[sp14 (S14), sp15 (S15)]   1.515274   3.858793
Omega1[sp15 (S15), sp15 (S15)]   4.573890  24.490485
Omega1[sp16 (S16), sp15 (S15)]   1.410171   1.987063
Omega1[sp17 (S17), sp15 (S15)]   1.248559   1.658844
Omega1[sp18 (S18), sp15 (S15)]   1.557974   2.707515
Omega1[sp19 (S19), sp15 (S15)]   1.776234   3.253162
Omega1[sp20 (S20), sp15 (S15)]   1.089468   1.260057
Omega1[sp1 (S1), sp16 (S16)]     1.043077   1.121326
Omega1[sp2 (S2), sp16 (S16)]     2.057829   3.322602
Omega1[sp3 (S3), sp16 (S16)]     1.467756   2.108053
Omega1[sp4 (S4), sp16 (S16)]     1.618682   3.222337
Omega1[sp5 (S5), sp16 (S16)]     1.036694   1.038282
Omega1[sp6 (S6), sp16 (S16)]     1.006542   1.007695
Omega1[sp7 (S7), sp16 (S16)]     1.225268   1.568044
Omega1[sp8 (S8), sp16 (S16)]     1.197246   1.525229
Omega1[sp9 (S9), sp16 (S16)]     1.157480   1.400326
Omega1[sp10 (S10), sp16 (S16)]   2.202570   8.011191
Omega1[sp11 (S11), sp16 (S16)]   2.082228   6.695805
Omega1[sp12 (S12), sp16 (S16)]   1.420419   2.005956
Omega1[sp13 (S13), sp16 (S16)]   1.093998   1.258677
Omega1[sp14 (S14), sp16 (S16)]   1.513683   2.312818
Omega1[sp15 (S15), sp16 (S16)]   1.410171   1.987063
Omega1[sp16 (S16), sp16 (S16)]   1.027769   1.067589
Omega1[sp17 (S17), sp16 (S16)]   1.039829   1.115800
Omega1[sp18 (S18), sp16 (S16)]   1.980004   3.127078
Omega1[sp19 (S19), sp16 (S16)]   2.105209   3.471378
Omega1[sp20 (S20), sp16 (S16)]   1.258151   1.639475
Omega1[sp1 (S1), sp17 (S17)]     1.094172   1.248602
Omega1[sp2 (S2), sp17 (S17)]     1.755395   2.718741
Omega1[sp3 (S3), sp17 (S17)]     1.299711   1.728927
Omega1[sp4 (S4), sp17 (S17)]     1.529137   2.742858
Omega1[sp5 (S5), sp17 (S17)]     1.060420   1.138031
Omega1[sp6 (S6), sp17 (S17)]     1.019085   1.021026
Omega1[sp7 (S7), sp17 (S17)]     1.037715   1.110113
Omega1[sp8 (S8), sp17 (S17)]     1.047342   1.128882
Omega1[sp9 (S9), sp17 (S17)]     1.178468   1.465912
Omega1[sp10 (S10), sp17 (S17)]   1.268834   1.317881
Omega1[sp11 (S11), sp17 (S17)]   2.564397  11.608945
Omega1[sp12 (S12), sp17 (S17)]   1.755644   3.033060
Omega1[sp13 (S13), sp17 (S17)]   1.244172   1.608766
Omega1[sp14 (S14), sp17 (S17)]   1.535215   2.524891
Omega1[sp15 (S15), sp17 (S17)]   1.248559   1.658844
Omega1[sp16 (S16), sp17 (S17)]   1.039829   1.115800
Omega1[sp17 (S17), sp17 (S17)]   1.063937   1.162702
Omega1[sp18 (S18), sp17 (S17)]   1.837554   2.931425
Omega1[sp19 (S19), sp17 (S17)]   2.148074   3.675870
Omega1[sp20 (S20), sp17 (S17)]   1.069366   1.201960
Omega1[sp1 (S1), sp18 (S18)]     2.113410   3.441301
Omega1[sp2 (S2), sp18 (S18)]     1.936293   3.157408
Omega1[sp3 (S3), sp18 (S18)]     2.254668   3.829278
Omega1[sp4 (S4), sp18 (S18)]     1.082042   1.217280
Omega1[sp5 (S5), sp18 (S18)]     2.635934   4.860072
Omega1[sp6 (S6), sp18 (S18)]     1.703217   2.600789
Omega1[sp7 (S7), sp18 (S18)]     1.968825   3.072294
Omega1[sp8 (S8), sp18 (S18)]     2.293735   4.094142
Omega1[sp9 (S9), sp18 (S18)]     1.113840   1.316262
Omega1[sp10 (S10), sp18 (S18)]   1.462826   4.627038
Omega1[sp11 (S11), sp18 (S18)]   1.355140   2.891993
Omega1[sp12 (S12), sp18 (S18)]   1.653418   2.704892
Omega1[sp13 (S13), sp18 (S18)]   1.945498   3.129905
Omega1[sp14 (S14), sp18 (S18)]   1.016173   1.049527
Omega1[sp15 (S15), sp18 (S18)]   1.557974   2.707515
Omega1[sp16 (S16), sp18 (S18)]   1.980004   3.127078
Omega1[sp17 (S17), sp18 (S18)]   1.837554   2.931425
Omega1[sp18 (S18), sp18 (S18)]   1.058854   1.170437
Omega1[sp19 (S19), sp18 (S18)]   1.014644   1.045528
Omega1[sp20 (S20), sp18 (S18)]   2.574496   4.650726
Omega1[sp1 (S1), sp19 (S19)]     2.064257   3.295152
Omega1[sp2 (S2), sp19 (S19)]     1.598978   2.457425
Omega1[sp3 (S3), sp19 (S19)]     1.843128   2.993526
Omega1[sp4 (S4), sp19 (S19)]     1.172287   1.498882
Omega1[sp5 (S5), sp19 (S19)]     2.177819   3.909196
Omega1[sp6 (S6), sp19 (S19)]     1.685327   2.550552
Omega1[sp7 (S7), sp19 (S19)]     2.276410   3.666242
Omega1[sp8 (S8), sp19 (S19)]     2.115426   3.705113
Omega1[sp9 (S9), sp19 (S19)]     1.075776   1.213683
Omega1[sp10 (S10), sp19 (S19)]   1.586783   5.706480
Omega1[sp11 (S11), sp19 (S19)]   1.418321   3.353283
Omega1[sp12 (S12), sp19 (S19)]   1.859729   3.346476
Omega1[sp13 (S13), sp19 (S19)]   2.222872   3.873230
Omega1[sp14 (S14), sp19 (S19)]   1.003182   1.007500
Omega1[sp15 (S15), sp19 (S19)]   1.776234   3.253162
Omega1[sp16 (S16), sp19 (S19)]   2.105209   3.471378
Omega1[sp17 (S17), sp19 (S19)]   2.148074   3.675870
Omega1[sp18 (S18), sp19 (S19)]   1.014644   1.045528
Omega1[sp19 (S19), sp19 (S19)]   1.001876   1.005376
Omega1[sp20 (S20), sp19 (S19)]   1.995629   3.431347
Omega1[sp1 (S1), sp20 (S20)]     2.044315   3.589435
Omega1[sp2 (S2), sp20 (S20)]     3.597732   6.921324
Omega1[sp3 (S3), sp20 (S20)]     4.484489   8.623411
Omega1[sp4 (S4), sp20 (S20)]     1.136485   1.403853
Omega1[sp5 (S5), sp20 (S20)]     3.326483   6.505529
Omega1[sp6 (S6), sp20 (S20)]     1.874537   3.205330
Omega1[sp7 (S7), sp20 (S20)]     2.076216   3.445866
Omega1[sp8 (S8), sp20 (S20)]     4.726713   9.203475
Omega1[sp9 (S9), sp20 (S20)]     1.398873   2.126414
Omega1[sp10 (S10), sp20 (S20)]   1.192098   1.398793
Omega1[sp11 (S11), sp20 (S20)]   1.016420   1.022556
Omega1[sp12 (S12), sp20 (S20)]   1.583294   2.333122
Omega1[sp13 (S13), sp20 (S20)]   1.115087   1.337931
Omega1[sp14 (S14), sp20 (S20)]   1.478544   2.298121
Omega1[sp15 (S15), sp20 (S20)]   1.089468   1.260057
Omega1[sp16 (S16), sp20 (S20)]   1.258151   1.639475
Omega1[sp17 (S17), sp20 (S20)]   1.069366   1.201960
Omega1[sp18 (S18), sp20 (S20)]   2.574496   4.650726
Omega1[sp19 (S19), sp20 (S20)]   1.995629   3.431347
Omega1[sp20 (S20), sp20 (S20)]   2.666599   4.952989

if (any(psrf_compt[, "Point est."] > 1.1)) {
  cat("::: {.callout-warning}\n")
  cat("**HMSC convergence concern:** One or more Omega parameters have Rhat > 1.1,")
  cat(" indicating incomplete chain mixing. Interpret these results cautiously.\n")
  cat(":::\n")
}

Warning

HMSC convergence concern: One or more Omega parameters have Rhat > 1.1, indicating incomplete chain mixing. Interpret these results cautiously.

assoc_compt <- computeAssociations(m_hmsc_compt)[[1]]
OmegaCor_compt <- assoc_compt$mean
OmegaCor_compt[assoc_compt$support > 0.05 & assoc_compt$support < 0.95] <- 0
hmsc_species_names_compt <- paste0("sp", rownames(mrIML_mat_compt))
hmsc_mat_compt <- OmegaCor_compt[hmsc_species_names_compt, hmsc_species_names_compt]
rownames(hmsc_mat_compt) <- rownames(mrIML_mat_compt)
colnames(hmsc_mat_compt) <- rownames(mrIML_mat_compt)

g_hmsc_compt <- graph_from_adjacency_matrix(
  hmsc_mat_compt,
  mode     = "undirected",
  weighted = TRUE,
  diag     = FALSE
)
edge_colors_hmsc_compt <- ifelse(E(g_hmsc_compt)$weight < 0, yes = "red", no = "blue")
edge_widths_hmsc_compt <- abs(E(g_hmsc_compt)$weight) * 10

par_old <- par(mfrow = c(1, 3))

plot(
  g_mrIML_compt,
  layout = layout1_compt,
  edge.color = edge_colors_mrIML_compt,
  edge.width = edge_widths_mrIML_compt,
  isolates = TRUE,
  main = "mrIML co-occurrence"
)

plot(
  g_hmsc_compt,
  layout = layout1_compt,
  edge.color = edge_colors_hmsc_compt,
  edge.width = edge_widths_hmsc_compt,
  isolates = TRUE,
  main = "HMSC co-occurrence"
)

plot(
  g_true_compt,
  layout = layout1_compt,
  edge.color = edge_colors_true_compt,
  isolates = TRUE,
  edge.width = edge_widths_true_compt,
  main = "Simulated (truth)"
)

par(par_old)

dist_mat_true_compt <- 1 - true_mat_compt
dist_mat_mrIML_compt <- 1 - mrIML_mat_compt
mantel_mrIML_compt <- vegan::mantel(
  as.dist(dist_mat_true_compt),
  as.dist(dist_mat_mrIML_compt),
  method = "spearman"
)
mantel_mrIML_compt

Mantel statistic based on Spearman's rank correlation rho 

Call:
vegan::mantel(xdis = as.dist(dist_mat_true_compt), ydis = as.dist(dist_mat_mrIML_compt),      method = "spearman") 

Mantel statistic r: 0.8534 
      Significance: 0.001 

Upper quantiles of permutations (null model):
   90%    95%  97.5%    99% 
0.0912 0.1173 0.1417 0.1633 
Permutation: free
Number of permutations: 999

dist_mat_hmsc_compt <- 1 - hmsc_mat_compt
mantel_hmsc_compt <- vegan::mantel(
  as.dist(dist_mat_true_compt),
  as.dist(dist_mat_hmsc_compt),
  method = "spearman"
)
mantel_hmsc_compt

Mantel statistic based on Spearman's rank correlation rho 

Call:
vegan::mantel(xdis = as.dist(dist_mat_true_compt), ydis = as.dist(dist_mat_hmsc_compt),      method = "spearman") 

Mantel statistic r: 0.3985 
      Significance: 0.001 

Upper quantiles of permutations (null model):
   90%    95%  97.5%    99% 
0.0947 0.1182 0.1432 0.1662 
Permutation: free
Number of permutations: 999

5.3 Mixed

simulated_mix_network <- generateM_specific_type(
  nn = network_size,
  k_ave = k_average,
  type.network = "random",
  type.interact = "mix",
  interact.str.max = 0.4,
  mix.compt.ratio = 0.5
)
true_mix_network <- simulated_mix_network[[1]]
M_mix <- simulated_mix_network[[2]]

data_mix <- generateDataSet(
  900,
  M_mix,
  count = network_size * 10000,
  mode = 4
) %>%
  t() %>%
  as.data.frame() %>%
  mutate(
    across(everything(), ~ ifelse(. < 1, yes = 0, no = 1))
  )

# Prepare data
Y_mix <- filterRareCommon(data_mix, lower = 0.01, higher = 0.99)
X1_mix <- Y_mix
# Fit mrIML obj
yhats_rf_sim_mix <- mrIMLpredicts(
  Y = Y_mix,
  X = NULL,
  X1 = X1_mix,
  Model = model_rf,
  prop = 0.7,
  k = 5,
  racing = FALSE
)

i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"
i Creating pre-processing data to finalize 1 unknown parameter: "mtry"

# Bootstrapping
bs_sim_mix <- mrBootstrap(yhats_rf_sim_mix)
# Extract network
assoc_net_mix <- mrCoOccurNet(bs_sim_mix)

mrIML_mat_mix <- assoc_net_mix %>%
  filter(mean_strength > 0.1) %>%
  mutate(
    mean_strength_dir = ifelse(
      direction == "negative",
      yes = -mean_strength,
      no = mean_strength
    ),
    across(contains("taxa"), ~ sub("^sp", "", .))
  ) %>%
  graph_from_data_frame(
    directed = FALSE,
    vertices = sub("^sp", "", names(Y_mix))
  ) %>%
  get.adjacency(attr = "mean_strength_dir", sparse = FALSE)

#plots
g_mrIML_mix <- graph_from_adjacency_matrix(
  mrIML_mat_mix,
  mode = 'undirected',
  weighted = T,
  diag = F
)

edge_colors_mrIML_mix <- ifelse(
  E(g_mrIML_mix)$weight < 0,
  yes = "red",
  no = "blue"
)
edge_widths_mrIML_mix <- abs(E(g_mrIML_mix)$weight) * 10 # Adjust the scaling factor as needed
matching_indices_mix <- as.numeric(rownames(mrIML_mat_mix))

true_mat_mix <- M_mix[matching_indices_mix, matching_indices_mix]
rownames(true_mat_mix) <- rownames(mrIML_mat_mix)
colnames(true_mat_mix) <- rownames(mrIML_mat_mix)
g_true_mix <- graph_from_adjacency_matrix(
  true_mat_mix,
  mode = 'undirected',
  weighted = T,
  diag = F
)
edge_colors_true_mix <- ifelse(
  E(g_true_mix)$weight < 0,
  yes = "red",
  no = "blue"
)
edge_widths_true_mix <- abs(E(g_true_mix)$weight) * 5 # Adjust the scaling factor as needed

layout1_mix <- layout_nicely(g_true_mix)

Warning: Non-positive edge weight found, ignoring all weights during graph
layout.

studyDesign_mix <- data.frame(sample = factor(1:nrow(Y_mix)))
rL_mix <- HmscRandomLevel(units = studyDesign_mix$sample)
m_hmsc_mix <- Hmsc(
  Y           = as.matrix(Y_mix),
  XFormula    = ~1,
  XData       = data.frame(intercept = rep(1, nrow(Y_mix))),
  distr       = "probit",
  studyDesign = studyDesign_mix,
  ranLevels   = list("sample" = rL_mix)
)
m_hmsc_mix <- sampleMcmc(
  m_hmsc_mix,
  thin      = 50,
  samples   = 1000,
  transient = 5000,
  nChains   = 4,
  verbose   = 0
)

mpost_mix <- convertToCodaObject(m_hmsc_mix)
psrf_mix <- gelman.diag(mpost_mix$Omega[[1]], multivariate = FALSE)$psrf
psrf_mix

                               Point est. Upper C.I.
Omega1[sp1 (S1), sp1 (S1)]      1.0257892  1.0571187
Omega1[sp2 (S2), sp1 (S1)]      1.2945092  1.7386810
Omega1[sp3 (S3), sp1 (S1)]      1.0694921  1.1670367
Omega1[sp4 (S4), sp1 (S1)]      1.0222942  1.0672694
Omega1[sp5 (S5), sp1 (S1)]      1.0130829  1.0369024
Omega1[sp6 (S6), sp1 (S1)]      1.0224422  1.0682251
Omega1[sp8 (S7), sp1 (S1)]      1.1447293  1.3891285
Omega1[sp9 (S8), sp1 (S1)]      1.0823786  1.2270406
Omega1[sp11 (S9), sp1 (S1)]     1.0120611  1.0303321
Omega1[sp12 (S10), sp1 (S1)]    1.0492889  1.1417067
Omega1[sp13 (S11), sp1 (S1)]    1.0624815  1.1337832
Omega1[sp15 (S12), sp1 (S1)]    1.1092255  1.3030465
Omega1[sp16 (S13), sp1 (S1)]    1.0926033  1.2557559
Omega1[sp17 (S14), sp1 (S1)]    1.0176373  1.0514458
Omega1[sp18 (S15), sp1 (S1)]    1.0507336  1.1500208
Omega1[sp19 (S16), sp1 (S1)]    1.0427451  1.0992912
Omega1[sp20 (S17), sp1 (S1)]    1.0160344  1.0497903
Omega1[sp1 (S1), sp2 (S2)]      1.2945092  1.7386810
Omega1[sp2 (S2), sp2 (S2)]      1.0482999  1.1099914
Omega1[sp3 (S3), sp2 (S2)]      1.2321819  1.6008927
Omega1[sp4 (S4), sp2 (S2)]      1.1429098  1.3964307
Omega1[sp5 (S5), sp2 (S2)]      1.2685403  1.6656836
Omega1[sp6 (S6), sp2 (S2)]      1.1919846  1.5124761
Omega1[sp8 (S7), sp2 (S2)]      1.0167739  1.0452938
Omega1[sp9 (S8), sp2 (S2)]      1.1524341  1.4337262
Omega1[sp11 (S9), sp2 (S2)]     1.3890030  1.9381467
Omega1[sp12 (S10), sp2 (S2)]    1.1364234  1.3806020
Omega1[sp13 (S11), sp2 (S2)]    1.1165436  1.3257207
Omega1[sp15 (S12), sp2 (S2)]    1.0702614  1.1848510
Omega1[sp16 (S13), sp2 (S2)]    1.4032908  1.9490793
Omega1[sp17 (S14), sp2 (S2)]    1.3663063  1.8979694
Omega1[sp18 (S15), sp2 (S2)]    1.1186717  1.3328516
Omega1[sp19 (S16), sp2 (S2)]    1.4364669  2.0665688
Omega1[sp20 (S17), sp2 (S2)]    1.2013885  1.5401285
Omega1[sp1 (S1), sp3 (S3)]      1.0694921  1.1670367
Omega1[sp2 (S2), sp3 (S3)]      1.2321819  1.6008927
Omega1[sp3 (S3), sp3 (S3)]      1.0215959  1.0660462
Omega1[sp4 (S4), sp3 (S3)]      1.0347669  1.1029241
Omega1[sp5 (S5), sp3 (S3)]      1.1146739  1.2892592
Omega1[sp6 (S6), sp3 (S3)]      1.0309477  1.0934494
Omega1[sp8 (S7), sp3 (S3)]      1.2480095  1.6107649
Omega1[sp9 (S8), sp3 (S3)]      1.1904654  1.5289871
Omega1[sp11 (S9), sp3 (S3)]     1.0219956  1.0242757
Omega1[sp12 (S10), sp3 (S3)]    1.0457945  1.1348252
Omega1[sp13 (S11), sp3 (S3)]    1.0544083  1.1576607
Omega1[sp15 (S12), sp3 (S3)]    1.3747740  1.8911888
Omega1[sp16 (S13), sp3 (S3)]    1.0339575  1.1023669
Omega1[sp17 (S14), sp3 (S3)]    1.0288303  1.0548919
Omega1[sp18 (S15), sp3 (S3)]    1.0075394  1.0243002
Omega1[sp19 (S16), sp3 (S3)]    1.0300076  1.0736514
Omega1[sp20 (S17), sp3 (S3)]    1.0634887  1.1837732
Omega1[sp1 (S1), sp4 (S4)]      1.0222942  1.0672694
Omega1[sp2 (S2), sp4 (S4)]      1.1429098  1.3964307
Omega1[sp3 (S3), sp4 (S4)]      1.0347669  1.1029241
Omega1[sp4 (S4), sp4 (S4)]      1.0027978  1.0102898
Omega1[sp5 (S5), sp4 (S4)]      1.0044611  1.0148079
Omega1[sp6 (S6), sp4 (S4)]      1.0084728  1.0277904
Omega1[sp8 (S7), sp4 (S4)]      1.1045903  1.2769616
Omega1[sp9 (S8), sp4 (S4)]      1.1950912  1.5400701
Omega1[sp11 (S9), sp4 (S4)]     1.2599502  1.6372530
Omega1[sp12 (S10), sp4 (S4)]    1.0157779  1.0484926
Omega1[sp13 (S11), sp4 (S4)]    1.0194660  1.0573504
Omega1[sp15 (S12), sp4 (S4)]    1.2126961  1.5429799
Omega1[sp16 (S13), sp4 (S4)]    1.0504450  1.1453536
Omega1[sp17 (S14), sp4 (S4)]    1.0805248  1.2295828
Omega1[sp18 (S15), sp4 (S4)]    1.0053179  1.0174728
Omega1[sp19 (S16), sp4 (S4)]    1.2460422  1.6060090
Omega1[sp20 (S17), sp4 (S4)]    1.0383131  1.1148268
Omega1[sp1 (S1), sp5 (S5)]      1.0130829  1.0369024
Omega1[sp2 (S2), sp5 (S5)]      1.2685403  1.6656836
Omega1[sp3 (S3), sp5 (S5)]      1.1146739  1.2892592
Omega1[sp4 (S4), sp5 (S5)]      1.0044611  1.0148079
Omega1[sp5 (S5), sp5 (S5)]      1.0140788  1.0350726
Omega1[sp6 (S6), sp5 (S5)]      1.0709812  1.1971272
Omega1[sp8 (S7), sp5 (S5)]      1.1135179  1.3181023
Omega1[sp9 (S8), sp5 (S5)]      1.0570520  1.1471153
Omega1[sp11 (S9), sp5 (S5)]     1.0041402  1.0111453
Omega1[sp12 (S10), sp5 (S5)]    1.0280540  1.0708060
Omega1[sp13 (S11), sp5 (S5)]    1.0577679  1.0972377
Omega1[sp15 (S12), sp5 (S5)]    1.0935552  1.2639160
Omega1[sp16 (S13), sp5 (S5)]    1.1367373  1.3673873
Omega1[sp17 (S14), sp5 (S5)]    1.0020225  1.0075812
Omega1[sp18 (S15), sp5 (S5)]    1.0225757  1.0696874
Omega1[sp19 (S16), sp5 (S5)]    1.0141366  1.0386859
Omega1[sp20 (S17), sp5 (S5)]    1.0652428  1.1867620
Omega1[sp1 (S1), sp6 (S6)]      1.0224422  1.0682251
Omega1[sp2 (S2), sp6 (S6)]      1.1919846  1.5124761
Omega1[sp3 (S3), sp6 (S6)]      1.0309477  1.0934494
Omega1[sp4 (S4), sp6 (S6)]      1.0084728  1.0277904
Omega1[sp5 (S5), sp6 (S6)]      1.0709812  1.1971272
Omega1[sp6 (S6), sp6 (S6)]      1.0013164  1.0050926
Omega1[sp8 (S7), sp6 (S6)]      1.1822931  1.4617803
Omega1[sp9 (S8), sp6 (S6)]      1.1820477  1.4934425
Omega1[sp11 (S9), sp6 (S6)]     1.0265060  1.0748354
Omega1[sp12 (S10), sp6 (S6)]    1.0165355  1.0512172
Omega1[sp13 (S11), sp6 (S6)]    1.0218321  1.0634240
Omega1[sp15 (S12), sp6 (S6)]    1.3201632  1.7753042
Omega1[sp16 (S13), sp6 (S6)]    1.0207517  1.0633952
Omega1[sp17 (S14), sp6 (S6)]    0.9997511  0.9999759
Omega1[sp18 (S15), sp6 (S6)]    1.0077390  1.0195252
Omega1[sp19 (S16), sp6 (S6)]    1.0677364  1.1852649
Omega1[sp20 (S17), sp6 (S6)]    1.0279120  1.0851364
Omega1[sp1 (S1), sp8 (S7)]      1.1447293  1.3891285
Omega1[sp2 (S2), sp8 (S7)]      1.0167739  1.0452938
Omega1[sp3 (S3), sp8 (S7)]      1.2480095  1.6107649
Omega1[sp4 (S4), sp8 (S7)]      1.1045903  1.2769616
Omega1[sp5 (S5), sp8 (S7)]      1.1135179  1.3181023
Omega1[sp6 (S6), sp8 (S7)]      1.1822931  1.4617803
Omega1[sp8 (S7), sp8 (S7)]      1.0066930  1.0191773
Omega1[sp9 (S8), sp8 (S7)]      1.1652407  1.4439572
Omega1[sp11 (S9), sp8 (S7)]     1.2778761  1.6975412
Omega1[sp12 (S10), sp8 (S7)]    1.0936101  1.2537958
Omega1[sp13 (S11), sp8 (S7)]    1.0853084  1.2266635
Omega1[sp15 (S12), sp8 (S7)]    1.0129651  1.0380399
Omega1[sp16 (S13), sp8 (S7)]    1.4503164  2.0408691
Omega1[sp17 (S14), sp8 (S7)]    1.2363844  1.6012829
Omega1[sp18 (S15), sp8 (S7)]    1.0567237  1.1485502
Omega1[sp19 (S16), sp8 (S7)]    1.3572158  1.8489834
Omega1[sp20 (S17), sp8 (S7)]    1.1973269  1.4972945
Omega1[sp1 (S1), sp9 (S8)]      1.0823786  1.2270406
Omega1[sp2 (S2), sp9 (S8)]      1.1524341  1.4337262
Omega1[sp3 (S3), sp9 (S8)]      1.1904654  1.5289871
Omega1[sp4 (S4), sp9 (S8)]      1.1950912  1.5400701
Omega1[sp5 (S5), sp9 (S8)]      1.0570520  1.1471153
Omega1[sp6 (S6), sp9 (S8)]      1.1820477  1.4934425
Omega1[sp8 (S7), sp9 (S8)]      1.1652407  1.4439572
Omega1[sp9 (S8), sp9 (S8)]      1.2879113  2.0372722
Omega1[sp11 (S9), sp9 (S8)]     1.3370293  1.9163929
Omega1[sp12 (S10), sp9 (S8)]    1.1578587  1.4167228
Omega1[sp13 (S11), sp9 (S8)]    1.1568204  1.4158089
Omega1[sp15 (S12), sp9 (S8)]    1.2251474  1.6028475
Omega1[sp16 (S13), sp9 (S8)]    1.2112185  1.6097625
Omega1[sp17 (S14), sp9 (S8)]    1.1682494  1.4601024
Omega1[sp18 (S15), sp9 (S8)]    1.1321299  1.3073987
Omega1[sp19 (S16), sp9 (S8)]    1.2961901  1.8234453
Omega1[sp20 (S17), sp9 (S8)]    1.2189916  1.6379443
Omega1[sp1 (S1), sp11 (S9)]     1.0120611  1.0303321
Omega1[sp2 (S2), sp11 (S9)]     1.3890030  1.9381467
Omega1[sp3 (S3), sp11 (S9)]     1.0219956  1.0242757
Omega1[sp4 (S4), sp11 (S9)]     1.2599502  1.6372530
Omega1[sp5 (S5), sp11 (S9)]     1.0041402  1.0111453
Omega1[sp6 (S6), sp11 (S9)]     1.0265060  1.0748354
Omega1[sp8 (S7), sp11 (S9)]     1.2778761  1.6975412
Omega1[sp9 (S8), sp11 (S9)]     1.3370293  1.9163929
Omega1[sp11 (S9), sp11 (S9)]    1.0090207  1.0241501
Omega1[sp12 (S10), sp11 (S9)]   1.2651504  1.6479522
Omega1[sp13 (S11), sp11 (S9)]   1.2054215  1.5301799
Omega1[sp15 (S12), sp11 (S9)]   1.2090602  1.5393191
Omega1[sp16 (S13), sp11 (S9)]   1.0312181  1.0942602
Omega1[sp17 (S14), sp11 (S9)]   1.0099207  1.0299654
Omega1[sp18 (S15), sp11 (S9)]   1.3349759  1.8010506
Omega1[sp19 (S16), sp11 (S9)]   1.0176030  1.0422582
Omega1[sp20 (S17), sp11 (S9)]   1.0427051  1.1268413
Omega1[sp1 (S1), sp12 (S10)]    1.0492889  1.1417067
Omega1[sp2 (S2), sp12 (S10)]    1.1364234  1.3806020
Omega1[sp3 (S3), sp12 (S10)]    1.0457945  1.1348252
Omega1[sp4 (S4), sp12 (S10)]    1.0157779  1.0484926
Omega1[sp5 (S5), sp12 (S10)]    1.0280540  1.0708060
Omega1[sp6 (S6), sp12 (S10)]    1.0165355  1.0512172
Omega1[sp8 (S7), sp12 (S10)]    1.0936101  1.2537958
Omega1[sp9 (S8), sp12 (S10)]    1.1578587  1.4167228
Omega1[sp11 (S9), sp12 (S10)]   1.2651504  1.6479522
Omega1[sp12 (S10), sp12 (S10)]  1.0136892  1.0428069
Omega1[sp13 (S11), sp12 (S10)]  1.0275638  1.0838574
Omega1[sp15 (S12), sp12 (S10)]  1.1864453  1.4823285
Omega1[sp16 (S13), sp12 (S10)]  1.0686218  1.1932571
Omega1[sp17 (S14), sp12 (S10)]  1.1145996  1.3075588
Omega1[sp18 (S15), sp12 (S10)]  1.0063411  1.0199913
Omega1[sp19 (S16), sp12 (S10)]  1.2451847  1.6080299
Omega1[sp20 (S17), sp12 (S10)]  1.0415962  1.1235358
Omega1[sp1 (S1), sp13 (S11)]    1.0624815  1.1337832
Omega1[sp2 (S2), sp13 (S11)]    1.1165436  1.3257207
Omega1[sp3 (S3), sp13 (S11)]    1.0544083  1.1576607
Omega1[sp4 (S4), sp13 (S11)]    1.0194660  1.0573504
Omega1[sp5 (S5), sp13 (S11)]    1.0577679  1.0972377
Omega1[sp6 (S6), sp13 (S11)]    1.0218321  1.0634240
Omega1[sp8 (S7), sp13 (S11)]    1.0853084  1.2266635
Omega1[sp9 (S8), sp13 (S11)]    1.1568204  1.4158089
Omega1[sp11 (S9), sp13 (S11)]   1.2054215  1.5301799
Omega1[sp12 (S10), sp13 (S11)]  1.0275638  1.0838574
Omega1[sp13 (S11), sp13 (S11)]  1.0086587  1.0274524
Omega1[sp15 (S12), sp13 (S11)]  1.1591999  1.4143051
Omega1[sp16 (S13), sp13 (S11)]  1.0803653  1.2248754
Omega1[sp17 (S14), sp13 (S11)]  1.1090872  1.2759723
Omega1[sp18 (S15), sp13 (S11)]  1.0039476  1.0070450
Omega1[sp19 (S16), sp13 (S11)]  1.2316696  1.5835294
Omega1[sp20 (S17), sp13 (S11)]  1.0515499  1.1503345
Omega1[sp1 (S1), sp15 (S12)]    1.1092255  1.3030465
Omega1[sp2 (S2), sp15 (S12)]    1.0702614  1.1848510
Omega1[sp3 (S3), sp15 (S12)]    1.3747740  1.8911888
Omega1[sp4 (S4), sp15 (S12)]    1.2126961  1.5429799
Omega1[sp5 (S5), sp15 (S12)]    1.0935552  1.2639160
Omega1[sp6 (S6), sp15 (S12)]    1.3201632  1.7753042
Omega1[sp8 (S7), sp15 (S12)]    1.0129651  1.0380399
Omega1[sp9 (S8), sp15 (S12)]    1.2251474  1.6028475
Omega1[sp11 (S9), sp15 (S12)]   1.2090602  1.5393191
Omega1[sp12 (S10), sp15 (S12)]  1.1864453  1.4823285
Omega1[sp13 (S11), sp15 (S12)]  1.1591999  1.4143051
Omega1[sp15 (S12), sp15 (S12)]  1.0114149  1.0331204
Omega1[sp16 (S13), sp15 (S12)]  1.5754150  2.2989660
Omega1[sp17 (S14), sp15 (S12)]  1.1754622  1.4587705
Omega1[sp18 (S15), sp15 (S12)]  1.1593917  1.4158850
Omega1[sp19 (S16), sp15 (S12)]  1.2688354  1.6671668
Omega1[sp20 (S17), sp15 (S12)]  1.3357229  1.8105796
Omega1[sp1 (S1), sp16 (S13)]    1.0926033  1.2557559
Omega1[sp2 (S2), sp16 (S13)]    1.4032908  1.9490793
Omega1[sp3 (S3), sp16 (S13)]    1.0339575  1.1023669
Omega1[sp4 (S4), sp16 (S13)]    1.0504450  1.1453536
Omega1[sp5 (S5), sp16 (S13)]    1.1367373  1.3673873
Omega1[sp6 (S6), sp16 (S13)]    1.0207517  1.0633952
Omega1[sp8 (S7), sp16 (S13)]    1.4503164  2.0408691
Omega1[sp9 (S8), sp16 (S13)]    1.2112185  1.6097625
Omega1[sp11 (S9), sp16 (S13)]   1.0312181  1.0942602
Omega1[sp12 (S10), sp16 (S13)]  1.0686218  1.1932571
Omega1[sp13 (S11), sp16 (S13)]  1.0803653  1.2248754
Omega1[sp15 (S12), sp16 (S13)]  1.5754150  2.2989660
Omega1[sp16 (S13), sp16 (S13)]  1.0003241  1.0024565
Omega1[sp17 (S14), sp16 (S13)]  1.0624011  1.1794480
Omega1[sp18 (S15), sp16 (S13)]  1.0137732  1.0433540
Omega1[sp19 (S16), sp16 (S13)]  1.0115066  1.0316549
Omega1[sp20 (S17), sp16 (S13)]  1.0522109  1.1532417
Omega1[sp1 (S1), sp17 (S14)]    1.0176373  1.0514458
Omega1[sp2 (S2), sp17 (S14)]    1.3663063  1.8979694
Omega1[sp3 (S3), sp17 (S14)]    1.0288303  1.0548919
Omega1[sp4 (S4), sp17 (S14)]    1.0805248  1.2295828
Omega1[sp5 (S5), sp17 (S14)]    1.0020225  1.0075812
Omega1[sp6 (S6), sp17 (S14)]    0.9997511  0.9999759
Omega1[sp8 (S7), sp17 (S14)]    1.2363844  1.6012829
Omega1[sp9 (S8), sp17 (S14)]    1.1682494  1.4601024
Omega1[sp11 (S9), sp17 (S14)]   1.0099207  1.0299654
Omega1[sp12 (S10), sp17 (S14)]  1.1145996  1.3075588
Omega1[sp13 (S11), sp17 (S14)]  1.1090872  1.2759723
Omega1[sp15 (S12), sp17 (S14)]  1.1754622  1.4587705
Omega1[sp16 (S13), sp17 (S14)]  1.0624011  1.1794480
Omega1[sp17 (S14), sp17 (S14)]  1.0102629  1.0323279
Omega1[sp18 (S15), sp17 (S14)]  1.1276162  1.3508692
Omega1[sp19 (S16), sp17 (S14)]  1.0304127  1.0873667
Omega1[sp20 (S17), sp17 (S14)]  1.0039407  1.0075576
Omega1[sp1 (S1), sp18 (S15)]    1.0507336  1.1500208
Omega1[sp2 (S2), sp18 (S15)]    1.1186717  1.3328516
Omega1[sp3 (S3), sp18 (S15)]    1.0075394  1.0243002
Omega1[sp4 (S4), sp18 (S15)]    1.0053179  1.0174728
Omega1[sp5 (S5), sp18 (S15)]    1.0225757  1.0696874
Omega1[sp6 (S6), sp18 (S15)]    1.0077390  1.0195252
Omega1[sp8 (S7), sp18 (S15)]    1.0567237  1.1485502
Omega1[sp9 (S8), sp18 (S15)]    1.1321299  1.3073987
Omega1[sp11 (S9), sp18 (S15)]   1.3349759  1.8010506
Omega1[sp12 (S10), sp18 (S15)]  1.0063411  1.0199913
Omega1[sp13 (S11), sp18 (S15)]  1.0039476  1.0070450
Omega1[sp15 (S12), sp18 (S15)]  1.1593917  1.4158850
Omega1[sp16 (S13), sp18 (S15)]  1.0137732  1.0433540
Omega1[sp17 (S14), sp18 (S15)]  1.1276162  1.3508692
Omega1[sp18 (S15), sp18 (S15)]  1.0241772  1.0660201
Omega1[sp19 (S16), sp18 (S15)]  1.2479745  1.6119658
Omega1[sp20 (S17), sp18 (S15)]  1.0025367  1.0097060
Omega1[sp1 (S1), sp19 (S16)]    1.0427451  1.0992912
Omega1[sp2 (S2), sp19 (S16)]    1.4364669  2.0665688
Omega1[sp3 (S3), sp19 (S16)]    1.0300076  1.0736514
Omega1[sp4 (S4), sp19 (S16)]    1.2460422  1.6060090
Omega1[sp5 (S5), sp19 (S16)]    1.0141366  1.0386859
Omega1[sp6 (S6), sp19 (S16)]    1.0677364  1.1852649
Omega1[sp8 (S7), sp19 (S16)]    1.3572158  1.8489834
Omega1[sp9 (S8), sp19 (S16)]    1.2961901  1.8234453
Omega1[sp11 (S9), sp19 (S16)]   1.0176030  1.0422582
Omega1[sp12 (S10), sp19 (S16)]  1.2451847  1.6080299
Omega1[sp13 (S11), sp19 (S16)]  1.2316696  1.5835294
Omega1[sp15 (S12), sp19 (S16)]  1.2688354  1.6671668
Omega1[sp16 (S13), sp19 (S16)]  1.0115066  1.0316549
Omega1[sp17 (S14), sp19 (S16)]  1.0304127  1.0873667
Omega1[sp18 (S15), sp19 (S16)]  1.2479745  1.6119658
Omega1[sp19 (S16), sp19 (S16)]  1.0354465  1.0757268
Omega1[sp20 (S17), sp19 (S16)]  1.0841176  1.2315633
Omega1[sp1 (S1), sp20 (S17)]    1.0160344  1.0497903
Omega1[sp2 (S2), sp20 (S17)]    1.2013885  1.5401285
Omega1[sp3 (S3), sp20 (S17)]    1.0634887  1.1837732
Omega1[sp4 (S4), sp20 (S17)]    1.0383131  1.1148268
Omega1[sp5 (S5), sp20 (S17)]    1.0652428  1.1867620
Omega1[sp6 (S6), sp20 (S17)]    1.0279120  1.0851364
Omega1[sp8 (S7), sp20 (S17)]    1.1973269  1.4972945
Omega1[sp9 (S8), sp20 (S17)]    1.2189916  1.6379443
Omega1[sp11 (S9), sp20 (S17)]   1.0427051  1.1268413
Omega1[sp12 (S10), sp20 (S17)]  1.0415962  1.1235358
Omega1[sp13 (S11), sp20 (S17)]  1.0515499  1.1503345
Omega1[sp15 (S12), sp20 (S17)]  1.3357229  1.8105796
Omega1[sp16 (S13), sp20 (S17)]  1.0522109  1.1532417
Omega1[sp17 (S14), sp20 (S17)]  1.0039407  1.0075576
Omega1[sp18 (S15), sp20 (S17)]  1.0025367  1.0097060
Omega1[sp19 (S16), sp20 (S17)]  1.0841176  1.2315633
Omega1[sp20 (S17), sp20 (S17)]  1.0369070  1.1099935

if (any(psrf_mix[, "Point est."] > 1.1)) {
  cat("::: {.callout-warning}\n")
  cat("**HMSC convergence concern:** One or more Omega parameters have Rhat > 1.1,")
  cat(" indicating incomplete chain mixing. Interpret these results cautiously.\n")
  cat(":::\n")
}

Warning

HMSC convergence concern: One or more Omega parameters have Rhat > 1.1, indicating incomplete chain mixing. Interpret these results cautiously.

assoc_mix <- computeAssociations(m_hmsc_mix)[[1]]
OmegaCor_mix <- assoc_mix$mean
OmegaCor_mix[assoc_mix$support > 0.05 & assoc_mix$support < 0.95] <- 0
hmsc_species_names_mix <- paste0("sp", rownames(mrIML_mat_mix))
hmsc_mat_mix <- OmegaCor_mix[hmsc_species_names_mix, hmsc_species_names_mix]
rownames(hmsc_mat_mix) <- rownames(mrIML_mat_mix)
colnames(hmsc_mat_mix) <- rownames(mrIML_mat_mix)

g_hmsc_mix <- graph_from_adjacency_matrix(
  hmsc_mat_mix,
  mode     = "undirected",
  weighted = TRUE,
  diag     = FALSE
)
edge_colors_hmsc_mix <- ifelse(E(g_hmsc_mix)$weight < 0, yes = "red", no = "blue")
edge_widths_hmsc_mix <- abs(E(g_hmsc_mix)$weight) * 10

par_old <- par(mfrow = c(1, 3))

plot(
  g_mrIML_mix,
  layout = layout1_mix,
  edge.color = edge_colors_mrIML_mix,
  edge.width = edge_widths_mrIML_mix,
  isolates = TRUE,
  main = "mrIML co-occurrence"
)

plot(
  g_hmsc_mix,
  layout = layout1_mix,
  edge.color = edge_colors_hmsc_mix,
  edge.width = edge_widths_hmsc_mix,
  isolates = TRUE,
  main = "HMSC co-occurrence"
)

plot(
  g_true_mix,
  layout = layout1_mix,
  edge.color = edge_colors_true_mix,
  isolates = TRUE,
  edge.width = edge_widths_true_mix,
  main = "Simulated (truth)"
)

par(par_old)

dist_mat_true_mix <- 1 - true_mat_mix
dist_mat_mrIML_mix <- 1 - mrIML_mat_mix
mantel_mrIML_mix <- vegan::mantel(
  as.dist(dist_mat_true_mix),
  as.dist(dist_mat_mrIML_mix),
  method = "spearman"
)
mantel_mrIML_mix

Mantel statistic based on Spearman's rank correlation rho 

Call:
vegan::mantel(xdis = as.dist(dist_mat_true_mix), ydis = as.dist(dist_mat_mrIML_mix),      method = "spearman") 

Mantel statistic r: 0.6892 
      Significance: 0.001 

Upper quantiles of permutations (null model):
  90%   95% 97.5%   99% 
0.119 0.157 0.181 0.219 
Permutation: free
Number of permutations: 999

dist_mat_hmsc_mix <- 1 - hmsc_mat_mix
mantel_hmsc_mix <- vegan::mantel(
  as.dist(dist_mat_true_mix),
  as.dist(dist_mat_hmsc_mix),
  method = "spearman"
)
mantel_hmsc_mix

Mantel statistic based on Spearman's rank correlation rho 

Call:
vegan::mantel(xdis = as.dist(dist_mat_true_mix), ydis = as.dist(dist_mat_hmsc_mix),      method = "spearman") 

Mantel statistic r: 0.5004 
      Significance: 0.001 

Upper quantiles of permutations (null model):
  90%   95% 97.5%   99% 
0.113 0.143 0.184 0.233 
Permutation: free
Number of permutations: 999

5.4 Summary

Comparison of Mantel r (Pearson) and p-values for mrIML and HMSC against the true interaction network across all simulation scenarios.

comparison_table <- tibble(
  Scenario = c("Random", "Mutualistic", "Competitive", "Mixed"),
  mrIML_r  = round(c(
    mantel_mrIML$statistic,
    mantel_mrIML_mutual$statistic,
    mantel_mrIML_compt$statistic,
    mantel_mrIML_mix$statistic
  ), 3),
  mrIML_p  = c(
    mantel_mrIML$signif,
    mantel_mrIML_mutual$signif,
    mantel_mrIML_compt$signif,
    mantel_mrIML_mix$signif
  ),
  HMSC_r   = round(c(
    mantel_hmsc$statistic,
    mantel_hmsc_mutual$statistic,
    mantel_hmsc_compt$statistic,
    mantel_hmsc_mix$statistic
  ), 3),
  HMSC_p   = c(
    mantel_hmsc$signif,
    mantel_hmsc_mutual$signif,
    mantel_hmsc_compt$signif,
    mantel_hmsc_mix$signif
  )
)

DT::datatable(
  comparison_table,
  colnames = c("Scenario", "mrIML r", "mrIML p", "HMSC r", "HMSC p"),
  options  = list(pageLength = 10, dom = "t")
)

6 Asymetric simulations

…