Step 0: Load the libraries

library("network")
library("sna")
library("RColorBrewer")
library("ggplot2")
library("GGally")
# Must load other packages first
library("sand")
library("intergraph")
library('scales')
library('mcclust')
library('ergm')
source('inv-logit.R')
source('mycugtest.R')
source('myqaptest.R')
source("lab5-utils.R")
library('ergm')
library('latticeExtra')

Loading bipartite product to user network

Samples created using Snowball approach

#mac
setwd("~/Documents/Master degree/Depaul/2016 Spring term/CSC 495 Social Network Analysis/Final Project/UxU/Final")
# windows
#setwd("~/Master degree/Depaul/2016 Spring term/CSC 495 Social Network Analysis/Final Project/UxU/Final")
#
#edge.dfB <- read.csv("Snowball2-Edges.csv", header = T)
edge.dfB <- read.csv("Snowball3-Edges.csv", header = T)
node.dfB <- read.csv("Snowball2-Nodes.csv", header = T)
AmazonNet <- graph.data.frame(edge.dfB,vertices = node.dfB, directed = F)

summary(AmazonNet)
## IGRAPH UN-B 3147 3602 -- 
## + attr: name (v/c), Rating (v/n), type (v/c), ERating (e/n)
is.bipartite(AmazonNet)
## [1] TRUE
u2u <- bipartite.projection(AmazonNet, which = "true")
## Warning in bipartite.projection(AmazonNet, which = "true"): vertex types
## converted to logical
summary(u2u)
## IGRAPH UNW- 3139 1266475 -- 
## + attr: name (v/c), Rating (v/n), weight (e/n)

Centrality

bet <- betweenness(u2u, directed = F, normalized = T)
clo <- closeness(u2u, mode = 'in' ,normalized = T)
save(bet,file = 'bet.Rdata')
save(clo,file = 'clo.Rdata')
deg <- degree(u2u, mode = 'in', normalized = T)
load('bet.Rdata')
load('clo.Rdata')
eig <- eigen_centrality(u2u,directed = F)
prk <- page.rank(u2u,directed = F) 
rank <- V(u2u)$Rating

u2u.centData <- data.frame(deg=deg,bet=bet,clo=clo,eig=eig$vector,prk=prk$vector, rating=rank)
head(u2u.centData)
##                      deg          bet       clo       eig          prk
## A1020G50MVAOLA 0.2660931 5.293107e-07 0.3581374 0.2123588 0.0003105653
## A103DZF5M6Y6C0 0.2660931 5.293107e-07 0.3581374 0.2123588 0.0003105653
## A103M039GIFX02 0.2989165 2.341633e-07 0.3623975 0.4953985 0.0003215187
## A1058L4DCEX548 0.2989165 2.341633e-07 0.3623975 0.4953985 0.0003215187
## A105FAAQ06FTWA 0.5079669 2.814998e-03 0.3901529 0.6247165 0.0005639858
## A10ASLX7DTTB6Z 0.2989165 2.341633e-07 0.3623975 0.4953985 0.0003215187
##                rating
## A1020G50MVAOLA      1
## A103DZF5M6Y6C0      5
## A103M039GIFX02      5
## A1058L4DCEX548      5
## A105FAAQ06FTWA      3
## A10ASLX7DTTB6Z      5
nrow(u2u.centData)
## [1] 3139
ggcorr(data=u2u.centData[,-6], label = T)

Setting indegree and betweenness attributes

j = 1
for (i in V(u2u)) {
  u2u <- set.vertex.attribute(u2u,'inDegree',i,deg[j])    
  j = j + 1
}
#save(g,file = 'AmazonNet.Rdata')
summary(u2u)

j = 1
for (i in V(u2u)) {
  u2u <- set.vertex.attribute(u2u,'betweenes',i,bet[j])    
  j = j + 1
}
save(u2u,file = 'u2u.Rdata')
summary(u2u)

Degree distribution

load('u2u.Rdata')
u2u.deg <- graph.strength(u2u)
ma <- max(u2u.deg)
mi <- min(u2u.deg)
p <- ggplot(data = data.frame(x=u2u.deg),aes(x=x))
p1 = p + geom_histogram(breaks =seq(mi,ma, by=200), fill='azure2',col='blue')
p2 <- p1 + scale_x_continuous('Weighted degree', breaks=seq(mi,ma,by=200))
p3 <- p2 + ggtitle('Weighted Degree Distribution of User to User Projection')
print(p3)

wts <-V(u2u)$betweenes
min(wts)
## [1] 0
max(wts)
## [1] 0.03387576
p <- ggplot(data=data.frame(x=wts), aes(x=x))
p <- p + geom_bar() #stat = 'count'
p <- p + scale_x_continuous('Weight', breaks=seq(0,.035, by=.005))
p <- p + ggtitle('Edge distribution for the network')
print(p)
## stat_bin: binwidth defaulted to range/30. Use 'binwidth = x' to adjust this.

max(bet)
## [1] 0.03387576
min(bet)
## [1] 0
p <- ggplot(data = data.frame(x=bet), aes(x=x))
p <- p + geom_histogram(breaks=seq(0,0.04, by= 0.005), fill='azure2', col='blue')
p <- p + scale_x_continuous('Betweenness', breaks=seq(0,0.04, by= 0.005))
p <- p + ggtitle('Betweenness for User to User projection')
print(p)

Rating distribution

wts <-V(u2u)$Rating
p <- ggplot(data=data.frame(x=wts), aes(x=x))
p <- p + geom_bar(stat='bin', binwidth=1, drop=T,fill='azure2', col='blue')#geom_bar(stat = 'count')
p <- p + scale_x_continuous('Rating', breaks=seq(0,5, by=1))
p <- p + ggtitle('Rating distribution fo sample network')
print(p)

