.Rhistory

C[rownames(P), colnames(P)] <- ifelse(is.na(C[rownames(P), colnames(P)]), NA, C[rownames(P), colnames(P)] + P[rownames(P), colnames(P)])
# C[C>0] <- 1 # Make binary
sum(is.na(C))
# might need to convert C into a binary matrix
### ---- transfer learning with SVD ----
# Define the grid of k and lambda values to search over
k_values <- c(2, 5, 10)            # Adjust as needed
lambda_values <- c(0.01, 0.05, 0.1)  # Adjust as needed
# Initialize variables to store the best results
results <- data.frame(k = integer(),
lambda = numeric(),
original_links = numeric(),
predicted_values = numeric())
# Loop over all combinations of k and lambda
for (k in k_values) {
for (lambda in lambda_values) {
# Apply softImpute
fit <- softImpute(C, rank.max = k, lambda = lambda, type = "svd", maxit = 600)
# Reconstruct the matrix
C_reconstructed <- softImpute::complete(C, fit)
# Extract the reconstructed P matrix from C_reconstructed
P_reconstructed <- C_reconstructed[rownames(P), colnames(P)]
# Combine indices of removed ones and zeros
test_indices <- rbind(
data.frame(row = remove_indices[, "row"], col = remove_indices[, "col"], label = rep(1, nrow(remove_indices))),
data.frame(row = zeros_to_remove_indices[, "row"], col = zeros_to_remove_indices[, "col"], label = rep(0, nrow(zeros_to_remove_indices)))
)
# Get the row and column names of the test links
test_rows <- rownames(P)[test_indices$row]  # These are the "node_to"
test_cols <- colnames(P)[test_indices$col]  # These are the "node_from"
# Actual labels and predictions
original_links <- P_original[cbind(test_rows, test_cols)]
predicted_values <- P_reconstructed[cbind(test_rows, test_cols)]
# Store the results with node information
results <- rbind(results, data.frame(
k = k,
lambda = lambda,
original_links = original_links,
predicted_values = predicted_values,
node_to = test_rows,
node_from = test_cols,
removed = 1 # mark these links as removed
))
# ---- (A) Enumerate ALL edges in P
all_edges <- expand.grid(
node_to = rownames(P),
node_from = colnames(P),
k = k,
lambda = lambda,
KEEP.OUT.ATTRS = FALSE,
stringsAsFactors = FALSE
)
# Fill in the original link values from P_original
all_edges$original_links <- mapply(
function(r, c) P_original[r, c],
all_edges$node_to,
all_edges$node_from
)
# Helper data frame of removed edges
removed_edges_idx <- data.frame(
node_to = test_rows,
node_from = test_cols,
stringsAsFactors = FALSE
)
# ---- (B) Subset edges NOT removed
not_removed <- all_edges[
!paste(all_edges$node_to, all_edges$node_from) %in%
paste(removed_edges_idx$node_to, removed_edges_idx$node_from),
]
not_removed$removed <- 0
not_removed$predicted_values <- NA
# ---- (C) Combine removed + not removed
not_removed$k <- k
not_removed$lambda <- lambda
complete_edges <- rbind(
results,     # removed edges with predictions
not_removed  # not-removed edges
)
print("before temp writting!")
write.csv(complete_edges,
file = "temp.csv",
row.names = FALSE)
}
}
### ---- save results for current combination ----
# After finishing the k/lambda loops, append the 'results' to 'combined_results'
# ---- (D) Append to combined_results
complete_edges_all <- read.csv("temp.csv")
combined_results <- rbind(
combined_results,
cbind(
data.frame(
emln_id = emln_id,
train_layer = layers_to_train,
test_layer = layer_to_predict,
prop_ones_removed = prop_ones_to_remove,
amount_of_removed_1 = num_1_to_remove,
prop_zeros_removed = prop_zeros_to_remove,
amount_of_removed_0 = num_0_to_remove
),
complete_edges_all
)
)
remove("temp.csv")
}
}
# Loop through all combinations of layers_to_train and layer_to_predict
for (layers_to_train in 1:num_layers) {
for (layer_to_predict in 1:num_layers) {
# Parameters for the current combination
emln_id <- 25
prop_ones_to_remove <- 0.2
prop_zeros_to_remove <- 0.2
# Load matrices
d <- load_emln(emln_id)
graph_list <- get_igraph(d, bipartite = TRUE, directed = FALSE)$layers_igraph
A_l <- d$extended
# Total number of layers
num_layers <- length(graph_list)
# Build the aggregated matrix A for training
A <- build_interaction_matrix(data = A_l, layers_to_filter = layers_to_train)
#A[A > 0] <- 1  # Make binary
# Build the layer to predict matrix P
P <- build_interaction_matrix(data = A_l, layers_to_filter = layer_to_predict)
#P[P > 0] <- 1  # Make binary
node_to <- rownames(P) # for the results
node_from <- colnames(P)
# # # crop to make toy example matrices
#A <- A[1:4, 1:3]
#P <- P[1:4, 1:3]
### ---- remove some links in P ----
P_original <- P # save it for later
num_1_to_remove <- floor(sum(P, na.rm = T)*prop_ones_to_remove)  # Number of links to remove
ones_in_P <- which(P > 0, arr.ind = TRUE)
remove_indices <- ones_in_P[sample(1:nrow(ones_in_P), num_1_to_remove), ]
P[remove_indices] <- NA  # Set removed links to NA
# Indices of zeros in P
num_0_to_remove <- floor((nrow(P)*ncol(P)-sum(P, na.rm = T))*prop_zeros_to_remove)  # Number of links to remove
zeros_in_P <- which(P == 0, arr.ind = TRUE)
# Randomly select zeros to remove
# set.seed(789)
zeros_to_remove_indices <- zeros_in_P[sample(1:nrow(zeros_in_P), num_0_to_remove), ]
# Set the selected zeros to NA
P[zeros_to_remove_indices] <- NA
### ---- creating a combined matrix C ----
# Combine A and P into a single matrix C with NAs representing missing data
all_row_ids <- unique(c(rownames(A), rownames(P)))
all_col_ids <- unique(c(colnames(A), colnames(P)))
C <- matrix(0, nrow = length(all_row_ids), ncol = length(all_col_ids),
dimnames = list(all_row_ids, all_col_ids))
# Place A into C
C[rownames(A), colnames(A)] <- A
# Place P into C
# Ensure that existing entries are not overwritten; sum overlapping entries
C[rownames(P), colnames(P)] <- ifelse(is.na(C[rownames(P), colnames(P)]), NA, C[rownames(P), colnames(P)] + P[rownames(P), colnames(P)])
# C[C>0] <- 1 # Make binary
sum(is.na(C))
# might need to convert C into a binary matrix
### ---- transfer learning with SVD ----
# Define the grid of k and lambda values to search over
k_values <- c(2, 5, 10)            # Adjust as needed
lambda_values <- c(0.01, 0.05, 0.1)  # Adjust as needed
# Initialize variables to store the best results
results <- data.frame(k = integer(),
lambda = numeric(),
original_links = numeric(),
predicted_values = numeric())
# Loop over all combinations of k and lambda
for (k in k_values) {
for (lambda in lambda_values) {
# Apply softImpute
fit <- softImpute(C, rank.max = k, lambda = lambda, type = "svd", maxit = 600)
# Reconstruct the matrix
C_reconstructed <- softImpute::complete(C, fit)
# Extract the reconstructed P matrix from C_reconstructed
P_reconstructed <- C_reconstructed[rownames(P), colnames(P)]
# Combine indices of removed ones and zeros
test_indices <- rbind(
data.frame(row = remove_indices[, "row"], col = remove_indices[, "col"], label = rep(1, nrow(remove_indices))),
data.frame(row = zeros_to_remove_indices[, "row"], col = zeros_to_remove_indices[, "col"], label = rep(0, nrow(zeros_to_remove_indices)))
)
# Get the row and column names of the test links
test_rows <- rownames(P)[test_indices$row]  # These are the "node_to"
test_cols <- colnames(P)[test_indices$col]  # These are the "node_from"
# Actual labels and predictions
original_links <- P_original[cbind(test_rows, test_cols)]
predicted_values <- P_reconstructed[cbind(test_rows, test_cols)]
# Store the results with node information
results <- rbind(results, data.frame(
k = k,
lambda = lambda,
original_links = original_links,
predicted_values = predicted_values,
node_to = test_rows,
node_from = test_cols,
removed = 1 # mark these links as removed
))
# ---- (A) Enumerate ALL edges in P
all_edges <- expand.grid(
node_to = rownames(P),
node_from = colnames(P),
k = k,
lambda = lambda,
KEEP.OUT.ATTRS = FALSE,
stringsAsFactors = FALSE
)
# Fill in the original link values from P_original
all_edges$original_links <- mapply(
function(r, c) P_original[r, c],
all_edges$node_to,
all_edges$node_from
)
# Helper data frame of removed edges
removed_edges_idx <- data.frame(
node_to = test_rows,
node_from = test_cols,
stringsAsFactors = FALSE
)
# ---- (B) Subset edges NOT removed
not_removed <- all_edges[
!paste(all_edges$node_to, all_edges$node_from) %in%
paste(removed_edges_idx$node_to, removed_edges_idx$node_from),
]
not_removed$removed <- 0
not_removed$predicted_values <- NA
# ---- (C) Combine removed + not removed
not_removed$k <- k
not_removed$lambda <- lambda
complete_edges <- rbind(
results,     # removed edges with predictions
not_removed  # not-removed edges
)
print("before temp writting!")
write.csv(complete_edges,
file = "temp.csv",
row.names = FALSE)
}
}
### ---- save results for current combination ----
# After finishing the k/lambda loops, append the 'results' to 'combined_results'
# ---- (D) Append to combined_results
complete_edges_all <- read.csv("temp.csv")
combined_results <- rbind(
combined_results,
cbind(
data.frame(
emln_id = emln_id,
train_layer = layers_to_train,
test_layer = layer_to_predict,
prop_ones_removed = prop_ones_to_remove,
amount_of_removed_1 = num_1_to_remove,
prop_zeros_removed = prop_zeros_to_remove,
amount_of_removed_0 = num_0_to_remove
),
complete_edges_all
)
)
#remove("temp.csv")
}
}
getwd()
## ---- functions ----
sigmoid <- function(x) {
1 / (1 + exp(-x))
}
d
## ---- load df ----
d <- read_csv('combined_results_0.2_rem_values_nonbinary_all_edges2.csv')
??read.csv
## ---- load libraries ----
library(tidyverse)
library(ggplot2)
## ---- load df ----
d <- read_csv('combined_results_0.2_rem_values_nonbinary_all_edges2.csv')
??softimpute
## --- Setup ------------------------------------------------------------
set.seed(42)  # for reproducibility
# Names (optional, just for readability)
pollinators <- paste0("Pol", 1:5)
plants      <- paste0("Plant", 1:6)
# Create a 5x6 weighted interaction matrix with strengths 0..9
A <- matrix(sample(0:9, size = length(pollinators) * length(plants), replace = TRUE),
nrow = length(pollinators), ncol = length(plants),
dimnames = list(pollinators, plants))
A
## --- Truncated SVD ----------------------------------------------------
# Choose the truncation rank k (e.g., 2 or 3)
k <- 2
# Full SVD (fine for small matrices). For big matrices, see 'irlba' note below.
svd_full <- svd(A)
U <- svd_full$u               # 5 x 5
d <- svd_full$d               # length 5 (since min(5,6) = 5)
V <- svd_full$v               # 6 x 5
U
# Keep only the top-k singular values/vectors
U_k <- U[, 1:k, drop = FALSE] # 5 x k
d_k <- d[1:k]                 # k
V_k <- V[, 1:k, drop = FALSE] # 6 x k
# Reconstruct the rank-k approximation: A_k = U_k * diag(d_k) * t(V_k)
A_k <- U_k %*% diag(d_k, nrow = k, ncol = k) %*% t(V_k)
A_k
set.seed(42)
mat <- matrix(sample(0:1, 15, replace = TRUE, prob = c(0.6, 0.4)),
nrow = 5, ncol = 3)
rownames(mat) <- paste0("C", 1:5) # Consumers
colnames(mat) <- paste0("P", 1:3) # Products
# 2. Convert to a network object
net <- network(mat, bipartite = TRUE, directed = FALSE)
# an example for visualization of a bipartite network with the R package ggnetwork.
# Load necessary libraries
library(network)
library(ggplot2)
library(ggnetwork)
library(dplyr)
# 1. Create a dummy bipartite adjacency matrix (5 Rows, 3 Columns)
# Rows could be "Consumers", Columns could be "Products"
set.seed(42)
mat <- matrix(sample(0:1, 15, replace = TRUE, prob = c(0.6, 0.4)),
nrow = 5, ncol = 3)
rownames(mat) <- paste0("C", 1:5) # Consumers
colnames(mat) <- paste0("P", 1:3) # Products
# 2. Convert to a network object
net <- network(mat, bipartite = TRUE, directed = FALSE)
# 3. Create the layout manually for the "row-front-of-row" effect
# We extract the number of nodes in each set
num_rows <- 5
num_cols <- 3
# Define coordinates
# Set 1 (Consumers): Y = 1, X spread from 1 to 5
# Set 2 (Products):  Y = 2, X spread from 1 to 3 (centered)
custom_layout <- matrix(0, nrow = 8, ncol = 2)
# Coordinates for Set 1
custom_layout[1:5, 1] <- 1:5        # x
custom_layout[1:5, 2] <- 1          # y
# Coordinates for Set 2 (centered relative to Set 1)
custom_layout[6:8, 1] <- seq(1.5, 4.5, length.out = 3) # x
custom_layout[6:8, 2] <- 2                             # y
# 4. Convert to ggnetwork data frame using the custom layout
df <- ggnetwork(net, layout = custom_layout)
# Add a "Type" column for better styling
df$type <- ifelse(df$y == 1, "Consumer", "Product")
# 5. Plot
ggplot(df, aes(x = x, y = y, xend = xend, yend = yend)) +
geom_edges(color = "grey70", size = 0.8, alpha = 0.6) +
geom_nodes(aes(color = type), size = 10) +
geom_nodetext(aes(label = vertex.names), color = "white", fontface = "bold") +
theme_blank() +
scale_color_manual(values = c("Consumer" = "#2C3E50", "Product" = "#E74C3C")) +
theme(legend.position = "bottom") +
labs(title = "Bipartite Network: Row-in-Front-of-Row",
subtitle = "Visualized with ggnetwork using custom coordinates")
data(memmott1999)
# memmott1999 is a matrix where:
# Rows = Pollinators (79 species)
# Cols = Plants (25 species)
# 2. Convert the matrix to a bipartite network object
# We transpose it if you prefer Plants on one side and Pollinators on the other
net_mat <- as.matrix(memmott1999)
# empirical example: use the memott data from the bipartite package as an example
# Load libraries
library(bipartite)
# 1. Load the data
data(memmott1999)
# 2. Convert the matrix to a bipartite network object
# We transpose it if you prefer Plants on one side and Pollinators on the other
net_mat <- as.matrix(memmott1999)
net <- network(net_mat, bipartite = TRUE, directed = FALSE)
# 3. Define the custom "Row-in-front-of-Row" Layout
# Set 1 (Pollinators): Indices 1 to 79
# Set 2 (Plants): Indices 80 to 104
num_pollinators <- nrow(net_mat)
num_plants <- ncol(net_mat)
custom_layout <- matrix(0, nrow = (num_pollinators + num_plants), ncol = 2)
# Coordinates for Pollinators (Top Row)
custom_layout[1:num_pollinators, 1] <- seq(1, 100, length.out = num_pollinators)
custom_layout[1:num_pollinators, 2] <- 2
# Coordinates for Plants (Bottom Row)
# We center the plants relative to the pollinator row for better aesthetics
custom_layout[(num_pollinators + 1):(num_pollinators + num_plants), 1] <- seq(10, 90, length.out = num_plants)
custom_layout[(num_pollinators + 1):(num_pollinators + num_plants), 2] <- 1
# 4. Convert to ggnetwork
df <- ggnetwork(net, layout = custom_layout)
# 5. Add metadata for styling
df$type <- ifelse(df$y == 2, "Pollinator", "Plant")
# 6. Plotting
ggplot(df, aes(x = x, y = y, xend = xend, yend = yend)) +
# Use thin lines and low alpha because there are many edges (415)
geom_edges(color = "grey80", size = 0.2, alpha = 0.4) +
# Points for species
geom_nodes(aes(color = type), size = 2) +
theme_blank() +
scale_color_manual(values = c("Pollinator" = "#3498DB", "Plant" = "#27AE60")) +
labs(title = "Memmott (1999) Bipartite Network",
subtitle = "Top Row: 79 Pollinators | Bottom Row: 25 Plants",
color = "Species Group") +
theme(legend.position = "bottom")
# 2. Add Random Edge Weights (1 to 10) to the matrix
# We only add weights where an interaction actually exists (non-zero)
set.seed(42)
weights <- mat_weight <- net_mat
weights[weights > 0] <- sample(1:10, sum(net_mat > 0), replace = TRUE)
# 3. Create the network object
# 'ignore.eval = FALSE' and 'names.eval' are required to keep the edge weights
net <- network(weights, bipartite = TRUE, directed = FALSE,
matrix.type = "bipartite", ignore.eval = FALSE, names.eval = "weight")
# 4. Generate and Add Random Abundances (Vertex Attributes)
# Total vertices = 79 pollinators + 25 plants = 104
num_nodes <- network.size(net)
abundances <- sample(10:100, num_nodes, replace = TRUE)
set.vertex.attribute(net, "abundance", abundances)
# 5. Define the "Row-in-front-of-Row" Layout
num_poll <- nrow(net_mat)
num_plan <- ncol(net_mat)
custom_layout <- matrix(0, nrow = num_nodes, ncol = 2)
# Pollinators (Top)
custom_layout[1:num_poll, 1] <- seq(1, 100, length.out = num_poll)
custom_layout[1:num_poll, 2] <- 2
# Plants (Bottom - Centered)
custom_layout[(num_poll + 1):num_nodes, 1] <- seq(15, 85, length.out = num_plan)
custom_layout[(num_poll + 1):num_nodes, 2] <- 1
# 6. Convert to ggnetwork
# We must specify 'weights = "weight"' to ensure the attribute is carried over
df <- ggnetwork(net, layout = custom_layout, weights = "weight")
# Add grouping for coloring
df$type <- ifelse(df$y == 2, "Pollinator", "Plant")
# 7. Final Plot
ggplot(df, aes(x = x, y = y, xend = xend, yend = yend)) +
# Map edge thickness (size) to 'weight'
geom_edges(aes(size = weight), color = "grey75", alpha = 0.3) +
# Map node size to 'abundance'
geom_nodes(aes(color = type, size = abundance)) +
geom_nodetext(aes(label = vertex.names), size = 1.5, repel = TRUE, max.overlaps = 10) +
scale_size_continuous(range = c(0.5, 4)) + # Control the spread of circle sizes
scale_color_manual(values = c("Pollinator" = "#2980B9", "Plant" = "#27AE60")) +
theme_blank() +
labs(title = "Memmott (1999) with Simulated Data",
subtitle = "Node size = Abundance | Edge width = Interaction Strength",
size = "Scale", color = "Group") +
theme(legend.position = "right")
# 1. "Standing" Projection Function
# x_data -> depth (skewed)
# y_data -> height (preserved)
# offset -> shift along the 'depth' axis
project_standing <- function(x, y, offset) {
# We use a 45-degree skew for depth.
# X is skewed by Y slightly to give perspective, then shifted by offset.
list(x = (x * 0.5) + offset, y = y + (x * 0.2))
}
set.seed(101)
offsets <- c(0, 1.5, 3) # The "Z" distance between plates
layer_names <- c("Time 1", "Time 2", "Time 3")
n_nodes <- 50
all_edges <- data.frame()
plane_data <- data.frame()
node_coords_list <- list()
# 2. Generate and Project "Standing" Layers
for (i in 1:3) {
adj <- matrix(sample(0:1, n_nodes^2, replace = TRUE, prob = c(0.96, 0.04)), n_nodes, n_nodes)
net <- network(adj, directed = FALSE)
# Using 'fruchtermanreingold' inside the "plates" for a more organic look
# since circle layout can look a bit rigid when tilted this way.
df <- ggnetwork(net, layout = "fruchtermanreingold")
# Apply the standing projection
proj_start <- project_standing(df$x, df$y, offsets[i])
proj_end   <- project_standing(df$xend, df$yend, offsets[i])
df$x <- proj_start$x; df$y <- proj_start$y
df$xend <- proj_end$x; df$yend <- proj_end$y
df$layer <- layer_names[i]
all_edges <- rbind(all_edges, df)
node_coords_list[[i]] <- df %>% select(x, y, vertex.names) %>% distinct()
# Standing Plane Background (The "Plate")
corners <- data.frame(px = c(0, 0, 1, 1), py = c(0, 1, 1, 0))
proj_c <- project_standing(corners$px, corners$py, offsets[i])
plane_data <- rbind(plane_data, data.frame(x = proj_c$x, y = proj_c$y, layer = layer_names[i]))
}
# 3. Generate 30 Random Inter-layer Links
interlink_indices <- data.frame(
from_node = sample(1:n_nodes, 30, replace = TRUE),
to_node = sample(1:n_nodes, 30, replace = TRUE),
connection = c(rep("1-2", 15), rep("2-3", 15))
)
interlayer_segments <- data.frame()
for (j in 1:nrow(interlink_indices)) {
idx <- interlink_indices[j, ]
l1 <- if(idx$connection == "1-2") 1 else 2
l2 <- if(idx$connection == "1-2") 2 else 3
interlayer_segments <- rbind(interlayer_segments, data.frame(
x = node_coords_list[[l1]]$x[idx$from_node],
y = node_coords_list[[l1]]$y[idx$from_node],
xend = node_coords_list[[l2]]$x[idx$to_node],
yend = node_coords_list[[l2]]$y[idx$to_node]
))
}
# 4. Final Visualization
ggplot() +
# Background Planes (Standing Plates)
geom_polygon(data = plane_data, aes(x = x, y = y, group = layer),
fill = "white", alpha = 0.3, color = "grey70") +
# INTER-layer links (Connecting across time/depth)
geom_segment(data = interlayer_segments, aes(x = x, y = y, xend = xend, yend = yend),
color = "steelblue", alpha = 0.4, size = 0.3, linetype = "dashed") +
# INTRA-layer edges
geom_edges(data = all_edges, aes(x = x, y = y, xend = xend, yend = yend),
color = "grey50", alpha = 0.2, size = 0.15) +
# Nodes
geom_nodes(data = all_edges, aes(x = x, y = y, color = layer), size = 1.3) +
scale_color_viridis_d(option = "plasma", end = 0.8) +
theme_void() +
coord_fixed(ratio = 0.8) + # Adjusted ratio to enhance the standing look
theme(legend.position = "none") +
annotate("text", x = offsets + 0.25, y = 1.2, label = layer_names, fontface = "bold")
setwd("~/Documents/github/figaro-")