diff --git a/multipartite_fungustree.R b/multipartite_fungustree.R index 45bcc37..870b537 100644 --- a/multipartite_fungustree.R +++ b/multipartite_fungustree.R @@ -5,8 +5,13 @@ library(sbm) data(fungusTreeNetwork) tree_genetic_dist_matrix <- fungusTreeNetwork$covar_tree$genetic_dist +tree_taxonomic_dist_matrix <- (fungusTreeNetwork$covar_tree$taxonomic_dist) + + tree_genet_sbm <- estimateSimpleSBM(netMat = tree_genetic_dist_matrix, model = "gaussian", estimOptions = list(plot = FALSE)) plot(tree_genet_sbm) +tree_taxonomic_sbm <- estimateSimpleSBM(netMat = tree_taxonomic_dist_matrix, model = "gaussian", estimOptions = list(plot = FALSE)) +plot(tree_taxonomic_sbm) tree_geo_sbm <- estimateSimpleSBM(netMat = fungusTreeNetwork$covar_tree$geographic_dist, model = "gaussian", estimOptions = list(plot = FALSE)) plot(tree_geo_sbm) @@ -22,6 +27,14 @@ namesFG <- c("Tree", "Fungus") givenclassif <- list(tree_geo_sbm$memberships, fungus_tree_sbm$memberships[["row"]]) geo_multi_noinit <- multipartiteBM(list_Net = list_Net, v_distrib = v_distrib, namesFG = namesFG, givenclassif = NULL, initBM = TRUE, verbose = TRUE) + +Tree_Geo_dist_sbm <- defineSBM(-fungusTreeNetwork$covar_tree$geographic_dist, model = "gaussian", type = "simple", dimLabels = c("Tree")) +Tree_Taxo_dist_sbm <- defineSBM(-fungusTreeNetwork$covar_tree$taxonomic_dist, model = "gaussian", type = "simple", dimLabels = c("Tree")) +Tree_Genetic_dist_sbm <- defineSBM(fungusTreeNetwork$covar_tree$genetic_dist, model = "gaussian", type = "simple", dimLabels = c("Tree")) +FungusTree_sbm <- defineSBM(fungusTreeNetwork$fungus_tree, dimLabels = c("Fungus", "Tree"), model = "bernoulli", type = "bipartite") + +geo_multi_noinit_sbm <- estimateMultipartiteSBM(list(Tree_Geo_dist_sbm, FungusTree_sbm), estimOptions = list(initBM = TRUE, )) + geo_multi_noinit_nobm <- multipartiteBM(list_Net = list_Net, v_distrib = v_distrib, namesFG = namesFG, givenclassif = NULL, initBM = FALSE, verbose = TRUE) # Les kmins donnent le même résultat ICL -449.21 geo_multi_noinit_nobm_kmin <- multipartiteBM(list_Net = list_Net, v_distrib = v_distrib, namesFG = namesFG, givenclassif = NULL, initBM = FALSE, verbose = TRUE, v_Kmin = c(7, 2)) @@ -32,4 +45,15 @@ geo_multi_init_kmin <- multipartiteBM(list_Net = list_Net, v_distrib = v_distrib # Le meilleur ICL est ici : -408.05 avec une classif de base issue de blockmodels indépendants geo_multi_init <- multipartiteBM(list_Net = list_Net, v_distrib = v_distrib, namesFG = namesFG, givenclassif = givenclassif) +library(knitr) + +kable(data.frame("Given classif" = c(FALSE, FALSE, FALSE, FALSE, TRUE, TRUE), "InitBM" = c(TRUE, FALSE, FALSE, TRUE, TRUE, TRUE), "Kmin given" = c(FALSE, FALSE, TRUE, TRUE, TRUE, FALSE), ICL = c( + geo_multi_noinit$fittedModel[[1]]$ICL, + geo_multi_noinit_nobm$fittedModel[[1]]$ICL, + geo_multi_noinit_nobm_kmin$fittedModel[[1]]$ICL, + geo_multi_noinit_kmin$fittedModel[[1]]$ICL, + geo_multi_init_kmin$fittedModel[[1]]$ICL, + geo_multi_init$fittedModel[[1]]$ICL +)), format = "markdown") + dataR6 <- formattingData(list_Net, v_distrib)