myprojects-hg-reborn

Fork

rm(list=ls())

library(tidyverse)
library(igraph)
library(viridis)
library(furrr)

source("/home/lorenzo/myprojects-hg/R-codes/stat_lib.R")

############################################################################
############################################################################
############################################################################
############################################################################
############################################################################

golden_ratio <- golden()

generate_data <- 1

cutoff <- 2.001

n_cores <- 10

k_min <- 3

kf_val <- 1.3
df_val <- 1.6



plan(multicore(workers=return_cores(n_cores)))



if (generate_data==1) {

file_path <- "./"



dat_files2 <- extract_file_list(file_path,".dat", full_path=1)


graph_list <- future_map2(dat_files2,cutoff, function(x,y) get_network_from_table(x,y))

print("I have the networks")

saveRDS(graph_list, "graph_list.RDS")


frag_dist <- future_map(graph_list, fragment_cluster)

print("I finished the fragmentation")

saveRDS(frag_dist, "frag_dist.RDS")

branch_size_dist <- future_map2(graph_list,k_min, function(x,y) remove_nodes_high_deg(x,y) )

print("I have the branch distribution")

saveRDS(branch_size_dist, "branch_size_dist.RDS")

## See https://kutt.it/5ejnu4

### or using the tilde notation
## tv <- map2(graph_list,3, ~ remove_nodes_high_deg(.x,.y) )


small_frag_exp <- future_map(frag_dist, function(x) mean(x$n2)) %>% flatten_dbl

saveRDS(small_frag_exp,"small_fragment_distribution.RDS")

## ## or using the tilde notation
## tt <- map(frag_dist, ~ mean(.x$n2)) %>% flatten_dbl

number_mon <- future_map(graph_list, vcount)  %>% flatten_int
saveRDS(number_mon, "aggregate_sizes.RDS")

diam_list <- future_map(graph_list, diameter) %>% flatten_dbl

saveRDS(diam_list, "diameters_list.RDS")

## var_branch <- map(branch_size_dist,function(x) l2_norm(x$value)) %>% flatten_dbl

} else{

    graph_list <- readRDS("graph_list.RDS")
    frag_dist <- readRDS("frag_dist.RDS")
    branch_size_dist <- readRDS("branch_size_dist.RDS")
    small_frag_exp <- readRDS("small_fragment_distribution.RDS")
    diam_list <- readRDS("diameters_list.RDS")
    number_mon <- readRDS("aggregate_sizes.RDS")
    


}


mm <- lm(small_frag_exp~diam_list)

real <- tibble(d=diam_list, s=small_frag_exp, type="numerics", df=df_val,
               kf=kf_val, n=number_mon)
simu <- tibble(d=diam_list, s=predict(mm), type="linear fit", df=df_val,
               kf=kf_val, n=number_mon)


my_pal <- viridis(3)[1:2]

df_plot <- bind_rows(real, simu) 

df_name <- as.character(df_val) %>%
       remove_special_characters    

    kf_name <- kf_val %>%
        as.character %>% 
       remove_special_characters    

fname <- paste("data_diameter_fragmentation_df_",df_name,
                   "_kf_", kf_name, ".RDS", sep="")

print("fname is, ")
print(fname)

saveRDS(real, fname)

gpl <- ggplot(df_plot, aes(x=d, y=s ,
                           color=type,
                          shape=type, linetype=type
                          )) +
    
geom_point(size=3, stroke=2) + 
    geom_line(size=1.2)+
## facet_wrap( ~ partner, nrow = 2, scales = "free_y"  )+
    my_ggplot_theme2("right")+
    ## theme(panel.spacing = unit(2, "lines"))+


 scale_color_manual(NULL, ## labels=c("Export Share","Import Share" ),
                    values=my_pal)+    

 scale_shape_manual(NULL, values=c(NA, 1)
                      )+ 

    scale_linetype_manual(NULL, values=c(1,NA)
                         )+ 

##     coord_cartesian(ylim = c(0, .5)) +

## scale_y_continuous(breaks=pretty_breaks(n=5),labels = mypercent)+
## scale_x_continuous(breaks=seq(1994, 2018, by=2))+

    
## labs(title="Share of World Imports and Exports under FTAs")+

xlab("Diameter")+
ylab("Small fragment size")



fname <- paste("small_fragment_vs_diameter.pdf")
ggsave(fname, gpl, width=7*golden_ratio,height=7)


print("So far so good")