Distribution of weight edge

wts <-E(u2u)$weight
min(wts)
## [1] 1
max(wts)
## [1] 6
p <- ggplot(data=data.frame(x=wts), aes(x=x))
p <- p + geom_histogram(breaks=seq(1,6, by= 1), fill='azure2', col='blue')#geom_bar(stat = 'bin',breaks=seq(0,5, by=1))
p <- p + scale_x_continuous('Weight', breaks=seq(1,6, by=1))
p <- p + ggtitle('Edge distribution for user to user projection')
print(p)

deg <- graph.strength(u2u)

max(deg)
## [1] 3594
min(deg)
## [1] 0
p <- ggplot(data = data.frame(x=deg),aes(x=x))
p <- p + geom_histogram(breaks =seq(0,3600, by=300), fill='azure2',col='blue')
p <- p + scale_x_continuous('Weighted degree', breaks=seq(0,3600, by=300))
p <- p + ggtitle('Weighted degree distribution for User to User projection')
print(p)

Filtering

u2uFil <- delete_edges(u2u, E(u2u)[(E(u2u)$weight<3)])
summary(u2u)
## IGRAPH UNW- 3139 1266475 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | weight (e/n)
summary(u2uFil)
## IGRAPH UNW- 3139 405 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | weight (e/n)
degF <- degree(u2uFil, normalized = T)
j = 1
for (i in V(u2uFil)) {
  u2uFil <- set.vertex.attribute(u2uFil,'degree',i,degF[j])    
  j = j + 1
}
u2uFil <- delete_vertices(u2uFil,V(u2uFil)[V(u2uFil)$degree==0])
summary(u2uFil)
## IGRAPH UNW- 53 405 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | degree (v/n), weight (e/n)

Community detection

ComBet <- cluster_edge_betweenness(u2u)
save(ComBet,file='ComBet.Rdata')
load('ComBet.Rdata')
ComBet$bridges

m = 0
for(i in 1:137){
  test = which(ComBet$membership==853)
  if (length(test) > m){
    m = i
  }
}
comunityBet = which(ComBet$membership==m)
comunityB <- induced.subgraph(u2u,comunityBet)
plot(comunityB, vertex.label='')
ComFG <- cluster_fast_greedy(u2u)
save(ComFG,file='ComFG.Rdata')

Fast greedy

load('ComFG.Rdata')
comunity = which(ComFG$membership==1)
length(comunity)/length(ComFG$membership)
## [1] 0.2720612
comunityFG <- induced.subgraph(u2u,comunity)
summary(comunityFG)
## IGRAPH UNW- 854 364231 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | weight (e/n)
comFG.deg <-degree(comunityFG, mode='in',normalized = T)
u2uComFG <- delete_vertices(comunityFG,V(comunityFG)[V(comunityFG)$inDegree<=0.3243*1.7])
summary(u2uComFG)
## IGRAPH UNW- 0 0 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | weight (e/n)

Walktrap

ComW10 <- cluster_walktrap(u2u,steps = 10)
save(ComW10,file='ComW10.Rdata')
load('ComW10.Rdata')
ComW10
## IGRAPH clustering walktrap, groups: 7, mod: 0.6
## + groups:
##   $`1`
##     [1] "A103M039GIFX02" "A1058L4DCEX548" "A105FAAQ06FTWA" "A10ASLX7DTTB6Z"
##     [5] "A10VOEBL5S337W" "A114IO16WXJ6GS" "A114ZZDRXMFGEN" "A11CF8XQA1LV06"
##     [9] "A11DBMDU3X5QLE" "A11DQXAB0OS2CF" "A11E010Z0NZCBQ" "A11F45IVW2WAK" 
##    [13] "A11JZYOPFMDGXW" "A11NGPI6U63B08" "A11RBPGTVNBIBS" "A11RKOUZ0NPCNO"
##    [17] "A1241U6QCSX5YJ" "A125ENS9P5BNJC" "A127R35BL1SGCY" "A12AN7TOF9TCG8"
##    [21] "A12QJOEW7AQE5R" "A12R7SL8IJZO4R" "A12UD2P4JC1HYN" "A12WR6UKHK3EZC"
##    [25] "A134PF2TZ3LHTW" "A1396PD1X1GPHX" "A13E99566QZ178" "A13S5XAZPGCAVC"
##    [29] "A13STKDYQO5DQ8" "A13T70B27VMBGD" "A13UNI8CGB0GC1" "A13VVN5MJGJRRM"
##    [33] "A13W69Z3K0ZX82" "A142W8NH9SRXS0" "A1467QB8PAVKB1" "A14B6XDB9MGKVZ"
##   + ... omitted several groups/vertices
comunityW10 = which(ComW10$membership==1)
length(comunityW10)/length(ComW10$membership)
## [1] 0.3029627
comunityW10 <- induced.subgraph(u2u,comunityW10)
summary(comunityW10)
## IGRAPH UNW- 951 441067 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | weight (e/n)
comFG.deg <-degree(comunityW10, mode='in',normalized = T)
u2uComFG <- delete_vertices(comunityW10,V(comunityW10)[V(comunityW10)$inDegree<=0.3243*1.7])
summary(u2uComFG)
## IGRAPH UNW- 35 595 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | weight (e/n)

Eigenvector

