In AnthonyEbert/StartNetwork: Calculates KL divergence between mechanistic and statistical network models

knitr::opts_chunk$set(cache = TRUE, echo = FALSE, message = FALSE, warning = FALSE)

We are interested in finding the closest mechanistic model $q(\cdot|\lambda)$, parameterised by $\lambda$ to a known statistical model $p(\cdot|\theta)$, where $\theta$ is known.

$$ \begin{aligned} KL[q(\cdot)||p(\cdot|\theta)] &= \sum_{y \in \mathcal{Y}} q(y) \log \frac{ q(y) }{p(y|\theta)} \ &= \log[z(\theta)] - \sum_{y \in \mathcal{Y}} q(y) \left[ \sum_i \theta_i s_i(y) \right] + \sum_{y \in \mathcal{Y}} q(y) \log q(y) \end{aligned} $$

The relative importance of the likelihood and entropy

set.seed(1)

library(StartNetwork)

library(parallel)

cl <- makeCluster(detectCores())

param_range_large <- rep(seq(0.025, 0.975, by = 0.025), 20)

x1 <- parSapply(cl, param_range_large, er_KL, nl = 5, pl = 0.3, include_entropy = TRUE, replicates = 1000)
x2 <- parSapply(cl, param_range_large, er_graphs_KL, nl = 5, pl = 0.3, include_entropy = TRUE, replicates = 1000)

library(ggplot2)


df <- data.frame(x = c(x1, x2), parameter = rep(param_range_large, 2), summary_statistic = c(rep("a", length(param_range_large)), rep("b", length(param_range_large))))

ggplot(df) + aes(x = parameter, y = x, group = parameter, col = summary_statistic) + geom_boxplot()

replicates <- 10

G_4 <- replicate(replicates,er_graphs_KL(0.3, pl = 0.3, nl = 4, replicates = 1000))
G_8 <- replicate(replicates,er_graphs_KL(0.3, pl = 0.3, nl = 8, replicates = 1000))
G_12 <- replicate(replicates,er_graphs_KL(0.3, pl = 0.3, nl = 12, replicates = 1000))
G_16 <- replicate(replicates,er_graphs_KL(0.3, pl = 0.3, nl = 16, replicates = 1000))

G <- data.frame(value = c(G_4, G_8, G_12, G_15), nodes = rep(c(4, 8, 12, 16), replicates), space = "G")


S2_4 <- replicate(replicates,er_dd_KL(0.3, pl = 0.3, nl = 4, replicates = 1000, sorted = FALSE))
S2_8 <- replicate(replicates,er_dd_KL(0.3, pl = 0.3, nl = 8, replicates = 1000, sorted = FALSE))
S2_12 <- replicate(replicates,er_dd_KL(0.3, pl = 0.3, nl = 12, replicates = 1000, sorted = FALSE))
S2_16 <- replicate(replicates,er_dd_KL(0.3, pl = 0.3, nl = 16, replicates = 1000, sorted = FALSE))

S2 <- data.frame(value = c(S2_4, S2_8, S2_12, S2_16), nodes = rep(c(4, 8, 12, 16), replicates), space = "S2")

S3_4 <- replicate(replicates,er_dd_KL(0.3, pl = 0.3, nl = 4, replicates = 1000, sorted = TRUE))
S3_8 <- replicate(replicates,er_dd_KL(0.3, pl = 0.3, nl = 8, replicates = 1000, sorted = TRUE))
S3_12 <- replicate(replicates,er_dd_KL(0.3, pl = 0.3, nl = 12, replicates = 1000, sorted = TRUE))
S3_16 <- replicate(replicates,er_dd_KL(0.3, pl = 0.3, nl = 16, replicates = 1000, sorted = TRUE))

S3 <- data.frame(value = c(S3_4, S3_8, S3_12, S3_16), nodes = rep(c(4, 8, 12, 16), replicates), space = "S3")

df <- dplyr::bind_rows(G, S2, S3)

ggplot(df) + aes(x = factor(nodes), y = value, col = space) + geom_boxplot() 

param_range_large <- rep(seq(0.2, 0.4, by = 0.05), 10)

g <- parSapply(cl, param_range_large, er_graphs_KL, nl = 8, pl = 0.3, include_entropy = TRUE, replicates = 2000)
s_2 <- parSapply(cl, param_range_large, er_dd_KL, nl = 8, pl = 0.3, replicates = 2000, sorted = FALSE)
s_3 <- parSapply(cl, param_range_large, er_dd_KL, nl = 8, pl = 0.3, replicates = 2000, sorted = TRUE)

df <- data.frame(x = c(g, s_2, s_3), parameter = rep(param_range_large, 3), space = c(rep("g", length(param_range_large)), rep("s_2", length(param_range_large)), rep("s_3", length(param_range_large)) ))

ggplot(df) + aes(x = parameter, y = x, group = parameter) + geom_boxplot() + facet_wrap(~space)

param_range_large <- rep(seq(1, 9), 10)
s_3 <- parSapply(cl, param_range_large, ba_dd_KL_m, nl = 15, pl = 1, replicates = 20, sorted = TRUE)
df <- data.frame(x = s_3, parameter = param_range_large)
ggplot(df) + aes(x = parameter, y = x, group = parameter) + geom_boxplot()

ba_dd_KL_m(2, nl = 30, pl = -15)

AnthonyEbert/StartNetwork documentation built on April 24, 2020, 3:28 a.m.