mia-rapport-2024/Rcodes/real_data/netclustering_dore_analyze.R
2024-06-28 10:49:49 +02:00

76 lines
3 KiB
R

require("ggplot2")
devtools::load_all("R/")
## Preparing data for analysis
interaction_data <- read.table(file = "real_data/data/interaction-data.txt", sep = "\t", header = TRUE)
insect_orders <- unique(interaction_data$insectorder)
plant_family <- unique(interaction_data$plantfamily)
# All results
iid_clustering <- readRDS("real_data/data/dore_collection_clustering_nb_run1_iid_123networks_24-05-23-21:40:42.Rds")
iid_best_partition <- extract_best_bipartite_partition(iid_clustering)
iid_unlist <- unlist(iid_best_partition)
rho_clustering <- readRDS("real_data/data/dore_collection_clustering_nb_run1_rho_123networks_25-05-23-13:58:30.Rds")
rho_best_partition <- extract_best_bipartite_partition(rho_clustering)
rho_unlist <- unlist(rho_best_partition)
pi_clustering <- readRDS("real_data/data/dore_collection_clustering_nb_run1_pi_123networks_25-05-23-17:31:25.Rds")
pi_best_partition <- extract_best_bipartite_partition(pi_clustering)
pi_unlist <- unlist(pi_best_partition)
pirho_clustering <- readRDS("real_data/data/dore_collection_clustering_nb_run1_pirho_123networks_26-05-23-19:22:55.Rds")
pirho_best_partition <- extract_best_bipartite_partition(pirho_clustering)
pirho_unlist <- unlist(pirho_best_partition)
### Matching taxonomy
taxonomy_in_clusters <- function(unlisted_model) {
lapply(seq_len(length(unlisted_model)), function(col_idx) {
# Per collection
# Empty init
insect_count <- t(sapply(insect_orders, function(order) {
out_count <- rep(0, unlisted_model[[col_idx]]$Q[2])
out_count
}))
plant_count <- t(sapply(plant_family, function(order) {
out_count <- rep(0, unlisted_model[[col_idx]]$Q[1])
out_count
}))
for (m in seq.int(unlisted_model[[col_idx]]$M)) {
#### Insect
insect_names <- names(unlisted_model[[col_idx]]$Z[[1]][[2]])
insect_count <- insect_count + t(sapply(insect_orders, function(order) {
out_count <- rep(0, unlisted_model[[col_idx]]$Q[2])
names(out_count) <- seq.int(unlisted_model[[col_idx]]$Q[2])
insect_count <- table(unlisted_model[[col_idx]]$Z[[m]][[2]][grep(order, insect_names)])
out_count[names(insect_count)] <- insect_count
out_count
}))
#### Plants
plant_names <- names(unlisted_model[[col_idx]]$Z[[1]][[1]])
plant_count <- t(sapply(plant_family, function(order) {
out_count <- rep(0, unlisted_model[[col_idx]]$Q[1])
names(out_count) <- seq.int(unlisted_model[[col_idx]]$Q[1])
plant_count <- table(unlisted_model[[col_idx]]$Z[[m]][[1]][grep(order, plant_names)])
out_count[names(plant_count)] <- plant_count
out_count
}))
}
return(list(insects = insect_count, plants = plant_count))
})
}
taxonomy_remove_empty <- function(taxonomy_collections_list) {
lapply(taxonomy_collections_list, function(collection) {
list(
insects = collection$insects[which(rowSums(collection$insects != 0) > 0), ],
plants = collection$plants[which(rowSums(collection$plants != 0) > 0), ]
)
})
}