ComEigen <- cluster_leading_eigen(u2u,weights = E(u2u)$weight)
save(ComEigen,file='ComEigen.Rdata')
load('ComEigen.Rdata')
ComEigen
## IGRAPH clustering leading eigenvector, groups: 6, mod: 0.62
## + groups:
##   $`1`
##     [1] "A1020G50MVAOLA" "A103DZF5M6Y6C0" "A10NTVNJWXWOWW" "A10OMQK3ZT67YC"
##     [5] "A10WON8F2TLFXY" "A114OOK94DZGY6" "A119COUQTI1L8P" "A11M8RJ3OVGT4X"
##     [9] "A11TYILTAFKPR3" "A11V8B7VYTBWYO" "A11YHTC4XZA7S7" "A11YOT86X3M4GU"
##    [13] "A12GQKJI1ARMB7" "A12UB9OEVFO2Y3" "A12WP9DPDKTA63" "A12WURMH84NQ58"
##    [17] "A12WZTC4YJ8ZEC" "A133HCUWZGFVRY" "A134FLH2KNA0JU" "A13LQ2RNZ1DCE" 
##    [21] "A13U7O5KGPUE9C" "A13YLV9JB3AQNH" "A13ZCL7UXEDF14" "A14228BBHSV633"
##    [25] "A142VB1X9G6MSQ" "A14IOCHMRFNL2Y" "A14PS3LPIO7UF8" "A15JZZEM0PSWXY"
##    [29] "A15LJVOOAF291U" "A15OYXDLQPOKGM" "A15U2WAT90V7DC" "A162CHKPA7CSXN"
##    [33] "A162GHNEKGDYVM" "A16H1VL03OZTDJ" "A16IAIV6KZHRI0" "A16KRCLLTHSB5F"
##   + ... omitted several groups/vertices
ComEigen$modularity
## [1] 0.6159343
comunityeig = which(ComEigen$membership==1)
length(comunityeig)/length(ComEigen$membership)
## [1] 0.2554954
comunityEig <- induced.subgraph(u2u,comunityeig)
summary(comunityEig)
## IGRAPH UNW- 802 321201 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | weight (e/n)
comFG.eig <-degree(comunityEig, mode='in',normalized = T)
u2uComEi <- delete_vertices(comunityEig,V(comunityEig)[V(comunityEig)$inDegree<=0.2862*1.6])
summary(u2uComEi)
## IGRAPH UNW- 28 378 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | weight (e/n)

Exporting User to User projection and communities to graphml file

write.graph(u2uFil,'u2uFil.graphml', format = c('graphml'))
write.graph(u2uComFG,'u2uComFG.graphml',format = c('graphml'))
write.graph(u2uComEi,'u2uComEi.graphml',format = c('graphml'))

CUG Test

Assortativity of the network by the Sex attribute

CUG does not like attributes = 0, so you will have to add 1 to the attribute when calling mycugtest.

assortativity(u2uFil,types1 = V(u2uFil)$Rating)
## [1] -0.02303069
u2uFilA <- mycugtest(u2uFil, assortativity, directed=F, cmode="edges", types1=V(u2uFil)$Rating)

Null hypothesis: assortativity of the projection with respect to the Rating attribute is similar to that of a random network of the same size and density.

The The null hypothesis is rejected

print.cug.test(u2uFilA)
## 
## Univariate Conditional Uniform Graph Test
## 
## Conditioning Method: edges 
## Graph Type: 
## Diagonal Used: FALSE 
## Replications: 1000 
## 
## Observed Value: -0.02303069 
## Pr(X>=Obs): 0.535 
## Pr(X<=Obs): 0.465
plot.cug.test(u2uFilA)

Null hypothesis: The degree assortativity of the usr to user network is similar to that of a random network of the same size and density. The null hypothesis is rejected. degree assortativity is much lower than a random network

assortativity.degree(u2uFil,directed = F)
## [1] -0.2881537
u2uFilDegree <- mycugtest(u2uFil, assortativity.degree, directed=F, cmode="edges")
print.cug.test(u2uFilDegree)
## 
## Univariate Conditional Uniform Graph Test
## 
## Conditioning Method: edges 
## Graph Type: 
## Diagonal Used: FALSE 
## Replications: 1000 
## 
## Observed Value: -0.2881537 
## Pr(X>=Obs): 1 
## Pr(X<=Obs): 0
plot.cug.test(u2uFilDegree)

transitivity of the network

Null hypothesis: The transitivity of the user to user network is similar to that of a random network of the same size and density. The null hypothesis is rejected. Transitivity is much higher than a random network.

transitivity(u2uFil,type = 'undirected')
## [1] 0.5882488
u2uFilT <- mycugtest(u2uFil, transitivity, directed=F, cmode="edges")
print.cug.test(u2uFilT)
## 
## Univariate Conditional Uniform Graph Test
## 
## Conditioning Method: edges 
## Graph Type: 
## Diagonal Used: FALSE 
## Replications: 1000 
## 
## Observed Value: 0.5882488 
## Pr(X>=Obs): 0 
## Pr(X<=Obs): 1
plot.cug.test(u2uFilT)

QAP | Assortativity of the network by degree

Null hypothesis: the ssortatitvity whoul be the same regardless the arrangement of the edges. We reject the null hypothesis since all the values are bellow

Null hypothesis is rejected. 100% of the probability mass of test networks has lower assortativity.

u2uFil.deg <- degree(u2uFil, mode = 'all',normalized = F)

