Homework 5

Week 5 Assignment: Brokerage and Power.

#Load the libraries
  library(dplyr)
  library(igraph)
  library(tidyr)
  library(corrr)
  library(ggplot2)

#Read the network into the environment. This will import it as an edgelist
  ##data_orig <- read.csv("C:/Users/Yifan Desktop/OneDrive/Social Network Analysis/COW_Trade_4.0/Dyadic_COW_4.0.csv", header=T,stringsAsFactors=FALSE)
data_orig <- read.csv("C:/Users/LYF-YOGA14s/OneDrive/Social Network Analysis/COW_Trade_4.0/Dyadic_COW_4.0.csv", header=T,stringsAsFactors=FALSE)
  ##flow 1 is Imports of Country A from Country B, in US millions of current dollars
  ##flow 2 is imports of country B from country A

  ##trade.nodes <- read.csv("C:/Users/Yifan Desktop/OneDrive/Social Network Analysis/COW_Trade_4.0/National_COW_4.0.csv", header = T, stringsAsFactors = F)
trade.nodes <- read.csv("C:/Users/LYF-YOGA14s/OneDrive/Social Network Analysis/COW_Trade_4.0/National_COW_4.0.csv", header = T, stringsAsFactors = F)

  trade2014 <- subset(data_orig, data_orig$year == 2014)
  
#Create weight column that holds the total trade on each edge. 
  temp <- trade2014 %>%
    mutate(from = importer1,
           to = importer2,
           weight = flow1) %>%
    select(from, to, weight, year)
  temp2 <- trade2014 %>%
    mutate(from = importer2,
           to = importer1,
           weight = flow2) %>%
    select(from, to, weight, year)

  temp_edgelist <- full_join(temp,temp2)
  network_edgelist <- subset(temp_edgelist, temp_edgelist$weight > 0)
  
  #Create igraph objects from the edgelists
  trade2014.ig <- graph.edgelist(as.matrix(network_edgelist[,c("from","to")]), directed = T)
  
  #Specify the weights values
  E(trade2014.ig)$weight <- network_edgelist$weight 

  #Add in long and short versions of country names for reports
    V(trade2014.ig)$long_name <- trade.nodes$statename[match(V(trade2014.ig)$name, trade.nodes$statename)]
    V(trade2014.ig)$short_name <- trade.nodes$stateabb[match(V(trade2014.ig)$name, trade.nodes$statename)]
    V(trade2014.ig)$ccode <- trade.nodes$ccode[match(V(trade2014.ig)$name, trade.nodes$statename)]
    V(trade2014.ig)$continent <- trade.nodes$continent[match(V(trade2014.ig)$name, trade.nodes$statename)]

  ##rm(data_orig)
    trade2014.nodes <- data.frame(name = V(trade2014.ig)$long_name,
                                  short_name = V(trade2014.ig)$short_name,
                                  ccode = V(trade2014.ig)$ccode)
    trade2014.nodes$continent <- V(trade2014.ig)$continent

#Find network features: 
is_bipartite(trade2014.ig)

[1] FALSE

is_directed(trade2014.ig)

[1] TRUE

is_weighted(trade2014.ig)

[1] TRUE

The original dataset is the trade dataset version 4 from the Correlates of War Project. In this subject I only use trade data of year 2014. The format is edgelist. The nodes are countries, and the ties are the trading relations between countries in 2014. The network is directed and weighted.

Let’s look at some basic descriptive facts.

#Original network size
vcount(trade2014.ig)

[1] 186

ecount(trade2014.ig)

[1] 22451

ecount(trade2014.ig) * 2 / vcount(trade2014.ig)

[1] 241.4086

ecount(trade2014.ig) / (vcount(trade2014.ig) * (vcount(trade2014.ig - 1)))

[1] 0.6524557

#summarize numeric network attribute
summary(E(trade2014.ig)$weight)

    Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
     0.0      0.2      4.9    813.1     78.5 472525.2 

The network in 2014 has 186 nodes, i.e. 186 countries involed in the trading network. There are 22451 egdes. Each country has 241 connected edges on average. On average, one country import 813 million dollars of good from another country each year. 65% of potential ties exit.

Let’s classify all dyads and traids in the network:

#Classify all dyads in the network:
igraph::dyad.census(trade2014.ig)

$mut
[1] 9933

$asym
[1] 2585

$null
[1] 4687

#Classify all triads in the network:
igraph::triad_census(trade2014.ig)

 [1]  75574  62420  92513  11021   9514  11233  60137  45246   5530
[10]    788 188669  14094   9006  15882 147383 306230

