first version of building the network

R/build_network.R

@@ -0,0 +1,31 @@
+source(here::here("R", "graphon_distribution.R"))
+
+compute_adj_matrix <- function(
+  X_matrix,
+  v,
+  a,
+  phi,
+  rho_n,
+  Fv = pnorm
+) {
+  xi <- vapply(v, FUN = function(x){pgraphon(x, a, Fv, X_matrix)}, numeric(1))
+  print(xi)
+  dim_adj <- length(v)
+  
+  rand_mat <- matrix(runif(dim_adj^2), dim_adj, dim_adj)
+  upper_mask <- upper.tri(rand_mat, diag=FALSE)
+  rand_mat[!upper_mask] <- 0
+  rand_mat <- rand_mat + t(rand_mat) + diag(nrow=dim_adj)
+  graphon_function <- function(x,y) { rho_n * phi(x,y)}
+  probabilities <- outer(xi, xi, FUN=graphon_function)
+  
+  adj_mat <- (rand_mat < probabilities) * 1L
+}
+
+set.seed(1)
+X <- matrix(rnorm(4), nrow = 4, ncol=1)
+a <- 0.5
+v <- rnorm(4)
+
+adj <- compute_adj_matrix(X, v, a, phi = function(x,y) {x * y}, 0.5)
+adj

scripts/plot_non_normal_X.qmd

@@ -348,3 +348,82 @@ results |>
        colour=latex2exp::TeX("$a$"),
        shape=latex2exp::TeX("$\\alpha$"))
 ```
+```{r k = n^alpha data generation, N(0,1)}
+#| cache: true
+#| echo: false
+#| collapse: true
+ns <- seq(100, 5000, 100)
+as <- seq(0, 20, 2)
+alphas <-  seq(0.1, 0.5, 0.1)
+
+set.seed(100)
+results <- data.frame(dim_n = integer(),
+                      dim_k = integer(),
+                      param_a = double(),
+                      param_alpha = double(),
+                      ssv = double())
+for (a in as) {
+  for (i in 1:length(ns)) {
+    for (j in 1:length(alphas)) {
+      n <- ns[i]
+      # HERE WE USE THE CEILING AND NOT FLOOR!
+      K <- ceiling(n^alphas[j])
+      if (!K > 0) next # skip if K is equal to zero
+      # use the default seed 1L
+      Q <- compute_matrix(seed=1L,
+                     a= a,
+                     n = n,
+                     K = K,
+                     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)
+      
+      ssv <- compute_minmax_sv(Q)[["smallest_singular_value"]]
+      
+      current_res <- data.frame(dim_n = n, dim_k = K, param_a = a, param_alpha=alphas[j], ssv =ssv)
+      results <- rbind(results, current_res)
+    }
+  }
+}
+```
+
+```{r k = n^alpha plotting, U[0,2]}
+results |>
+  filter(param_a %in% c(0, 10, 20)) |>
+  mutate(param_a = as.factor(param_a),
+         param_alpha = as.factor(param_alpha)) |>
+  group_by(param_a, param_alpha) |>
+  ggplot(aes(dim_n, ssv * dim_k, col=param_a, shape=param_alpha, interaction(param_a, param_alpha))) +
+  geom_point(size=1.5) +
+  geom_line() +
+  geom_function(fun = function(x) {x^(0.5)}, colour="black") +
+  #scale_y_log10() +
+  theme_bw() +
+  labs(x=latex2exp::TeX("$n$"),
+       y=latex2exp::TeX("Smallest singular value of $Q$"),
+       title=latex2exp::TeX("Smallest singular value of $Q$ with respect to $a$."),
+       subtitle = latex2exp::TeX(("Hyperparameter $k = n^{\\alpha}$. Black line is $\\sqrt{n}$, and $X \\sim N(0,1) $")),
+       colour=latex2exp::TeX("$a$"),
+       shape=latex2exp::TeX("$\\alpha$"))
+```
+
+```{r k = n^alpha plotting, U[0,2]}
+results |>
+  filter(param_a %in% c(0, 10, 20)) |>
+  mutate(param_a = as.factor(param_a),
+         param_alpha = as.factor(param_alpha)) |>
+  group_by(param_a, param_alpha) |>
+  ggplot(aes(dim_n, ssv / sqrt(dim_n) * dim_k, col=param_a, shape=param_alpha, interaction(param_a, param_alpha))) +
+  geom_point(size=1.5) +
+  geom_line() +
+ # geom_function(fun = function(x) {x^(0.5)}, colour="black") +
+  #scale_y_log10() +
+  theme_bw() +
+  labs(x=latex2exp::TeX("$n$"),
+       y=latex2exp::TeX("Smallest singular value of $Q$ / sqrt(n)"),
+       title=latex2exp::TeX("Smallest singular value of $Q$ with respect to $a$."),
+       subtitle = latex2exp::TeX(("Hyperparameter $k = n^{\\alpha}$. Black line is $\\sqrt{n}$, and $X \\sim N(0,1) $")),
+       colour=latex2exp::TeX("$a$"),
+       shape=latex2exp::TeX("$\\alpha$"))
+```