j = 1
for (i in V(u2uFil)) {
  u2uFil <- set.vertex.attribute(u2uFil,'Degree',i,u2uFil.deg[j])    
  j = j + 1
}
summary(u2uFil)
## IGRAPH UNW- 53 405 -- 
## + attr: name (v/c), Rating (v/n), inDegree (v/n), betweenes (v/n),
## | degree (v/n), Degree (v/n), weight (e/n)
u2uFil_qap <- myqaptest(u2uFil,assortativity.nominal, types=V(u2uFil)$Degree, directed=F)
summary.qaptest(u2uFil_qap)
## 
## QAP Test Results
## 
## Estimated p-values:
##  p(f(perm) >= f(d)): 0 
##  p(f(perm) <= f(d)): 1 
## 
## Test Diagnostics:
##  Test Value (f(d)): 0.1611313 
##  Replications: 1000 
##  Distribution Summary:
##      Min:     -0.07192299 
##      1stQ:    -0.03824288 
##      Med:     -0.02986611 
##      Mean:    -0.02916207 
##      3rdQ:    -0.02106631 
##      Max:     0.03273163
print.qaptest(u2uFil_qap)
## 
## QAP Test Results
## 
## Estimated p-values:
##  p(f(perm) >= f(d)): 0 
##  p(f(perm) <= f(d)): 1
plot.qaptest(u2uFil_qap)

u2uNet <- asNetwork(u2uFil)
u2uM1 <- ergm(u2uNet ~ edges+nodemix('Degree',base=c(1,2)))
## Observed statistic(s) mix.Degree.6.6, mix.Degree.1.7, mix.Degree.6.7, mix.Degree.1.9, mix.Degree.6.9, mix.Degree.7.9, mix.Degree.1.10, mix.Degree.6.10, mix.Degree.7.10, mix.Degree.9.10, mix.Degree.1.11, mix.Degree.6.11, mix.Degree.7.11, mix.Degree.9.11, mix.Degree.10.11, mix.Degree.1.12, mix.Degree.7.12, mix.Degree.9.12, mix.Degree.10.12, mix.Degree.1.13, mix.Degree.6.13, mix.Degree.9.13, mix.Degree.11.13, mix.Degree.12.13, mix.Degree.1.15, mix.Degree.6.15, mix.Degree.9.15, mix.Degree.10.15, mix.Degree.15.15, mix.Degree.1.17, mix.Degree.6.17, mix.Degree.7.17, mix.Degree.9.17, mix.Degree.10.17, mix.Degree.11.17, mix.Degree.12.17, mix.Degree.13.17, mix.Degree.15.17, mix.Degree.1.30, mix.Degree.6.30, mix.Degree.7.30, mix.Degree.10.30, mix.Degree.12.30, mix.Degree.13.30, mix.Degree.15.30, mix.Degree.1.50, and mix.Degree.52.52 are at their smallest attainable values. Their coefficients will be fixed at -Inf.
## Evaluating log-likelihood at the estimate.
summary(u2uM1)
## 
## ==========================
## Summary of model fit
## ==========================
## 
## Formula:   u2uNet ~ edges + nodemix("Degree", base = c(1, 2))
## 
## Iterations:  17 out of 20 
## 
## Monte Carlo MLE Results:
##                  Estimate Std. Error MCMC % p-value    
## edges              -19.16    5066.66      0   0.997    
## mix.Degree.6.6       -Inf       0.00      0  <1e-04 ***
## mix.Degree.1.7       -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.7       -Inf       0.00      0  <1e-04 ***
## mix.Degree.7.7      37.73    8259.29      0   0.996    
## mix.Degree.1.9       -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.9       -Inf       0.00      0  <1e-04 ***
## mix.Degree.7.9       -Inf       0.00      0  <1e-04 ***
## mix.Degree.9.9      38.32    7165.34      0   0.996    
## mix.Degree.1.10      -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.10      -Inf       0.00      0  <1e-04 ***
## mix.Degree.7.10      -Inf       0.00      0  <1e-04 ***
## mix.Degree.9.10      -Inf       0.00      0  <1e-04 ***
## mix.Degree.10.10    39.64    6707.01      0   0.995    
## mix.Degree.1.11      -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.11      -Inf       0.00      0  <1e-04 ***
## mix.Degree.7.11      -Inf       0.00      0  <1e-04 ***
## mix.Degree.9.11      -Inf       0.00      0  <1e-04 ***
## mix.Degree.10.11     -Inf       0.00      0  <1e-04 ***
## mix.Degree.11.11    18.24    5066.66      0   0.997    
## mix.Degree.1.12      -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.12     38.32    7165.34      0   0.996    
## mix.Degree.7.12      -Inf       0.00      0  <1e-04 ***
## mix.Degree.9.12      -Inf       0.00      0  <1e-04 ***
## mix.Degree.10.12     -Inf       0.00      0  <1e-04 ***
## mix.Degree.11.12    18.47    5066.66      0   0.997    
## mix.Degree.12.12    38.32    7165.34      0   0.996    
## mix.Degree.1.13      -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.13      -Inf       0.00      0  <1e-04 ***
## mix.Degree.7.13     38.52    7023.16      0   0.996    
## mix.Degree.9.13      -Inf       0.00      0  <1e-04 ***
## mix.Degree.10.13    39.44    6735.91      0   0.995    
## mix.Degree.11.13     -Inf       0.00      0  <1e-04 ***
## mix.Degree.12.13     -Inf       0.00      0  <1e-04 ***
## mix.Degree.13.13    37.73    8259.29      0   0.996    
## mix.Degree.1.15      -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.15      -Inf       0.00      0  <1e-04 ***
## mix.Degree.7.15     38.06    7446.35      0   0.996    
## mix.Degree.9.15      -Inf       0.00      0  <1e-04 ***
## mix.Degree.10.15     -Inf       0.00      0  <1e-04 ***
## mix.Degree.11.15    18.47    5066.66      0   0.997    
## mix.Degree.12.15    38.32    7165.34      0   0.996    
## mix.Degree.13.15    38.06    7446.35      0   0.996    
## mix.Degree.15.15     -Inf       0.00      0  <1e-04 ***
## mix.Degree.1.17      -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.17      -Inf       0.00      0  <1e-04 ***
## mix.Degree.7.17      -Inf       0.00      0  <1e-04 ***
## mix.Degree.9.17      -Inf       0.00      0  <1e-04 ***
## mix.Degree.10.17     -Inf       0.00      0  <1e-04 ***
## mix.Degree.11.17     -Inf       0.00      0  <1e-04 ***
## mix.Degree.12.17     -Inf       0.00      0  <1e-04 ***
## mix.Degree.13.17     -Inf       0.00      0  <1e-04 ***
## mix.Degree.15.17     -Inf       0.00      0  <1e-04 ***
## mix.Degree.17.17    40.88    6622.79      0   0.995    
## mix.Degree.1.30      -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.30      -Inf       0.00      0  <1e-04 ***
## mix.Degree.7.30      -Inf       0.00      0  <1e-04 ***
## mix.Degree.9.30     39.44    6735.91      0   0.995    
## mix.Degree.10.30     -Inf       0.00      0  <1e-04 ***
## mix.Degree.11.30    19.85    5066.66      0   0.997    
## mix.Degree.12.30     -Inf       0.00      0  <1e-04 ***
## mix.Degree.13.30     -Inf       0.00      0  <1e-04 ***
## mix.Degree.15.30     -Inf       0.00      0  <1e-04 ***
## mix.Degree.17.30    40.67    6630.69      0   0.995    
## mix.Degree.30.30    38.84    6879.97      0   0.995    
## mix.Degree.1.50      -Inf       0.00      0  <1e-04 ***
## mix.Degree.6.50     38.06    7446.35      0   0.996    
## mix.Degree.7.50     38.52    7023.16      0   0.996    
## mix.Degree.9.50     38.84    6879.97      0   0.995    
## mix.Degree.10.50    39.44    6735.91      0   0.995    
## mix.Degree.11.50    40.30    6649.13      0   0.995    
## mix.Degree.12.50    38.84    6879.97      0   0.995    
## mix.Degree.13.50    38.52    7023.16      0   0.996    
## mix.Degree.15.50    38.06    7446.35      0   0.996    
## mix.Degree.17.50    40.00    6670.19      0   0.995    
## mix.Degree.30.50    39.08    6808.04      0   0.995    
## mix.Degree.50.50    37.73    8259.29      0   0.996    
## mix.Degree.1.52     38.06    7446.35      0   0.996    
## mix.Degree.6.52     37.73    8259.29      0   0.996    
## mix.Degree.7.52     38.06    7446.35      0   0.996    
## mix.Degree.9.52     38.32    7165.34      0   0.996    
## mix.Degree.10.52    38.84    6879.97      0   0.995    
## mix.Degree.11.52    39.64    6707.01      0   0.995    
## mix.Degree.12.52    38.32    7165.34      0   0.996    
## mix.Degree.13.52    38.06    7446.35      0   0.996    
## mix.Degree.15.52    37.73    8259.29      0   0.996    
## mix.Degree.17.52    39.36    6749.04      0   0.995    
## mix.Degree.30.52    38.52    7023.16      0   0.996    
## mix.Degree.50.52    38.06    7446.35      0   0.996    
## mix.Degree.52.52     -Inf       0.00      0  <1e-04 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##      Null Deviance: 1910  on 1378  degrees of freedom
##  Residual Deviance: 1436  on 1288  degrees of freedom
##  
## AIC: 1616    BIC: 2087    (Smaller is better.) 
## 
##  Warning: The following terms have infinite coefficient estimates:
##   mix.Degree.6.6 mix.Degree.1.7 mix.Degree.6.7 mix.Degree.1.9 mix.Degree.6.9 mix.Degree.7.9 mix.Degree.1.10 mix.Degree.6.10 mix.Degree.7.10 mix.Degree.9.10 mix.Degree.1.11 mix.Degree.6.11 mix.Degree.7.11 mix.Degree.9.11 mix.Degree.10.11 mix.Degree.1.12 mix.Degree.7.12 mix.Degree.9.12 mix.Degree.10.12 mix.Degree.1.13 mix.Degree.6.13 mix.Degree.9.13 mix.Degree.11.13 mix.Degree.12.13 mix.Degree.1.15 mix.Degree.6.15 mix.Degree.9.15 mix.Degree.10.15 mix.Degree.15.15 mix.Degree.1.17 mix.Degree.6.17 mix.Degree.7.17 mix.Degree.9.17 mix.Degree.10.17 mix.Degree.11.17 mix.Degree.12.17 mix.Degree.13.17 mix.Degree.15.17 mix.Degree.1.30 mix.Degree.6.30 mix.Degree.7.30 mix.Degree.10.30 mix.Degree.12.30 mix.Degree.13.30 mix.Degree.15.30 mix.Degree.1.50 mix.Degree.52.52