75574 / sum(igraph::triad_census(trade2014.ig))

[1] 0.07161783

(75574+62420+92513+11021+9514+11233+60137+45246) / sum(igraph::triad_census(trade2014.ig))

[1] 0.3484117

igraph::components(trade2014.ig)$no

[1] 1

There are 9933 mutual bilateral trade relations, 2585 unilateral trade relations, and 4687 pairs of coutries have no trading in 2014.

7% of the triads are empty, and 65% of them are triangle. The network is quite dense.

Let’s look at the distribution of nodes centrality.

trade2014.nodes$closeness <- centr_clo(trade2014.ig)$res

trade2014.nodes$betweenness <- centr_betw(trade2014.ig, directed = TRUE)$res

trade2014.nodes$eigen <- centr_eigen(trade2014.ig, directed=TRUE, scale=FALSE)$vector

trade2014.nodes$bonpow <- power_centrality(trade2014.ig)


mat2014<-as.matrix(as_adjacency_matrix(trade2014.ig, attr="weight"))

#square the adjacency matrix
mat2014sq<-t(mat2014) %*% mat2014

#Calculate the proportion of reflected centrality.
trade2014.nodes$rc<-diag(mat2014sq)/rowSums(mat2014sq)
#replace missing values with 0
trade2014.nodes$rc<-ifelse(is.nan(trade2014.nodes$rc),0,trade2014.nodes$rc)

#Calculate received eigenvalue centrality
trade2014.nodes$eigen.rc<-trade2014.nodes$eigen*trade2014.nodes$rc

#Calculate the proportion of derived centrality.
  trade2014.nodes$dc<-1-diag(mat2014sq)/rowSums(mat2014sq)
#replace missing values with 0
  trade2014.nodes$dc<-ifelse(is.nan(trade2014.nodes$dc),1,trade2014.nodes$dc)
#Calculate received eigenvalue centrality
trade2014.nodes$eigen.dc<-trade2014.nodes$eigen*trade2014.nodes$dc

---

Following last week’s analysis, I continued to add brokerage.

To calculate brokerage scores with the trade dataset, I try to create a statnet network (whose import script is not provided.) I’m not sure what I did wrong. The betweenness is much different from what I expect (and very different from results igraph). Then I check the degree. The statnet gives a totally different outcome from igraph too. The two packages shouldn’t give different results for such a simple parameter.

I also failed to calculate brokerage scores. I only got the total scores while all roles scores are NA.

My codes show below. What did I do wrongly?

##I need the statnet package and network data to calculate brokerage.
library(statnet)
trade2014.stat <- network(as.matrix(network_edgelist[1:2], matrix.type = 'edgelist'))

print(trade2014.stat)

 Network attributes:
  vertices = 186 
  directed = TRUE 
  hyper = FALSE 
  loops = FALSE 
  multiple = FALSE 
  bipartite = FALSE 
  total edges= 22451 
    missing edges= 0 
    non-missing edges= 22451 

 Vertex attribute names: 
    vertex.names 

 Edge attribute names not shown 

trade2014.stat%e%"weight"      <- network_edgelist$weight
trade2014.stat%e%"year"      <- network_edgelist$year

trade2014.stat%v%"long_name" <- trade.nodes$statename[match(trade2014.stat%v%"vertex.names", trade.nodes$statename)]
trade2014.stat%v%"short_name" <- trade.nodes$stateabb[match(trade2014.stat%v%"vertex.names", trade.nodes$statename)]
trade2014.stat%v%"ccode" <- trade.nodes$ccode[match(trade2014.stat%v%"vertex.names", trade.nodes$statename)]
trade2014.stat%v%"continent" <- trade.nodes$continent[match(trade2014.stat%v%"vertex.names", trade.nodes$statename)]

network::list.vertex.attributes(trade2014.stat)

[1] "ccode"        "continent"    "long_name"    "na"          
[5] "short_name"   "vertex.names"

trade2014.nodes$between_stat <- sna::betweenness(trade2014.stat)
trade2014.nodes$degree <- igraph::degree(trade2014.ig)
trade2014.nodes$degree_stat <- sna::degree(trade2014.stat)

brokerage_stat <- data.frame(brokerage(trade2014.stat, cl = trade2014.nodes$continent)$z.nli)
trade2014.nodes<-trade2014.nodes %>%
  mutate(broker.tot = brokerage_stat$t,
         broker.coord = brokerage_stat$w_I,
         broker.itin = brokerage_stat$w_O,
         broker.rep = brokerage_stat$b_IO,
         broker.gate = brokerage_stat$b_OI,
         broker.lia = brokerage_stat$b_O)