rm(list = ls())
library("mcPAFit")
alpha <- 1
net <- simple_net(time_step = 10000,
p_const = 0.01, alpha = alpha)
stats <- get_my_statistics(net$net)
result_perfect <- empirical_estimate_perfect(net$net)
result <- empirical_estimate(stats$final_deg)
final_deg <- stats$final_deg
result_square <- empirical_estimate_square(stats$final_deg)
sqrt(result$variance) / result$A
y_min <- min(result$A,result_square$A, result_perfect$A[as.character(result$k)], result_square$A_1,
result_square$A_2, result_square$A_3, result$lower_A)
y_max <- max(result$A,result_square$A, result_perfect$A[as.character(result$k)], result_square$A_1,
result_square$A_2, result_square$A_3, result$upper_A)
plot(result$k, result$A, log = "xy", pch = 20,
#ylim = c(y_min,y_max),
xlab = "Degree", ylab = "Attachment rate")
polygon(c(result$k,rev(result$k)) ,
c(result$lower_A[as.character(result$k)],
rev(result$upper_A[as.character(result$k)])),
col = rgb(1,0,0,0.1),border = NA)
lines(result$k,result$k^alpha)
points(result_square$k[as.character(result$k)], result_square$A[as.character(result$k)], pch = 20,col = "blue")
points(result_perfect$k[as.character(result$k)],
result_perfect$A[as.character(result$k)], pch = 20,col = "red")
points(result_square$k, result_square$A_1, pch = 20, col = "brown")
points(result_square$k, result_square$A_2, pch = 20, col = "yellow")
points(result_square$k, result_square$A_3, pch = 20, col = "green")
points(result_square$k, result_square$A_4, pch = 20, col = "purple")
points(result_square$k, result_square$A_5, pch = 20, col = "pink")
deg_vec <- table(stats$final_deg)
degree <- as.integer(names(deg_vec))
k_1 <- quantile(stats$final_deg,95/100)
abline(v = k_1)
result$alpha
result$ci
result_perfect$alpha
result_perfect$ci
#alpha_timeline
N_k_t <- result$N_k_t
N_k_t <- N_k_t / N_k_t[1] * result_square$W_k[1]
y_min <- min(c(deg_vec,result_perfect$T_k, N_k_t))
y_max <- max(c(deg_vec,result_perfect$T_k,N_k_t))
plot(names(result_perfect$T_k),result_perfect$T_k, pch = 20,
ylim = c(y_min, y_max), col = "black", log = "xy",
xlab = "Degree", ylab = "Num. of appearances")
points(names(result_square$W_k), result_square$W_k, pch = 20, col = "red")
points(names(N_k_t),N_k_t, col = "blue")
points(as.integer(names(deg_vec[-1])), deg_vec[-1],
col = "green")
#points(degree,deg_vec/PA, col = "green")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.