Plottin a simulated version of the network

gaplot(asIgraph(simulate(u2uM1)),names=F)

Plotting the parameter fit

Looking at the plot, we see that the model does not capture the connection 1-2 (male to unkown) as well as the other statistics. However the p-values for the sex feature are all good, but not good for the edges.

u2uM1gofm <-gof(u2uM1,GOF=~model)
summary(u2uM1gofm)
## 
## Goodness-of-fit for model statistics 
## 
##                  obs min   mean max MC p-value
## edges            405 386 403.95 421       0.84
## mix.Degree.7.7     1   1   1.00   1       1.00
## mix.Degree.9.9     3   3   3.00   3       1.00
## mix.Degree.10.10  15  15  15.00  15       1.00
## mix.Degree.11.11  30  18  29.38  41       0.92
## mix.Degree.6.12    3   3   3.00   3       1.00
## mix.Degree.11.12  15   7  14.79  22       1.00
## mix.Degree.12.12   3   3   3.00   3       1.00
## mix.Degree.7.13    4   4   4.00   4       1.00
## mix.Degree.10.13  12  12  12.00  12       1.00
## mix.Degree.13.13   1   1   1.00   1       1.00
## mix.Degree.7.15    2   2   2.00   2       1.00
## mix.Degree.11.15   5   2   4.89   8       1.00
## mix.Degree.12.15   3   3   3.00   3       1.00
## mix.Degree.13.15   2   2   2.00   2       1.00
## mix.Degree.17.17  55  55  55.00  55       1.00
## mix.Degree.9.30   12  12  12.00  12       1.00
## mix.Degree.11.30  40  32  39.89  49       1.00
## mix.Degree.17.30  44  44  44.00  44       1.00
## mix.Degree.30.30   6   6   6.00   6       1.00
## mix.Degree.6.50    2   2   2.00   2       1.00
## mix.Degree.7.50    4   4   4.00   4       1.00
## mix.Degree.9.50    6   6   6.00   6       1.00
## mix.Degree.10.50  12  12  12.00  12       1.00
## mix.Degree.11.50  30  30  30.00  30       1.00
## mix.Degree.12.50   6   6   6.00   6       1.00
## mix.Degree.13.50   4   4   4.00   4       1.00
## mix.Degree.15.50   2   2   2.00   2       1.00
## mix.Degree.17.50  22  22  22.00  22       1.00
## mix.Degree.30.50   8   8   8.00   8       1.00
## mix.Degree.50.50   1   1   1.00   1       1.00
## mix.Degree.1.52    2   2   2.00   2       1.00
## mix.Degree.6.52    1   1   1.00   1       1.00
## mix.Degree.7.52    2   2   2.00   2       1.00
## mix.Degree.9.52    3   3   3.00   3       1.00
## mix.Degree.10.52   6   6   6.00   6       1.00
## mix.Degree.11.52  15  15  15.00  15       1.00
## mix.Degree.12.52   3   3   3.00   3       1.00
## mix.Degree.13.52   2   2   2.00   2       1.00
## mix.Degree.15.52   1   1   1.00   1       1.00
## mix.Degree.17.52  11  11  11.00  11       1.00
## mix.Degree.30.52   4   4   4.00   4       1.00
## mix.Degree.50.52   2   2   2.00   2       1.00
plot(u2uM1gofm)

