library(igraph)
library(ineq)
# Random model
g1 = sample_gnp(n=100, p=5/100)
coords = layout_with_fr(g1)
plot(g1, layout=coords, vertex.size = 3, vertex.label=NA)

# Preferential attachment model
g2 = sample_pa(n=100, m=3, directed=FALSE)
coords = layout_with_fr(g2)
plot(g2, layout=coords, vertex.size = 3, vertex.label=NA)

# random model
d1 = degree(g1)
dist = table(d1)
barplot(dist / sum(dist), xlab="Degree", ylab="Frequency");

# mean and median
summary(d1)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1.00    3.00    4.50    4.64    6.00   10.00
# skewness
skewness = function(x) mean( ((x - mean(x)) / sd(x))^3 )
skewness(d1)
## [1] 0.4666584
# Lorenz curve and Gini coefficient
Gini(d1)
## [1] 0.2375862
lc = Lc(d1, plot=TRUE)

# preferential attachment model
d2 = degree(g2)
dist = table(d2)
barplot(dist / sum(dist), xlab="Degree", ylab="Frequency");

# mean and median
summary(d2)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    3.00    3.00    4.00    5.88    6.00   23.00
# skewness
skewness(d2)
## [1] 2.248101
# Lorenz curve and Gini coefficient
Gini(d2)
## [1] 0.330034
lc = Lc(d2, plot=TRUE)

# Cumulative distribution
ccdf = function(d) {
  n = length(d)
  max = max(d)
  P = rep(0, max);
  for (i in 1:length(P)) {
    P[i] = length(d[d >= i]) / n;
  } 
  return(P)
}

# random model
P = ccdf(d1)
plot(1:max(d1), P, log="xy", type = "l", xlab="Degree", ylab="CCDF");

# preferential attachment model
P = ccdf(d2)
plot(1:max(d2), P, log="xy", type = "l", xlab="Degree", ylab="CCDF");