# load local files
source(here::here("R", "singular_values.R"))
source(here::here("R", "graphon_distribution.R"))
source(here::here("R","singular_value_plot.R"))

# load libaries for data handling
library(ggplot2)
library(tidyr)
library(dplyr)

# Line plots -------------------------------------------------------------------
# Create a grid of a‑values
a_grid <- seq(-20, 20, length.out = 200)

# function which takes only a to compute Q_c
make_matrix <- function(a) { compute_matrix(seed=4L,
                              a= a,
                              n = 2,
                              K = 2,
                              sample_X_fn = function(n) {matrix(rnorm(n), ncol = 1L)},
                              fv = function(x) {dnorm(x, mean=0, sd=1)},
                              Fv = function(x) {pnorm(x, mean=0, sd=1)},
                              guard = 1e-12)}

# Compute the matrices and reshape to long format
df_entries <- tibble(a = a_grid) %>%
  mutate(
    M = purrr::map(a, make_matrix),                     # list‑column of matrices
    m11 = purrr::map_dbl(M, ~ .x[1, 1]),
    m12 = purrr::map_dbl(M, ~ .x[1, 2]),
    m21 = purrr::map_dbl(M, ~ .x[2, 1]),
    m22 = purrr::map_dbl(M, ~ .x[2, 2])
  ) %>%
  select(a, m11, m12, m21, m22) %>%
  pivot_longer(-a, names_to = "entry", values_to = "value")

# Plot
ggplot(df_entries, aes(x = a, y = value, colour = entry, linetype = entry)) +
  geom_line(linewidth = 1) +
  labs(
    title = "Matrix entries as a function of the parameter `a`",
    x = "a",
    y = "Matrix entry value",
    colour = "Entry"
  ) +
  theme_minimal()

# Heat map for a single larger matrix ------------------------------------------

# Choose a value of a
a0 <- -10
M0 <-  compute_matrix(seed=1L,
                      a= a0,
                      n = 50,
                      K = 50,
                      sample_X_fn = function(n) {matrix(rnorm(n), ncol = 1L)},
                      fv = function(x) {dnorm(x, mean=0, sd=1)},
                      Fv = function(x) {pnorm(x, mean=0, sd=1)},
                      guard = 1e-12)

# Convert to a tidy data frame for ggplot
df_heat <- as.data.frame(M0) %>%
  mutate(row = row_number()) %>%
  pivot_longer(-row, names_to = "col", values_to = "value") %>%
  mutate(
    col = as.integer(gsub("V", "", col)),   # turn "V1","V2" into 1,2
    row = factor(row, levels = rev(unique(row)))   # reverse y‑axis for proper orientation
  )

ggplot(df_heat, aes(x = col, y = row, fill = value)) +
  geom_tile(colour = "grey30", size = 0.5) +
  scale_fill_gradient2(
    low = "steelblue", mid = "white", high = "tomato",
    midpoint = 0.5, limits = c(min(df_heat$value), max(df_heat$value))
  ) +
  labs(
    title = sprintf("Heat‑map of the matrix at a = %.2f", a0),
    x = "Column", y = "Row", fill = "Value"
  ) +
  theme_minimal() +
  theme(axis.text = element_blank(),
        axis.ticks = element_blank())