Plotting the goodness of fit

u2uM1gof <- gof(u2uM1)
plot(u2uM1gof)

## Warning in max((1:(n - 1))[x$pval.dist[1:(n - 1), "MC p-value"] < 1]): no
## non-missing arguments to max; returning -Inf

u2uM2 <- ergm(u2uNet ~ edges+nodemix('Degree',base=c(1,2))+degree(9))
## Observed statistic(s) mix.Degree.6.6, mix.Degree.1.7, mix.Degree.6.7, mix.Degree.1.9, mix.Degree.6.9, mix.Degree.7.9, mix.Degree.1.10, mix.Degree.6.10, mix.Degree.7.10, mix.Degree.9.10, mix.Degree.1.11, mix.Degree.6.11, mix.Degree.7.11, mix.Degree.9.11, mix.Degree.10.11, mix.Degree.1.12, mix.Degree.7.12, mix.Degree.9.12, mix.Degree.10.12, mix.Degree.1.13, mix.Degree.6.13, mix.Degree.9.13, mix.Degree.11.13, mix.Degree.12.13, mix.Degree.1.15, mix.Degree.6.15, mix.Degree.9.15, mix.Degree.10.15, mix.Degree.15.15, mix.Degree.1.17, mix.Degree.6.17, mix.Degree.7.17, mix.Degree.9.17, mix.Degree.10.17, mix.Degree.11.17, mix.Degree.12.17, mix.Degree.13.17, mix.Degree.15.17, mix.Degree.1.30, mix.Degree.6.30, mix.Degree.7.30, mix.Degree.10.30, mix.Degree.12.30, mix.Degree.13.30, mix.Degree.15.30, mix.Degree.1.50, and mix.Degree.52.52 are at their smallest attainable values. Their coefficients will be fixed at -Inf.
## Starting maximum likelihood estimation via MCMLE:
## Iteration 1 of at most 20: 
## The log-likelihood improved by 1.237 
## Iteration 2 of at most 20: 
## The log-likelihood improved by 0.3803 
## Iteration 3 of at most 20: 
## The log-likelihood improved by 0.1583 
## Iteration 4 of at most 20: 
## The log-likelihood improved by 0.05655 
## Iteration 5 of at most 20: 
## The log-likelihood improved by 0.03513 
## Iteration 6 of at most 20: 
## The log-likelihood improved by 0.0136 
## Iteration 7 of at most 20: 
## The log-likelihood improved by 0.01572 
## Iteration 8 of at most 20: 
## The log-likelihood improved by 0.8226 
## Iteration 9 of at most 20: 
## The log-likelihood improved by 0.0214 
## Iteration 10 of at most 20: 
## The log-likelihood improved by 0.8759 
## Iteration 11 of at most 20: 
## The log-likelihood improved by 0.6214 
## Iteration 12 of at most 20: 
## The log-likelihood improved by 0.9331 
## Iteration 13 of at most 20: 
## The log-likelihood improved by 0.03958 
## Iteration 14 of at most 20: 
## The log-likelihood improved by 0.9064 
## Iteration 15 of at most 20: 
## The log-likelihood improved by 4.539 
## Iteration 16 of at most 20: 
## The log-likelihood improved by 0.4141 
## Iteration 17 of at most 20: 
## The log-likelihood improved by 3.701 
## Iteration 18 of at most 20: 
## The log-likelihood improved by 1.646 
## Iteration 19 of at most 20: 
## The log-likelihood improved by 0.195 
## Iteration 20 of at most 20:
## Warning in ergm.MCMCse.lognormal(theta = theta, init = init, statsmatrix =
## statsmatrix0, : Approximate Hessian matrix is singular. Standard errors due
## to MCMC approximation of the likelihood cannot be evaluated. This is likely
## due to insufficient MCMC sample size or highly correlated model terms.
## The log-likelihood improved by 0.006558 
## Step length converged once. Increasing MCMC sample size.
## MCMLE estimation did not converge after 20 iterations. The estimated coefficients may not be accurate. Estimation may be resumed by passing the coefficients as initial values; see 'init' under ?control.ergm for details.
## Evaluating log-likelihood at the estimate. Using 20 bridges: 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 .
## 
## This model was fit using MCMC.  To examine model diagnostics and check for degeneracy, use the mcmc.diagnostics() function.
save(u2uM2,file=('dolM2.Rdata'))

gaplot(asIgraph(simulate(u2uM2)),names=F)

u2uM2gofm <-gof(u2uM2,GOF=~model)
summary(u2uM2gofm)
## 
## Goodness-of-fit for model statistics 
## 
##                  obs min   mean max MC p-value
## edges            405 398 408.22 418       0.44
## mix.Degree.7.7     1   1   1.00   1       1.00
## mix.Degree.9.9     3   3   3.00   3       1.00
## mix.Degree.10.10  15  15  15.00  15       1.00
## mix.Degree.11.11  30  22  31.48  39       0.72
## mix.Degree.6.12    3   3   3.00   3       1.00
## mix.Degree.11.12  15  10  15.92  24       0.86
## mix.Degree.12.12   3   3   3.00   3       1.00
## mix.Degree.7.13    4   4   4.00   4       1.00
## mix.Degree.10.13  12  12  12.00  12       1.00
## mix.Degree.13.13   1   1   1.00   1       1.00
## mix.Degree.7.15    2   2   2.00   2       1.00
## mix.Degree.11.15   5   0   5.08  10       1.00
## mix.Degree.12.15   3   3   3.00   3       1.00
## mix.Degree.13.15   2   2   2.00   2       1.00
## mix.Degree.17.17  55  55  55.00  55       1.00
## mix.Degree.9.30   12  12  12.00  12       1.00
## mix.Degree.11.30  40  32  40.74  51       0.90
## mix.Degree.17.30  44  44  44.00  44       1.00
## mix.Degree.30.30   6   6   6.00   6       1.00
## mix.Degree.6.50    2   2   2.00   2       1.00
## mix.Degree.7.50    4   4   4.00   4       1.00
## mix.Degree.9.50    6   6   6.00   6       1.00
## mix.Degree.10.50  12  12  12.00  12       1.00
## mix.Degree.11.50  30  30  30.00  30       1.00
## mix.Degree.12.50   6   6   6.00   6       1.00
## mix.Degree.13.50   4   4   4.00   4       1.00
## mix.Degree.15.50   2   2   2.00   2       1.00
## mix.Degree.17.50  22  22  22.00  22       1.00
## mix.Degree.30.50   8   8   8.00   8       1.00
## mix.Degree.50.50   1   1   1.00   1       1.00
## mix.Degree.1.52    2   2   2.00   2       1.00
## mix.Degree.6.52    1   1   1.00   1       1.00
## mix.Degree.7.52    2   2   2.00   2       1.00
## mix.Degree.9.52    3   3   3.00   3       1.00
## mix.Degree.10.52   6   6   6.00   6       1.00
## mix.Degree.11.52  15  15  15.00  15       1.00
## mix.Degree.12.52   3   3   3.00   3       1.00
## mix.Degree.13.52   2   2   2.00   2       1.00
## mix.Degree.15.52   1   1   1.00   1       1.00
## mix.Degree.17.52  11  11  11.00  11       1.00
## mix.Degree.30.52   4   4   4.00   4       1.00
## mix.Degree.50.52   2   2   2.00   2       1.00
## degree9            3   3   3.00   3       1.00
plot(u2uM2gofm)

u2uM2gof <- gof(u2uM2)
plot(u2uM2gof)

## Warning in max((1:(n - 1))[x$pval.dist[1:(n - 1), "MC p-value"] < 1]): no
## non-missing arguments to max; returning -Inf

mcmc.diagnostics(u2uM2, center=FALSE)
## Sample statistics summary:
## 
## Iterations = 16384:1063936
## Thinning interval = 1024 
## Number of chains = 1 
## Sample size per chain = 1024 
## 
## 1. Empirical mean and standard deviation for each variable,
##    plus standard error of the mean:
## 
##                     Mean    SD Naive SE Time-series SE
## edges            405.419 3.882  0.12132        0.16246
## mix.Degree.7.7     1.000 0.000  0.00000        0.00000
## mix.Degree.9.9     3.000 0.000  0.00000        0.00000
## mix.Degree.10.10  15.000 0.000  0.00000        0.00000
## mix.Degree.11.11  30.081 3.093  0.09664        0.13172
## mix.Degree.6.12    3.000 0.000  0.00000        0.00000
## mix.Degree.11.12  15.047 2.671  0.08347        0.10922
## mix.Degree.12.12   3.000 0.000  0.00000        0.00000
## mix.Degree.7.13    4.000 0.000  0.00000        0.00000
## mix.Degree.10.13  12.000 0.000  0.00000        0.00000
## mix.Degree.13.13   1.000 0.000  0.00000        0.00000
## mix.Degree.7.15    2.000 0.000  0.00000        0.00000
## mix.Degree.11.15   5.011 1.789  0.05590        0.08007
## mix.Degree.12.15   3.000 0.000  0.00000        0.00000
## mix.Degree.13.15   2.000 0.000  0.00000        0.00000
## mix.Degree.17.17  55.000 0.000  0.00000        0.00000
## mix.Degree.9.30   12.000 0.000  0.00000        0.00000
## mix.Degree.11.30  40.280 3.655  0.11420        0.22885
## mix.Degree.17.30  44.000 0.000  0.00000        0.00000
## mix.Degree.30.30   6.000 0.000  0.00000        0.00000
## mix.Degree.6.50    2.000 0.000  0.00000        0.00000
## mix.Degree.7.50    4.000 0.000  0.00000        0.00000
## mix.Degree.9.50    6.000 0.000  0.00000        0.00000
## mix.Degree.10.50  12.000 0.000  0.00000        0.00000
## mix.Degree.11.50  30.000 0.000  0.00000        0.00000
## mix.Degree.12.50   6.000 0.000  0.00000        0.00000
## mix.Degree.13.50   4.000 0.000  0.00000        0.00000
## mix.Degree.15.50   2.000 0.000  0.00000        0.00000
## mix.Degree.17.50  22.000 0.000  0.00000        0.00000
## mix.Degree.30.50   8.000 0.000  0.00000        0.00000
## mix.Degree.50.50   1.000 0.000  0.00000        0.00000
## mix.Degree.1.52    2.000 0.000  0.00000        0.00000
## mix.Degree.6.52    1.000 0.000  0.00000        0.00000
## mix.Degree.7.52    2.000 0.000  0.00000        0.00000
## mix.Degree.9.52    3.000 0.000  0.00000        0.00000
## mix.Degree.10.52   6.000 0.000  0.00000        0.00000
## mix.Degree.11.52  15.000 0.000  0.00000        0.00000
## mix.Degree.12.52   3.000 0.000  0.00000        0.00000
## mix.Degree.13.52   2.000 0.000  0.00000        0.00000
## mix.Degree.15.52   1.000 0.000  0.00000        0.00000
## mix.Degree.17.52  11.000 0.000  0.00000        0.00000
## mix.Degree.30.52   4.000 0.000  0.00000        0.00000
## mix.Degree.50.52   2.000 0.000  0.00000        0.00000
## degree9            3.000 0.000  0.00000        0.00000
## 
## 2. Quantiles for each variable:
## 
##                  2.5% 25%   50% 75%  97.5%
## edges             398 403 405.5 408 413.00
## mix.Degree.7.7      1   1   1.0   1   1.00
## mix.Degree.9.9      3   3   3.0   3   3.00
## mix.Degree.10.10   15  15  15.0  15  15.00
## mix.Degree.11.11   24  28  30.0  32  36.00
## mix.Degree.6.12     3   3   3.0   3   3.00
## mix.Degree.11.12   10  13  15.0  17  20.42
## mix.Degree.12.12    3   3   3.0   3   3.00
## mix.Degree.7.13     4   4   4.0   4   4.00
## mix.Degree.10.13   12  12  12.0  12  12.00
## mix.Degree.13.13    1   1   1.0   1   1.00
## mix.Degree.7.15     2   2   2.0   2   2.00
## mix.Degree.11.15    2   4   5.0   6   8.00
## mix.Degree.12.15    3   3   3.0   3   3.00
## mix.Degree.13.15    2   2   2.0   2   2.00
## mix.Degree.17.17   55  55  55.0  55  55.00
## mix.Degree.9.30    12  12  12.0  12  12.00
## mix.Degree.11.30   33  38  40.0  43  47.00
## mix.Degree.17.30   44  44  44.0  44  44.00
## mix.Degree.30.30    6   6   6.0   6   6.00
## mix.Degree.6.50     2   2   2.0   2   2.00
## mix.Degree.7.50     4   4   4.0   4   4.00
## mix.Degree.9.50     6   6   6.0   6   6.00
## mix.Degree.10.50   12  12  12.0  12  12.00
## mix.Degree.11.50   30  30  30.0  30  30.00
## mix.Degree.12.50    6   6   6.0   6   6.00
## mix.Degree.13.50    4   4   4.0   4   4.00
## mix.Degree.15.50    2   2   2.0   2   2.00
## mix.Degree.17.50   22  22  22.0  22  22.00
## mix.Degree.30.50    8   8   8.0   8   8.00
## mix.Degree.50.50    1   1   1.0   1   1.00
## mix.Degree.1.52     2   2   2.0   2   2.00
## mix.Degree.6.52     1   1   1.0   1   1.00
## mix.Degree.7.52     2   2   2.0   2   2.00
## mix.Degree.9.52     3   3   3.0   3   3.00
## mix.Degree.10.52    6   6   6.0   6   6.00
## mix.Degree.11.52   15  15  15.0  15  15.00
## mix.Degree.12.52    3   3   3.0   3   3.00
## mix.Degree.13.52    2   2   2.0   2   2.00
## mix.Degree.15.52    1   1   1.0   1   1.00
## mix.Degree.17.52   11  11  11.0  11  11.00
## mix.Degree.30.52    4   4   4.0   4   4.00
## mix.Degree.50.52    2   2   2.0   2   2.00
## degree9             3   3   3.0   3   3.00
## 
## 
## Sample statistics cross-correlations:

## 
## MCMC diagnostics shown here are from the last round of simulation, prior to computation of final parameter estimates. Because the final estimates are refinements of those used for this simulation run, these diagnostics may understate model performance. To directly assess the performance of the final model on in-model statistics, please use the GOF command: gof(ergmFitObject, GOF=~model).