ergm.count package
krhee
1/31/2021
# Method 1: From a socio matrix
## In many cases, a valued sociomatrix is available (or can easily be constructed). In this case, we will build one from the Sampson data.
library(ergm.count)## Loading required package: ergm## Loading required package: network## network: Classes for Relational Data
## Version 1.16.1 created on 2020-10-06.
## copyright (c) 2005, Carter T. Butts, University of California-Irvine
## Mark S. Handcock, University of California -- Los Angeles
## David R. Hunter, Penn State University
## Martina Morris, University of Washington
## Skye Bender-deMoll, University of Washington
## For citation information, type citation("network").
## Type help("network-package") to get started.##
## ergm: version 3.11.0, created on 2020-10-14
## Copyright (c) 2020, Mark S. Handcock, University of California -- Los Angeles
## David R. Hunter, Penn State University
## Carter T. Butts, University of California -- Irvine
## Steven M. Goodreau, University of Washington
## Pavel N. Krivitsky, UNSW Sydney
## Martina Morris, University of Washington
## with contributions from
## Li Wang
## Kirk Li, University of Washington
## Skye Bender-deMoll, University of Washington
## Chad Klumb
## Michał Bojanowski, Kozminski University
## Ben Bolker
## Based on "statnet" project software (statnet.org).
## For license and citation information see statnet.org/attribution
## or type citation("ergm").## NOTE: Versions before 3.6.1 had a bug in the implementation of the bd()
## constraint which distorted the sampled distribution somewhat. In
## addition, Sampson's Monks datasets had mislabeled vertices. See the
## NEWS and the documentation for more details.## NOTE: Some common term arguments pertaining to vertex attribute and
## level selection have changed in 3.10.0. See terms help for more
## details. Use 'options(ergm.term=list(version="3.9.4"))' to use old
## behavior.##
## ergm.count: version 3.4.0, created on 2019-05-15
## Copyright (c) 2019, Pavel N. Krivitsky, University of Wollongong
## with contributions from
## Mark S. Handcock, University of California -- Los Angeles
## David R. Hunter, Penn State University
## Based on "statnet" project software (statnet.org).
## For license and citation information see statnet.org/attribution
## or type citation("ergm.count").## NOTE: The form of the term 'CMP' has been changed in version 3.2 of
## 'ergm.count'. See the news or help('CMP') for more information.data(samplk)
ls()## [1] "samplk1" "samplk2" "samplk3"# Create a sociomatrix totaling the nominations
samplk.tot.m <- as.matrix(samplk1) + as.matrix(samplk2) + as.matrix(samplk3)
samplk.tot.m[1:5, 1:5]## John Bosco Gregory Basil Peter Bonaventure
## John Bosco 0 1 2 0 2
## Gregory 3 0 0 0 0
## Basil 3 1 0 0 0
## Peter 0 0 0 0 3
## Bonaventure 1 0 0 3 0# Create a network where the number of nominations becomes an attribute of an edge
samplk.tot <- as.network(samplk.tot.m, directed=TRUE, matrix.type="a", ignore.eval=FALSE, names.eval="nominations") # Important!
# Add vertex attributes (Note that the naems were already imported)
samplk.tot %v% "group" <- samplk1 %v% "group" # Groups identified by Sampson
samplk.tot %v% "group"## [1] "Turks" "Turks" "Outcasts" "Loyal" "Loyal" "Loyal"
## [7] "Turks" "Waverers" "Loyal" "Waverers" "Loyal" "Turks"
## [13] "Waverers" "Turks" "Turks" "Turks" "Outcasts" "Outcasts"# We can view the attribute as a sociomatrix
as.matrix(samplk.tot, attrname="nominations")[1:5, 1:5]## John Bosco Gregory Basil Peter Bonaventure
## John Bosco 0 1 2 0 2
## Gregory 3 0 0 0 0
## Basil 3 1 0 0 0
## Peter 0 0 0 0 3
## Bonaventure 1 0 0 3 0# Also, note that sampl.tot now has an edge if i nominated j *at least once*.
as.matrix(samplk.tot)[1:5, 1:5]## John Bosco Gregory Basil Peter Bonaventure
## John Bosco 0 1 1 0 1
## Gregory 1 0 0 0 0
## Basil 1 1 0 0 0
## Peter 0 0 0 0 1
## Bonaventure 1 0 0 1 0# Method 2: Form an edgelist
## Sociomatrices are simple to work with, but not very convenient for large, sparse networks. In the latter case, edgelists are often preferred. For our present case, suppose that instead of a sociomtrix we have an edgelist with values:
samplk.tot.el <- as.matrix(samplk.tot, attrname="nominations", matrix.type="edgelist")
samplk.tot.el[1:5,]## [,1] [,2] [,3]
## [1,] 2 1 3
## [2,] 3 1 3
## [3,] 5 1 1
## [4,] 6 1 2
## [5,] 7 1 1# and an initial empty network
samplk.tot2 <- samplk1 # Copy samplk1
delete.edges(samplk.tot2, seq_along(samplk.tot2$mel)) # Empty it out
samplk.tot2 # We could also have used network.initialize(18)## Network attributes:
## vertices = 18
## directed = TRUE
## hyper = FALSE
## loops = FALSE
## multiple = FALSE
## bipartite = FALSE
## total edges= 0
## missing edges= 0
## non-missing edges= 0
##
## Vertex attribute names:
## cloisterville group vertex.names
##
## No edge attributessamplk.tot2[samplk.tot.el[, 1:2], names.eval="nominations", add.edges=TRUE] <- samplk.tot.el[, 3]
as.matrix(samplk.tot2, attrname="nominations")[1:5, 1:5]## John Bosco Gregory Basil Peter Bonaventure
## John Bosco 0 1 2 0 2
## Gregory 3 0 0 0 0
## Basil 3 1 0 0 0
## Peter 0 0 0 0 3
## Bonaventure 1 0 0 3 0In general, the construction net[i, j, names.eval="attrname", add.edges=TRUE] <- value can be used to modify individual edge values for attribute "attrname". This way, we can also add more than one edge attribute to a network. Note that network objects can suppor tan almost unlimited number of vertex, edge, or network attributes, and that these attributes can contain any data type.
We can fit the equivalent of logistic regression on the probability of nomination, with every ordered pair of monks observed 3 times. We will look at differential homophily on group.
# summary(samplk.tot ~ sum)This produces an error because no such term has been implemented for binary mode, while
summary(samplk.tot ~ sum, response="nominations")## sum
## 168This gives the summary statistics. We will introduce more statistics.
New Concept: A Reference Distribution With binary ERGMs, we only concern ourselves with presence and absence of ties among actors – who is connected with whom? If we want to model values as well, we need to think about who is connected with whom and how strong or intense these connections are. In particular, we need to think about how the values for connections we measure are distributed. The reference distribution (a reference measure, for the mathematically inclined) specifies the model for the data ‘before’ we add any ERGM terms, and is the first step in modeling these values. The reference distribution is specified using a one-sided formula as a reference argument to an ergm or simulate call. Running help("ergm-references") will list the choices implemented in the various packages, and are given as a one-sided formula.
Conceptually, it has two ingredients: the sample space and the baseline distribution (h(y)). An ERGM that “borrows” these from a distribution X for which we have a name is called an “X-reference ERGM”.
The Sample Space: For binary ERGMs, the sample space ‘y’ – the set of possible networks that can occur – is usually some subset of 2^Y, the set of all possible ways in which relationships among the actors may occur.
For the sample space of valued ERGMs, we need to define S, the set of possible values each relationship may take. For example, for count data, that’s S = {0, 1, …, s} if the maximum count is s, and {0, 1, …} if there is no a priori upper bound. Having specified that, ‘y’ is defined as some subset of S^Y: the set of possible ways to assign to each relationship a value.
As with binary ERGMs, other constraints like degree distribution may be imposed on ‘y’.
h(y): The baseline distribution, what difference does it make? Suppose that we have a sample space with S = {0, 1, 2, 3} (e.g., number of monk-monk nominations) and let’s have one ERGM term: the sum of values of all relations: There are two values for h(y) that might be familiar: h(y) = 1 (or any constant) == Uniform or truncated geometric h(y) = m! / ((y_i,j)!(m-(y_i,j))!) == Binomial
Suppose that we don’t have an a priori upper bound on the counts (S = {0, 1, …}). Then there are two familiar reference distributions: h(y) = 1 (or any constant) == Geometric h(y) == Poisson
Let’s simulate these distributions
y <- network.initialize(2, directed=FALSE) # A network with one dyad
sim.geom <- simulate(y ~ sum, coef=log(2/3), reference=~Geometric, response="w",output="stats",nsim=1000)
sim.pois <- simulate(y ~ sum, coef=log(2), reference=~Poisson, response="w", output="stats", nsim=1000)
mean(sim.geom)## [1] 1.883mean(sim.pois)## [1] 2.015# Similar means, but what do they mean?par(mfrow = c(1, 2))
hist(sim.geom, breaks=diff(range(sim.geom))*4)
hist(sim.pois, breaks=diff(range(sim.pois))*4)
But using reference=~Geometric can be dangerous. This issue only arises with ERGMs that have an infinite sample space.
par(mfrow=c(1,1))
sim.geo0 <- simulate(y ~ sum, coef=0, reference=~Geometric, response="w", output="stats", nsim=100, control=control.simulate(MCMC.burnin=0, MCMC.interval=1))
mean(sim.geo0) #45446.96## [1] 302162.1plot(c(sim.geo0, xlab="MCMC iteration", ylab="Value of the tie"))## Warning in xy.coords(x, y, xlabel, ylabel, log): NAs introduced by coercion
Valued ERGM Terms GLM-style terms: Many of the familiar ERGM effects can be modeled using the very same terms in the valued case, but applied a little differently.
If y_i,j is allowed to have values other than 0 and 1, then simply using such a term in a Poisson-reference ERGM creates the familiar log-linear effect. Similarly, in a Binomial-reference ERGM, such terms produce an effect on log-odds of a success.
The good news is that almost every dyad-independent ergm term has been reimplemented to allow this. It is involed by specifying form="sum" argument for one of the terms inherited from binary ERGMs, though this not required, as it’s the default. Also, note that for valued ERGMs, the “intercept” term is sum, not edges.
help("ergm-terms) has the complete list across all the loaded package. In particular, the one in package ‘ergm’ has each term be tagged with whether it’s binary or valued.
Example: Sampson’s Monks We will look at differential homophily on group. That is, Y_i,j ~ IID (Binomial(3, pi_i,j))
samplk.tot.nm <- ergm(samplk.tot ~ sum + nodematch("group", diff=TRUE, form="sum"), response="nominations", reference=~Binomial(3))## Starting contrastive divergence estimation via CD-MCMLE:## Iteration 1 of at most 60:## Convergence test P-value:0e+00## Optimizing with step length 0.981645215654761.## The log-likelihood improved by 2.157.## Iteration 2 of at most 60:## Convergence test P-value:6.7e-249## Optimizing with step length 1.## The log-likelihood improved by 1.041.## Iteration 3 of at most 60:## Convergence test P-value:9.7e-82## Optimizing with step length 1.## The log-likelihood improved by 0.2265.## Iteration 4 of at most 60:## Convergence test P-value:5.9e-11## Optimizing with step length 1.## The log-likelihood improved by 0.02831.## Iteration 5 of at most 60:## Convergence test P-value:2.2e-05## Optimizing with step length 1.## The log-likelihood improved by 0.01441.## Iteration 6 of at most 60:## Convergence test P-value:2.9e-02## Optimizing with step length 1.## The log-likelihood improved by 0.006149.## Iteration 7 of at most 60:## Convergence test P-value:6.3e-02## Optimizing with step length 1.## The log-likelihood improved by 0.00533.## Iteration 8 of at most 60:## Convergence test P-value:1.2e-02## Optimizing with step length 1.## The log-likelihood improved by 0.007229.## Iteration 9 of at most 60:## Convergence test P-value:1.4e-01## Optimizing with step length 1.## The log-likelihood improved by 0.004116.## Iteration 10 of at most 60:## Convergence test P-value:1.1e-01## Optimizing with step length 1.## The log-likelihood improved by 0.004529.## Iteration 11 of at most 60:## Convergence test P-value:2e-02## Optimizing with step length 1.## The log-likelihood improved by 0.006655.## Iteration 12 of at most 60:## Convergence test P-value:2.4e-01## Optimizing with step length 1.## The log-likelihood improved by 0.003378.## Iteration 13 of at most 60:## Convergence test P-value:6.2e-04## Optimizing with step length 1.## The log-likelihood improved by 0.01102.## Iteration 14 of at most 60:## Convergence test P-value:2.9e-01## Optimizing with step length 1.## The log-likelihood improved by 0.003062.## Iteration 15 of at most 60:## Convergence test P-value:5e-02## Optimizing with step length 1.## The log-likelihood improved by 0.005441.## Iteration 16 of at most 60:## Convergence test P-value:8.5e-01## Convergence detected. Stopping.## Optimizing with step length 1.## The log-likelihood improved by 0.0009824.## Starting Monte Carlo maximum likelihood estimation (MCMLE):## Iteration 1 of at most 20:## Optimizing with step length 1.## The log-likelihood improved by 0.4308.## Step length converged once. Increasing MCMC sample size.## Iteration 2 of at most 20:## Optimizing with step length 1.## The log-likelihood improved by 0.01948.## Step length converged twice. Stopping.## Finished MCMLE.## Note: Null model likelihood calculation is not implemented for valued
## ERGMs at this time. This means that all likelihood-based inference
## (LRT, Analysis of Deviance, AIC, BIC, etc.) is only valid between
## models with the same reference distribution and constraints.## 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.mcmc.diagnostics(samplk.tot.nm)## Sample statistics summary:
##
## Iterations = 16384:4209664
## Thinning interval = 1024
## Number of chains = 1
## Sample size per chain = 4096
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## sum 1.809326 10.431 0.16299 0.19274
## nodematch.sum.group.Loyal -0.009277 3.826 0.05978 0.06373
## nodematch.sum.group.Outcasts -0.019775 1.788 0.02793 0.03275
## nodematch.sum.group.Turks 0.568359 5.599 0.08749 0.09578
## nodematch.sum.group.Waverers -0.010742 1.729 0.02701 0.03315
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## sum -18 -5 2 9 22
## nodematch.sum.group.Loyal -7 -3 0 3 8
## nodematch.sum.group.Outcasts -4 -1 0 1 3
## nodematch.sum.group.Turks -11 -3 1 4 11
## nodematch.sum.group.Waverers -3 -1 0 1 4
##
##
## Sample statistics cross-correlations:
## sum nodematch.sum.group.Loyal
## sum 1.0000000 0.344646966
## nodematch.sum.group.Loyal 0.3446470 1.000000000
## nodematch.sum.group.Outcasts 0.1970078 0.026255034
## nodematch.sum.group.Turks 0.5323071 -0.005168988
## nodematch.sum.group.Waverers 0.1799000 0.016122639
## nodematch.sum.group.Outcasts
## sum 0.19700780
## nodematch.sum.group.Loyal 0.02625503
## nodematch.sum.group.Outcasts 1.00000000
## nodematch.sum.group.Turks 0.02191042
## nodematch.sum.group.Waverers 0.01384062
## nodematch.sum.group.Turks
## sum 0.532307115
## nodematch.sum.group.Loyal -0.005168988
## nodematch.sum.group.Outcasts 0.021910418
## nodematch.sum.group.Turks 1.000000000
## nodematch.sum.group.Waverers 0.010748703
## nodematch.sum.group.Waverers
## sum 0.17990001
## nodematch.sum.group.Loyal 0.01612264
## nodematch.sum.group.Outcasts 0.01384062
## nodematch.sum.group.Turks 0.01074870
## nodematch.sum.group.Waverers 1.00000000
##
## Sample statistics auto-correlation:
## Chain 1
## sum nodematch.sum.group.Loyal nodematch.sum.group.Outcasts
## Lag 0 1.000000000 1.0000000000 1.00000000
## Lag 1024 0.165964968 0.0638165140 0.15776347
## Lag 2048 0.036220900 0.0024458254 0.03189946
## Lag 3072 0.012999808 -0.0293581564 0.02976116
## Lag 4096 -0.004267446 0.0137736932 0.02043776
## Lag 5120 0.006606039 -0.0005915909 0.01799074
## nodematch.sum.group.Turks nodematch.sum.group.Waverers
## Lag 0 1.0000000000 1.00000000
## Lag 1024 0.0901811719 0.17428977
## Lag 2048 -0.0020400387 0.05790912
## Lag 3072 0.0006293503 0.01982681
## Lag 4096 -0.0059560416 0.01091568
## Lag 5120 0.0224468065 -0.02675798
##
## Sample statistics burn-in diagnostic (Geweke):
## Chain 1
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## sum nodematch.sum.group.Loyal
## 0.36137 0.39484
## nodematch.sum.group.Outcasts nodematch.sum.group.Turks
## -0.37043 -0.06829
## nodematch.sum.group.Waverers
## -0.28131
##
## Individual P-values (lower = worse):
## sum nodematch.sum.group.Loyal
## 0.7178265 0.6929638
## nodematch.sum.group.Outcasts nodematch.sum.group.Turks
## 0.7110617 0.9455520
## nodematch.sum.group.Waverers
## 0.7784757
## Joint P-value (lower = worse): 0.9852957 .
##
## 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).Note that it looks like fitting the model twice. This is because the first run is using an approximation technique called ‘contrastive divergence’ to find a good starting value for the MLE fit.
summary(samplk.tot.nm)## Call:
## ergm(formula = samplk.tot ~ sum + nodematch("group", diff = TRUE,
## form = "sum"), response = "nominations", reference = ~Binomial(3))
##
## Iterations: 2 out of 20
##
## Monte Carlo MLE Results:
## Estimate Std. Error MCMC % z value Pr(>|z|)
## sum -2.3321 0.1308 0 -17.827 <1e-04 ***
## nodematch.sum.group.Loyal 2.2579 0.2883 0 7.831 <1e-04 ***
## nodematch.sum.group.Outcasts 3.5957 0.5755 0 6.248 <1e-04 ***
## nodematch.sum.group.Turks 2.2072 0.2203 0 10.021 <1e-04 ***
## nodematch.sum.group.Waverers 1.0750 0.5935 0 1.811 0.0701 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Null Deviance: 0.0 on 306 degrees of freedom
## Residual Deviance: -559.3 on 301 degrees of freedom
##
## Note that the null model likelihood and deviance are defined to be 0.
## This means that all likelihood-based inference (LRT, Analysis of
## Deviance, AIC, BIC, etc.) is only valid between models with the same
## reference distribution and constraints.
##
## AIC: -549.3 BIC: -530.7 (Smaller is better.)Based on this, we can say that the odds of a monk nominating another monk not in the same group during a given time step are exp(beta1) = exp(-2.3274) = 0.0975, that the odds of a Loyal Opposition monk nominating another Loyal Opposition monk are exp(beta2) = exp(2.259) = 9.5738 times higher etc…
Let’s do another example, with Zachary’s Karate Club data. We will use a Poisson log-linear model for the number of contexts in which each pair of individuals interacted, as a function of whether this individual is a faction leader (Mr. Hi or John A.) That is, Y_i,j ~ IID (Poisson(u_i,j)).
We will do this by constructing a dummy variable, a vertex attribute "leader":
data(zach)
unique(zach %v% "role")## [1] "Instructor" "Member" "President"# Vertex attr. 'leader' is TRUE for Hi and John, FALSE for others
zach %v% "leader" <- zach %v% "role" != "Member"
zach.lead <- ergm(zach ~ sum + nodefactor("leader"), response="contexts", reference=~Poisson)## Starting contrastive divergence estimation via CD-MCMLE:## Iteration 1 of at most 60:## Convergence test P-value:0e+00## Optimizing with step length 1.## The log-likelihood improved by 2.23.## Iteration 2 of at most 60:## Convergence test P-value:7.4e-107## Optimizing with step length 1.## The log-likelihood improved by 0.2751.## Iteration 3 of at most 60:## Convergence test P-value:6.3e-53## Optimizing with step length 1.## The log-likelihood improved by 0.1206.## Iteration 4 of at most 60:## Convergence test P-value:1.7e-29## Optimizing with step length 1.## The log-likelihood improved by 0.06513.## Iteration 5 of at most 60:## Convergence test P-value:4.6e-24## Optimizing with step length 1.## The log-likelihood improved by 0.05251.## Iteration 6 of at most 60:## Convergence test P-value:1.2e-13## Optimizing with step length 1.## The log-likelihood improved by 0.02885.## Iteration 7 of at most 60:## Convergence test P-value:1.6e-08## Optimizing with step length 1.## The log-likelihood improved by 0.0173.## Iteration 8 of at most 60:## Convergence test P-value:5.7e-03## Optimizing with step length 1.## The log-likelihood improved by 0.005009.## Iteration 9 of at most 60:## Convergence test P-value:1.1e-02## Optimizing with step length 1.## The log-likelihood improved by 0.004362.## Iteration 10 of at most 60:## Convergence test P-value:6e-01## Convergence detected. Stopping.## Optimizing with step length 1.## The log-likelihood improved by 0.0004981.## Starting Monte Carlo maximum likelihood estimation (MCMLE):## Iteration 1 of at most 20:## Optimizing with step length 0.784770980730095.## The log-likelihood improved by 4.861.## Iteration 2 of at most 20:## Optimizing with step length 1.## The log-likelihood improved by 0.1129.## Step length converged once. Increasing MCMC sample size.## Iteration 3 of at most 20:## Optimizing with step length 1.## The log-likelihood improved by 0.004734.## Step length converged twice. Stopping.## Finished MCMLE.## Note: Null model likelihood calculation is not implemented for valued
## ERGMs at this time. This means that all likelihood-based inference
## (LRT, Analysis of Deviance, AIC, BIC, etc.) is only valid between
## models with the same reference distribution and constraints.## 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.mcmc.diagnostics(zach.lead)## Sample statistics summary:
##
## Iterations = 16384:4209664
## Thinning interval = 1024
## Number of chains = 1
## Sample size per chain = 4096
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## sum 0.4353 15.09 0.2358 0.3489
## nodefactor.sum.leader.TRUE 0.9446 10.32 0.1613 0.3224
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## sum -28.62 -10 0 11.00 31
## nodefactor.sum.leader.TRUE -18.00 -6 1 7.25 22
##
##
## Sample statistics cross-correlations:
## sum nodefactor.sum.leader.TRUE
## sum 1.0000000 0.6034392
## nodefactor.sum.leader.TRUE 0.6034392 1.0000000
##
## Sample statistics auto-correlation:
## Chain 1
## sum nodefactor.sum.leader.TRUE
## Lag 0 1.00000000 1.00000000
## Lag 1024 0.27563218 0.54472974
## Lag 2048 0.14239063 0.32308868
## Lag 3072 0.08958532 0.19159214
## Lag 4096 0.03539085 0.14251610
## Lag 5120 0.01640211 0.08734523
##
## Sample statistics burn-in diagnostic (Geweke):
## Chain 1
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## sum nodefactor.sum.leader.TRUE
## -0.9293 -0.4680
##
## Individual P-values (lower = worse):
## sum nodefactor.sum.leader.TRUE
## 0.3527471 0.6398056
## Joint P-value (lower = worse): 0.6471615 .
##
## 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).summary(zach.lead)## Call:
## ergm(formula = zach ~ sum + nodefactor("leader"), response = "contexts",
## reference = ~Poisson)
##
## Iterations: 3 out of 20
##
## Monte Carlo MLE Results:
## Estimate Std. Error MCMC % z value Pr(>|z|)
## sum -1.22097 0.08311 0 -14.69 <1e-04 ***
## nodefactor.sum.leader.TRUE 1.44021 0.12153 0 11.85 <1e-04 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Null Deviance: 0.0 on 561 degrees of freedom
## Residual Deviance: -355.1 on 559 degrees of freedom
##
## Note that the null model likelihood and deviance are defined to be 0.
## This means that all likelihood-based inference (LRT, Analysis of
## Deviance, AIC, BIC, etc.) is only valid between models with the same
## reference distribution and constraints.
##
## AIC: -351.1 BIC: -342.5 (Smaller is better.)Based on this, we can say that the expected number of contexts of interaction between two non-leaders is exp(beta1) = exp(-1.2231) = 0.2943, that the xpected number of contexts of interaction between a leader and a non-leader is exp(beta2) = exp(1.4457) = 4.245 times higher, and that the expected number of contexts of interaction between the two leaders is exp(2beta2) = exp(2 1.4457) = 18.0201 times higher than that between two non-leaders. (Because the leaders were hostile to each other, this may not be a very good prediction).
Spartsity and zero-modification: It is often the case that in networks of counts, the network is sparse, yet if two actors do interact, their interaction count is relatively high. This amounts to zero-inflation.
We can model this using the binary-ERGM-based terms with the term nonzero and GLM-style terms with argument form="nonzero". For example,
samplk.tot.nm.nz <- ergm(samplk.tot ~ sum + nonzero + nodematch("group", diff=TRUE, form="sum"), response="nominations", reference=~Binomial(3))## Starting contrastive divergence estimation via CD-MCMLE:## Iteration 1 of at most 60:## Convergence test P-value:0e+00## Optimizing with step length 0.62040703403172.## The log-likelihood improved by 1.848.## Iteration 2 of at most 60:## Convergence test P-value:2.8e-258## Optimizing with step length 1.## The log-likelihood improved by 1.083.## Iteration 3 of at most 60:## Convergence test P-value:5.2e-62## Optimizing with step length 1.## The log-likelihood improved by 0.1786.## Iteration 4 of at most 60:## Convergence test P-value:7.1e-08## Optimizing with step length 1.## The log-likelihood improved by 0.02218.## Iteration 5 of at most 60:## Convergence test P-value:3.6e-06## Optimizing with step length 1.## The log-likelihood improved by 0.01767.## Iteration 6 of at most 60:## Convergence test P-value:3.9e-01## Optimizing with step length 1.## The log-likelihood improved by 0.003098.## Iteration 7 of at most 60:## Convergence test P-value:6.7e-01## Convergence detected. Stopping.## Optimizing with step length 1.## The log-likelihood improved by 0.001991.## Starting Monte Carlo maximum likelihood estimation (MCMLE):## Iteration 1 of at most 20:## Optimizing with step length 1.## The log-likelihood improved by 1.31.## Step length converged once. Increasing MCMC sample size.## Iteration 2 of at most 20:## Optimizing with step length 1.## The log-likelihood improved by 0.02304.## Step length converged twice. Stopping.## Finished MCMLE.## Note: Null model likelihood calculation is not implemented for valued
## ERGMs at this time. This means that all likelihood-based inference
## (LRT, Analysis of Deviance, AIC, BIC, etc.) is only valid between
## models with the same reference distribution and constraints.## 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.mcmc.diagnostics(samplk.tot.nm.nz)## Sample statistics summary:
##
## Iterations = 16384:4209664
## Thinning interval = 1024
## Number of chains = 1
## Sample size per chain = 4096
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## sum 0.45801 13.916 0.21744 0.33494
## nonzero 0.42212 6.949 0.10858 0.16646
## nodematch.sum.group.Loyal -0.30591 5.200 0.08124 0.12032
## nodematch.sum.group.Outcasts 0.06543 2.221 0.03470 0.05058
## nodematch.sum.group.Turks -0.55786 7.624 0.11912 0.17873
## nodematch.sum.group.Waverers 0.41284 2.480 0.03875 0.06164
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## sum -26 -9 0 10 27
## nonzero -13 -4 0 5 14
## nodematch.sum.group.Loyal -11 -4 0 3 10
## nodematch.sum.group.Outcasts -5 -1 0 2 4
## nodematch.sum.group.Turks -15 -6 -1 5 14
## nodematch.sum.group.Waverers -4 -1 0 2 5
##
##
## Sample statistics cross-correlations:
## sum nonzero nodematch.sum.group.Loyal
## sum 1.0000000 0.87779563 0.386404854
## nonzero 0.8777956 1.00000000 0.263044846
## nodematch.sum.group.Loyal 0.3864049 0.26304485 1.000000000
## nodematch.sum.group.Outcasts 0.1459605 0.06579631 -0.004610990
## nodematch.sum.group.Turks 0.5544391 0.40113521 0.014753965
## nodematch.sum.group.Waverers 0.1486642 0.13344602 0.003508859
## nodematch.sum.group.Outcasts
## sum 0.14596051
## nonzero 0.06579631
## nodematch.sum.group.Loyal -0.00461099
## nodematch.sum.group.Outcasts 1.00000000
## nodematch.sum.group.Turks -0.00756570
## nodematch.sum.group.Waverers -0.04689825
## nodematch.sum.group.Turks
## sum 0.55443909
## nonzero 0.40113521
## nodematch.sum.group.Loyal 0.01475397
## nodematch.sum.group.Outcasts -0.00756570
## nodematch.sum.group.Turks 1.00000000
## nodematch.sum.group.Waverers -0.02222342
## nodematch.sum.group.Waverers
## sum 0.148664239
## nonzero 0.133446023
## nodematch.sum.group.Loyal 0.003508859
## nodematch.sum.group.Outcasts -0.046898249
## nodematch.sum.group.Turks -0.022223421
## nodematch.sum.group.Waverers 1.000000000
##
## Sample statistics auto-correlation:
## Chain 1
## sum nonzero nodematch.sum.group.Loyal
## Lag 0 1.000000000 1.000000000 1.000000000
## Lag 1024 0.386746026 0.376427961 0.373561546
## Lag 2048 0.169842054 0.168381204 0.146519338
## Lag 3072 0.068800835 0.068135039 0.046619911
## Lag 4096 0.022062851 0.029405267 -0.006721694
## Lag 5120 -0.009772093 -0.007917907 -0.018647947
## nodematch.sum.group.Outcasts nodematch.sum.group.Turks
## Lag 0 1.00000000 1.00000000
## Lag 1024 0.35989981 0.36197849
## Lag 2048 0.11777958 0.15391963
## Lag 3072 0.05518176 0.07568681
## Lag 4096 0.04125574 0.04947828
## Lag 5120 0.02312996 0.01595116
## nodematch.sum.group.Waverers
## Lag 0 1.00000000
## Lag 1024 0.41333870
## Lag 2048 0.19101731
## Lag 3072 0.08826120
## Lag 4096 0.03857880
## Lag 5120 0.02185699
##
## Sample statistics burn-in diagnostic (Geweke):
## Chain 1
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## sum nonzero
## -1.6923 -1.4225
## nodematch.sum.group.Loyal nodematch.sum.group.Outcasts
## -0.7243 -1.5018
## nodematch.sum.group.Turks nodematch.sum.group.Waverers
## 0.4785 -0.3824
##
## Individual P-values (lower = worse):
## sum nonzero
## 0.09059308 0.15486914
## nodematch.sum.group.Loyal nodematch.sum.group.Outcasts
## 0.46885918 0.13314561
## nodematch.sum.group.Turks nodematch.sum.group.Waverers
## 0.63226370 0.70215465
## Joint P-value (lower = worse): 0.3290667 .
##
## 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).summary(samplk.tot.nm.nz)## Call:
## ergm(formula = samplk.tot ~ sum + nonzero + nodematch("group",
## diff = TRUE, form = "sum"), response = "nominations", reference = ~Binomial(3))
##
## Iterations: 2 out of 20
##
## Monte Carlo MLE Results:
## Estimate Std. Error MCMC % z value Pr(>|z|)
## sum -0.3346 0.1978 0 -1.691 0.0908 .
## nonzero -2.9828 0.3236 0 -9.218 <1e-04 ***
## nodematch.sum.group.Loyal 1.2181 0.2265 0 5.377 <1e-04 ***
## nodematch.sum.group.Outcasts 1.9802 0.4714 0 4.201 <1e-04 ***
## nodematch.sum.group.Turks 1.1950 0.1750 0 6.827 <1e-04 ***
## nodematch.sum.group.Waverers 0.6076 0.4159 0 1.461 0.1441
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Null Deviance: 0.0 on 306 degrees of freedom
## Residual Deviance: -652.4 on 300 degrees of freedom
##
## Note that the null model likelihood and deviance are defined to be 0.
## This means that all likelihood-based inference (LRT, Analysis of
## Deviance, AIC, BIC, etc.) is only valid between models with the same
## reference distribution and constraints.
##
## AIC: -640.4 BIC: -618 (Smaller is better.)This fits a zero-modified Binomial model, with a coefficient on the number of non-zero relations -2.9799 is negative and highly significant, indicating that there is an excess of zeros in the data relative to the binomial distribution, and gives the rest of the model.
More examples: Sampson’s Monks Suppose that we want to fit a model with a zero-modified Biomodal baseline, mutuality, transitive (hierarchical) triads, and cyclical (antihierarchical) triads, to this dataset:
samplk.tot.ergm <- ergm(samplk.tot ~ sum + nonzero + mutual("min") + transitiveweights("min", "max", "min") + cyclicalweights("min", "max", "min"), reference=~Binomial(3), response="nominations")## Starting contrastive divergence estimation via CD-MCMLE:## Iteration 1 of at most 60:## Convergence test P-value:0e+00## Optimizing with step length 0.838434268741269.## The log-likelihood improved by 2.25.## Iteration 2 of at most 60:## Convergence test P-value:2.9e-276## Optimizing with step length 1.## The log-likelihood improved by 1.09.## Iteration 3 of at most 60:## Convergence test P-value:8.9e-104## Optimizing with step length 1.## The log-likelihood improved by 0.3285.## Iteration 4 of at most 60:## Convergence test P-value:2.4e-43## Optimizing with step length 1.## The log-likelihood improved by 0.1187.## Iteration 5 of at most 60:## Convergence test P-value:4.2e-10## Optimizing with step length 1.## The log-likelihood improved by 0.0269.## Iteration 6 of at most 60:## Convergence test P-value:2.4e-04## Optimizing with step length 1.## The log-likelihood improved by 0.01186.## Iteration 7 of at most 60:## Convergence test P-value:9.5e-01## Convergence detected. Stopping.## Optimizing with step length 1.## The log-likelihood improved by 0.0005811.## Starting Monte Carlo maximum likelihood estimation (MCMLE):## Iteration 1 of at most 20:## Optimizing with step length 1.## The log-likelihood improved by 2.891.## Step length converged once. Increasing MCMC sample size.## Iteration 2 of at most 20:## Optimizing with step length 1.## The log-likelihood improved by 0.07167.## Step length converged twice. Stopping.## Finished MCMLE.## Note: Null model likelihood calculation is not implemented for valued
## ERGMs at this time. This means that all likelihood-based inference
## (LRT, Analysis of Deviance, AIC, BIC, etc.) is only valid between
## models with the same reference distribution and constraints.## 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.mcmc.diagnostics(samplk.tot.ergm)## Sample statistics summary:
##
## Iterations = 16384:4209664
## Thinning interval = 1024
## Number of chains = 1
## Sample size per chain = 4096
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## sum 4.401 19.604 0.3063 0.6471
## nonzero 2.964 9.294 0.1452 0.3168
## mutual.min 1.454 9.007 0.1407 0.2795
## transitiveweights.min.max.min 6.930 24.848 0.3882 0.7882
## cyclicalweights.min.max.min 7.137 26.257 0.4103 0.8250
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## sum -33.00 -9.00 4 18 43.00
## nonzero -15.00 -3.00 3 9 21.00
## mutual.min -15.00 -5.00 1 8 19.62
## transitiveweights.min.max.min -39.62 -10.25 6 24 58.00
## cyclicalweights.min.max.min -41.00 -11.00 6 25 60.00
##
##
## Sample statistics cross-correlations:
## sum nonzero mutual.min
## sum 1.0000000 0.9245203 0.8626262
## nonzero 0.9245203 1.0000000 0.7779440
## mutual.min 0.8626262 0.7779440 1.0000000
## transitiveweights.min.max.min 0.9089147 0.8969098 0.7648867
## cyclicalweights.min.max.min 0.8847698 0.8791953 0.7973502
## transitiveweights.min.max.min
## sum 0.9089147
## nonzero 0.8969098
## mutual.min 0.7648867
## transitiveweights.min.max.min 1.0000000
## cyclicalweights.min.max.min 0.9570465
## cyclicalweights.min.max.min
## sum 0.8847698
## nonzero 0.8791953
## mutual.min 0.7973502
## transitiveweights.min.max.min 0.9570465
## cyclicalweights.min.max.min 1.0000000
##
## Sample statistics auto-correlation:
## Chain 1
## sum nonzero mutual.min transitiveweights.min.max.min
## Lag 0 1.0000000 1.0000000 1.0000000 1.0000000
## Lag 1024 0.6162129 0.6059622 0.5052799 0.5722101
## Lag 2048 0.3975894 0.3984843 0.3016641 0.3657304
## Lag 3072 0.2680860 0.2770034 0.1917940 0.2394811
## Lag 4096 0.1764258 0.1854697 0.1305061 0.1559500
## Lag 5120 0.1210870 0.1212495 0.1040810 0.1135902
## cyclicalweights.min.max.min
## Lag 0 1.0000000
## Lag 1024 0.5719283
## Lag 2048 0.3593614
## Lag 3072 0.2338994
## Lag 4096 0.1596769
## Lag 5120 0.1197892
##
## Sample statistics burn-in diagnostic (Geweke):
## Chain 1
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## sum nonzero
## -0.28653 -0.04773
## mutual.min transitiveweights.min.max.min
## -0.42507 -0.18840
## cyclicalweights.min.max.min
## 0.07114
##
## Individual P-values (lower = worse):
## sum nonzero
## 0.7744760 0.9619316
## mutual.min transitiveweights.min.max.min
## 0.6707826 0.8505642
## cyclicalweights.min.max.min
## 0.9432878
## Joint P-value (lower = worse): 0.4988715 .
##
## 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).summary(samplk.tot.ergm)## Call:
## ergm(formula = samplk.tot ~ sum + nonzero + mutual("min") + transitiveweights("min",
## "max", "min") + cyclicalweights("min", "max", "min"), response = "nominations",
## reference = ~Binomial(3))
##
## Iterations: 2 out of 20
##
## Monte Carlo MLE Results:
## Estimate Std. Error MCMC % z value Pr(>|z|)
## sum -0.04365 0.18772 0 -0.233 0.8161
## nonzero -3.47124 0.30637 0 -11.330 <1e-04 ***
## mutual.min 1.38593 0.23726 0 5.841 <1e-04 ***
## transitiveweights.min.max.min 0.21782 0.16677 0 1.306 0.1915
## cyclicalweights.min.max.min -0.25025 0.14299 0 -1.750 0.0801 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Null Deviance: 0.0 on 306 degrees of freedom
## Residual Deviance: -609.4 on 301 degrees of freedom
##
## Note that the null model likelihood and deviance are defined to be 0.
## This means that all likelihood-based inference (LRT, Analysis of
## Deviance, AIC, BIC, etc.) is only valid between models with the same
## reference distribution and constraints.
##
## AIC: -599.4 BIC: -580.8 (Smaller is better.)What does it tell us? The negative coefficient on nonzero suggests zero-inflation, there is strong evidence of mutuality, and the positive coefficient on transitive weights and negative on the cyclical weights suggests hierarchy, but they are not significant.
Example2: Zachary’s Karate Club Now let’s try using Poisson to model the Zachary Karate Club data: a zero-modified Poisson, with potentially different levels of activity for the faction leaders, heterogeneity in actor activity level overall, and an effect of difference in faction membership, a model that looks like this:
summary(zach ~ sum + nonzero + nodefactor("leader") + absdiffcat("faction.id") + nodesqrtcovar(TRUE), response="contexts")## sum nonzero
## 231.00000 78.00000
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 90.00000 74.00000
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 12.00000 11.00000
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 10.00000 16.17248A few other notes: - By default, Poisson-reference ERGMs use a variant of the TNT proposal to make sampling from spares (or 0-inflated) networks more efficient. There is a tuning parameter MCMC.prop.args = lsit(p0 = ...) that can be used to control how much “zero-inflation” there is. It has a sensible default. - Informally, a Poisson random variable contains more “information” than a Bernoulli random variable, which means that large changes in likelihood are not necessarily symptomatic of a problem. Thus, it often helps to set MCMLE.trustregion, which is normally 20 to something higher.
Now, for the fit and the diagnostics:
zach.pois <- ergm(zach ~ sum + nonzero + nodefactor("leader") + absdiffcat("faction.id") + nodesqrtcovar(TRUE) , response="contexts", reference=~Poisson, control=control.ergm(MCMLE.trustregion=100, MCMLE.maxit=50), verbose=TRUE)## Evaluating network in model.## Initializing Metropolis-Hastings proposal(s): ergm.count:MH_PoissonTNT
## Initializing model.
## Using initial method 'CD'.
## Fitting initial model.
## Starting contrastive divergence estimation via CD-MCMLE:
## Iteration 1 of at most 60:
## Convergence test P-value:0e+00
## Optimizing with step length 0.611245110020712.
## The log-likelihood improved by 1.812.
## Iteration 2 of at most 60:
## Convergence test P-value:0e+00
## Optimizing with step length 0.912862030647809.
## The log-likelihood improved by 1.763.
## Iteration 3 of at most 60:
## Convergence test P-value:7.9e-176
## Optimizing with step length 1.
## The log-likelihood improved by 0.6215.
## Iteration 4 of at most 60:
## Convergence test P-value:9.3e-54
## Optimizing with step length 1.
## The log-likelihood improved by 0.1654.
## Iteration 5 of at most 60:
## Convergence test P-value:7.3e-13
## Optimizing with step length 1.
## The log-likelihood improved by 0.03799.
## Iteration 6 of at most 60:
## Convergence test P-value:7.9e-04
## Optimizing with step length 1.
## The log-likelihood improved by 0.01338.
## Iteration 7 of at most 60:
## Convergence test P-value:4.5e-01
## Optimizing with step length 1.
## The log-likelihood improved by 0.003854.
## Iteration 8 of at most 60:
## Convergence test P-value:5.1e-01
## Convergence detected. Stopping.
## Optimizing with step length 1.
## The log-likelihood improved by 0.003572.
## Starting Monte Carlo maximum likelihood estimation (MCMLE):
## Density guard set to 10000 from an initial count of 78 edges.
##
## Iteration 1 of at most 50 with parameter:
## sum nonzero
## 1.06012423 -4.32173053
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -0.04846859 0.12640124
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.25332102 -0.91051794
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.92364092 2.80787170
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 385.38477 19.03516
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -54.71777 210.08398
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 191.00195 52.40332
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -8.52832 195.62193
## Starting MCMLE Optimization...
## Optimizing with step length 0.123466198773696.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.853.
##
## Iteration 2 of at most 50 with parameter:
## sum nonzero
## 1.19403083 -4.47310634
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.22993012 0.01231704
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.47709187 -0.88066711
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.75014503 2.65821337
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 600.54590 54.44824
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 176.80078 251.67676
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 175.31445 90.75977
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 43.17578 255.69528
## Starting MCMLE Optimization...
## Optimizing with step length 0.102912959419001.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 4.488.
##
## Iteration 3 of at most 50 with parameter:
## sum nonzero
## 1.35942885 -5.31423039
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.27565686 -0.07318498
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.53162920 -0.84130670
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.70925515 2.17682473
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 461.75977 29.87109
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 219.91504 146.83301
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 149.10059 56.93652
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 57.77734 211.95933
## Starting MCMLE Optimization...
## Optimizing with step length 0.12468918335365.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 2.979.
##
## Iteration 4 of at most 50 with parameter:
## sum nonzero
## 1.4402976 -5.7371195
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.2330784 -0.1278777
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.5733584 -0.7717300
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.7576397 1.9407619
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 342.243164 5.242188
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 391.564453 47.046875
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 34.724609 49.681641
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 108.610352 176.148398
## Starting MCMLE Optimization...
## Optimizing with step length 0.166335874696315.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.104.
##
## Iteration 5 of at most 50 with parameter:
## sum nonzero
## 1.4112998 -5.6793203
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.2670265 -0.1009915
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.5083182 -0.7720514
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.7580388 1.6883462
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 295.680664 4.518555
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 384.864258 27.591797
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 29.171875 37.333984
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 106.858398 154.430696
## Starting MCMLE Optimization...
## Optimizing with step length 0.132270004432296.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 2.839.
##
## Iteration 6 of at most 50 with parameter:
## sum nonzero
## 1.55987805 -6.12157693
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -0.07794834 -0.09988690
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.49691924 -0.78497569
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.77465272 2.01080927
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 256.765625 3.730469
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 330.493164 22.414062
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 25.517578 30.867188
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 91.029297 131.481853
## Starting MCMLE Optimization...
## Optimizing with step length 0.21735141208215.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 4.973.
##
## Iteration 7 of at most 50 with parameter:
## sum nonzero
## 1.65577600 -6.41753624
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -0.54081664 -0.07725244
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.50213467 -0.80217877
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.77334251 2.64766295
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 683.68652 22.39941
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -46.58301 189.90332
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 72.30957 126.92188
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 151.65137 311.92353
## Starting MCMLE Optimization...
## Optimizing with step length 0.186903977552815.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.073.
##
## Iteration 8 of at most 50 with parameter:
## sum nonzero
## 1.57937140 -6.55095705
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -0.63408047 -0.06479942
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.43788335 -0.73159318
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.73055397 2.63337716
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 643.19141 21.34473
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -54.26367 173.92090
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 73.57129 126.69727
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 150.00586 303.74648
## Starting MCMLE Optimization...
## Optimizing with step length 0.0484736368297892.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.1.
##
## Iteration 9 of at most 50 with parameter:
## sum nonzero
## 1.80617722 -8.19341297
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -0.38501058 -0.06340425
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.42218706 -0.70026865
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.69702936 1.85773727
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 433.422852 -1.510742
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -62.356445 138.925781
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 53.293945 100.631836
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 82.903320 223.620270
## Starting MCMLE Optimization...
## Optimizing with step length 0.110713313522399.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.567.
##
## Iteration 10 of at most 50 with parameter:
## sum nonzero
## 1.80117305 -7.95523930
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -0.02259040 -0.09675795
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.41829480 -0.70758007
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.68072925 1.69902446
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 362.204102 -4.994141
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -51.679688 120.494141
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 44.367188 88.219727
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 63.299805 189.243904
## Starting MCMLE Optimization...
## Optimizing with step length 0.149768434307909.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.704.
##
## Iteration 11 of at most 50 with parameter:
## sum nonzero
## 1.7438984 -7.6258406
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.1286752 -0.1174595
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.4789251 -0.7379670
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.7137420 1.6498629
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 232.72949 -18.83008
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 268.76074 30.87500
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 23.32910 29.60742
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 82.68457 130.83317
## Starting MCMLE Optimization...
## Optimizing with step length 0.171442753563562.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.245.
##
## Iteration 12 of at most 50 with parameter:
## sum nonzero
## 1.6103447 -6.7220170
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.1482553 -0.1306844
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.4876581 -0.7443429
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.7428642 1.7125008
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 329.5302734 -0.9189453
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 427.4765625 26.4785156
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 32.0605469 41.9335938
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 116.5410156 172.8196808
## Starting MCMLE Optimization...
## Optimizing with step length 0.144523245045526.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.779.
##
## Iteration 13 of at most 50 with parameter:
## sum nonzero
## 1.77818875 -7.29659492
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -0.26009942 -0.09910257
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.49853506 -0.72136789
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.78613173 2.11768798
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 438.7919922 -0.1865234
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 5.8417969 139.2539062
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 53.5224609 87.6699219
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 60.3330078 211.5976443
## Starting MCMLE Optimization...
## Optimizing with step length 0.171415483049749.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 2.679.
##
## Iteration 14 of at most 50 with parameter:
## sum nonzero
## 1.6922854 -6.7687628
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -0.2946200 -0.1259875
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.4761426 -0.7118592
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.7778209 2.1552663
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 405.226562 -1.329102
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -59.910156 134.635742
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 55.925781 93.069336
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 55.954102 210.706174
## Starting MCMLE Optimization...
## Optimizing with step length 0.0917683931673208.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.542.
##
## Iteration 15 of at most 50 with parameter:
## sum nonzero
## 1.77807683 -7.49683975
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.04927294 -0.10252732
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.41165866 -0.68827915
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.74049784 1.66463155
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 362.8964844 -0.6298828
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -43.8632812 130.3867188
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 50.6367188 77.0878906
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 37.6435547 171.9432248
## Starting MCMLE Optimization...
## Optimizing with step length 0.163908095127206.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 2.899.
##
## Iteration 16 of at most 50 with parameter:
## sum nonzero
## 1.6997336 -7.0718126
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.1300363 -0.1329366
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.4201784 -0.6966249
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.6871549 1.6640456
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 320.459961 -2.359375
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 147.162109 106.303711
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 44.612305 66.678711
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 39.216797 154.129234
## Starting MCMLE Optimization...
## Optimizing with step length 0.196819638256654.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 2.857.
##
## Iteration 17 of at most 50 with parameter:
## sum nonzero
## 1.5713756 -6.4793892
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.1568670 -0.1203021
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.3997759 -0.6928432
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.6953927 1.7118526
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 322.611328 3.064453
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 210.080078 96.754883
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 45.650391 66.444336
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 55.393555 163.609808
## Starting MCMLE Optimization...
## Optimizing with step length 0.189417420719721.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 2.953.
##
## Iteration 18 of at most 50 with parameter:
## sum nonzero
## 1.5034226 -6.1638696
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.1807291 -0.1281155
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.3992656 -0.6831119
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.7074475 1.6233173
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 296.737305 3.243164
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 379.260742 24.503906
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 33.181641 40.063477
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 104.834961 153.483584
## Starting MCMLE Optimization...
## Optimizing with step length 0.156534088033794.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.291.
##
## Iteration 19 of at most 50 with parameter:
## sum nonzero
## 1.53074763 -6.30989732
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.05413629 -0.15449921
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.45344095 -0.71956733
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.70168966 1.63390645
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 91.08887 -20.04199
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 144.99707 -14.85156
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 11.91113 14.03809
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 42.04883 65.50605
## Starting MCMLE Optimization...
## Optimizing with step length 0.272052113471358.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.944.
##
## Iteration 20 of at most 50 with parameter:
## sum nonzero
## 1.37723070 -5.43985841
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.03106387 -0.11182222
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.51700455 -0.73102228
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.77142374 1.79052037
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 76.962891 -11.710938
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 104.120117 -4.709961
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 9.140625 12.530273
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 32.552734 51.304453
## Starting MCMLE Optimization...
## Optimizing with step length 0.351108385516244.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.018.
##
## Iteration 21 of at most 50 with parameter:
## sum nonzero
## 1.25833295 -4.85880515
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.02259589 -0.11560444
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.54674375 -0.76460069
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -0.86435539 1.85578366
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 51.487305 -5.799805
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -52.432617 35.960938
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 8.758789 14.704102
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -2.586914 26.883003
## Starting MCMLE Optimization...
## Optimizing with step length 0.570726995740676.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.121.
##
## Iteration 22 of at most 50 with parameter:
## sum nonzero
## 1.12060584 -4.32999483
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.09289602 -0.09093094
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.58098183 -0.76362331
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -1.05256653 1.92338704
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 46.1669922 0.9267578
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -43.8623047 42.7138672
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 8.0947266 15.2792969
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -4.6787109 24.8956024
## Starting MCMLE Optimization...
## Optimizing with step length 0.934030023608374.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.351.
##
## Iteration 23 of at most 50 with parameter:
## sum nonzero
## 1.01594097 -4.07912931
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.19336055 -0.08206791
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.55729528 -0.73659932
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -1.02676082 1.96231215
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 57.8554688 3.2744141
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 99.7001953 -0.5341797
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 9.2568359 14.0888672
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 22.3095703 38.6012653
## Starting MCMLE Optimization...
## Optimizing with step length 0.981593402854804.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 3.958.
##
## Iteration 24 of at most 50 with parameter:
## sum nonzero
## 1.01621914 -3.96506330
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.18347192 -0.05972676
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.63857878 -0.81021347
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -1.01942338 1.81128260
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 11.305664 1.420898
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 8.061523 6.836914
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 1.078125 2.161133
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 3.420898 3.555497
## Starting MCMLE Optimization...
## Optimizing with step length 1.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## The log-likelihood improved by 0.3377.
## Step length converged once. Increasing MCMC sample size.
##
## Iteration 25 of at most 50 with parameter:
## sum nonzero
## 0.99773620 -3.86261765
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 0.19728166 -0.09423273
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.68070933 -0.88068915
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -1.13696500 1.84288723
## Starting unconstrained MCMC...
## Back from unconstrained MCMC.
## Average estimating function values:
## sum nonzero
## 6.83471680 1.06030273
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 4.72729492 1.92236328
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -0.01513672 0.57275391
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 0.79833984 2.72164487
## Starting MCMLE Optimization...
## Optimizing with step length 1.
## Using lognormal metric (see control.ergm function).
## Using log-normal approx (no optim)
## Starting MCMC s.e. computation.
## The log-likelihood improved by 0.05352.
## Step length converged twice. Stopping.
## Finished MCMLE.
## Note: Null model likelihood calculation is not implemented for valued
## ERGMs at this time. This means that all likelihood-based inference
## (LRT, Analysis of Deviance, AIC, BIC, etc.) is only valid between
## models with the same reference distribution and constraints.
## Evaluating log-likelihood at the estimate. Using 20 bridges: Running theta=[ 0.023687434,-0.091923294, 0.004838928,-0.002290386,-0.015798061,-0.020821652,-0.026972007, 0.042825189].
## Running theta=[ 0.073555715,-0.285446019, 0.015026146,-0.007112252,-0.049057135,-0.064656710,-0.083755181, 0.132983482].
## Running theta=[ 0.12342400,-0.47896874, 0.02521336,-0.01193412,-0.08231621,-0.10849177,-0.14053835, 0.22314177].
## Running theta=[ 0.17329228,-0.67249147, 0.03540058,-0.01675598,-0.11557528,-0.15232683,-0.19732153, 0.31330007].
## Running theta=[ 0.22316056,-0.86601419, 0.04558780,-0.02157785,-0.14883436,-0.19616188,-0.25410470, 0.40345836].
## Running theta=[ 0.27302884,-1.05953692, 0.05577502,-0.02639972,-0.18209343,-0.23999694,-0.31088787, 0.49361665].
## Running theta=[ 0.32289712,-1.25305964, 0.06596223,-0.03122158,-0.21535251,-0.28383200,-0.36767105, 0.58377495].
## Running theta=[ 0.37276540,-1.44658237, 0.07614945,-0.03604345,-0.24861158,-0.32766706,-0.42445422, 0.67393324].
## Running theta=[ 0.42263368,-1.64010509, 0.08633667,-0.04086531,-0.28187066,-0.37150211,-0.48123739, 0.76409153].
## Running theta=[ 0.47250197,-1.83362781, 0.09652388,-0.04568718,-0.31512973,-0.41533717,-0.53802057, 0.85424983].
## Running theta=[ 0.52237025,-2.02715054, 0.10671110,-0.05050905,-0.34838881,-0.45917223,-0.59480374, 0.94440812].
## Running theta=[ 0.57223853,-2.22067326, 0.11689832,-0.05533091,-0.38164788,-0.50300729,-0.65158691, 1.03456641].
## Running theta=[ 0.62210681,-2.41419599, 0.12708554,-0.06015278,-0.41490696,-0.54684234,-0.70837009, 1.12472470].
## Running theta=[ 0.67197509,-2.60771871, 0.13727275,-0.06497465,-0.44816603,-0.59067740,-0.76515326, 1.21488300].
## Running theta=[ 0.72184337,-2.80124144, 0.14745997,-0.06979651,-0.48142511,-0.63451246,-0.82193643, 1.30504129].
## Running theta=[ 0.77171165,-2.99476416, 0.15764719,-0.07461838,-0.51468418,-0.67834752,-0.87871961, 1.39519958].
## Running theta=[ 0.82157993,-3.18828688, 0.16783441,-0.07944024,-0.54794326,-0.72218258,-0.93550278, 1.48535788].
## Running theta=[ 0.87144822,-3.38180961, 0.17802162,-0.08426211,-0.58120233,-0.76601763,-0.99228595, 1.57551617].
## Running theta=[ 0.92131650,-3.57533233, 0.18820884,-0.08908398,-0.61446141,-0.80985269,-1.04906913, 1.66567446].
## Running theta=[ 0.97118478,-3.76885506, 0.19839606,-0.09390584,-0.64772048,-0.85368775,-1.10585230, 1.75583275].
## .
## This model was fit using MCMC. To examine model diagnostics and check
## for degeneracy, use the mcmc.diagnostics() function.mcmc.diagnostics(zach.pois)## Sample statistics summary:
##
## Iterations = 16384:4209664
## Thinning interval = 1024
## Number of chains = 1
## Sample size per chain = 4096
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## sum 6.83472 30.850 0.48202 1.8652
## nonzero 1.06030 8.089 0.12640 0.3264
## nodefactor.sum.leader.TRUE 4.72729 29.862 0.46659 2.7825
## absdiff.sum.faction.id.1 1.92236 16.519 0.25811 0.6993
## absdiff.sum.faction.id.2 -0.01514 5.610 0.08766 0.1889
## absdiff.sum.faction.id.3 0.57275 5.655 0.08835 0.2264
## absdiff.sum.faction.id.4 0.79834 6.139 0.09593 0.3549
## nodesqrtcovar.centered 2.72164 9.757 0.15246 0.8575
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## sum -52.00 -14.000 7.000 28.000 69.00
## nonzero -15.00 -4.000 1.000 7.000 17.00
## nodefactor.sum.leader.TRUE -48.00 -17.000 3.000 25.000 65.00
## absdiff.sum.faction.id.1 -27.00 -10.000 1.000 12.000 37.00
## absdiff.sum.faction.id.2 -9.00 -4.000 -1.000 4.000 12.00
## absdiff.sum.faction.id.3 -8.00 -4.000 0.000 4.000 13.00
## absdiff.sum.faction.id.4 -8.00 -4.000 0.000 4.000 15.00
## nodesqrtcovar.centered -12.94 -4.785 2.291 8.777 24.21
##
##
## Sample statistics cross-correlations:
## sum nonzero nodefactor.sum.leader.TRUE
## sum 1.0000000 0.8132389 0.56059956
## nonzero 0.8132389 1.0000000 0.29577039
## nodefactor.sum.leader.TRUE 0.5605996 0.2957704 1.00000000
## absdiff.sum.faction.id.1 0.5353617 0.4988607 0.03979358
## absdiff.sum.faction.id.2 0.2883610 0.2827363 0.30779395
## absdiff.sum.faction.id.3 0.3658751 0.3142282 0.45563951
## absdiff.sum.faction.id.4 0.3794053 0.2893082 0.64963244
## nodesqrtcovar.centered 0.6438022 0.3617255 0.81515711
## absdiff.sum.faction.id.1 absdiff.sum.faction.id.2
## sum 0.53536167 0.28836104
## nonzero 0.49886069 0.28273625
## nodefactor.sum.leader.TRUE 0.03979358 0.30779395
## absdiff.sum.faction.id.1 1.00000000 0.06793945
## absdiff.sum.faction.id.2 0.06793945 1.00000000
## absdiff.sum.faction.id.3 0.06384909 0.16617840
## absdiff.sum.faction.id.4 -0.08537683 0.19770516
## nodesqrtcovar.centered 0.13380776 0.31472676
## absdiff.sum.faction.id.3 absdiff.sum.faction.id.4
## sum 0.36587515 0.37940531
## nonzero 0.31422821 0.28930816
## nodefactor.sum.leader.TRUE 0.45563951 0.64963244
## absdiff.sum.faction.id.1 0.06384909 -0.08537683
## absdiff.sum.faction.id.2 0.16617840 0.19770516
## absdiff.sum.faction.id.3 1.00000000 0.28610114
## absdiff.sum.faction.id.4 0.28610114 1.00000000
## nodesqrtcovar.centered 0.48850090 0.63390863
## nodesqrtcovar.centered
## sum 0.6438022
## nonzero 0.3617255
## nodefactor.sum.leader.TRUE 0.8151571
## absdiff.sum.faction.id.1 0.1338078
## absdiff.sum.faction.id.2 0.3147268
## absdiff.sum.faction.id.3 0.4885009
## absdiff.sum.faction.id.4 0.6339086
## nodesqrtcovar.centered 1.0000000
##
## Sample statistics auto-correlation:
## Chain 1
## sum nonzero nodefactor.sum.leader.TRUE
## Lag 0 1.0000000 1.0000000 1.0000000
## Lag 1024 0.7275382 0.3872732 0.9208029
## Lag 2048 0.5941120 0.2629616 0.8605541
## Lag 3072 0.5227698 0.2118173 0.8100151
## Lag 4096 0.4642496 0.1733856 0.7645480
## Lag 5120 0.4223528 0.1771779 0.7267307
## absdiff.sum.faction.id.1 absdiff.sum.faction.id.2
## Lag 0 1.0000000 1.0000000
## Lag 1024 0.6758913 0.3737886
## Lag 2048 0.5277069 0.2478362
## Lag 3072 0.4215962 0.1869893
## Lag 4096 0.3335231 0.1671493
## Lag 5120 0.2629816 0.1344209
## absdiff.sum.faction.id.3 absdiff.sum.faction.id.4
## Lag 0 1.0000000 1.0000000
## Lag 1024 0.4509826 0.5569002
## Lag 2048 0.3379687 0.4621891
## Lag 3072 0.2939022 0.3750311
## Lag 4096 0.2482221 0.3271915
## Lag 5120 0.2305064 0.2852206
## nodesqrtcovar.centered
## Lag 0 1.0000000
## Lag 1024 0.8845695
## Lag 2048 0.8190385
## Lag 3072 0.7628534
## Lag 4096 0.7160669
## Lag 5120 0.6777946
##
## Sample statistics burn-in diagnostic (Geweke):
## Chain 1
##
## Fraction in 1st window = 0.1
## Fraction in 2nd window = 0.5
##
## sum nonzero
## -2.1699 -1.9417
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## -2.9571 0.4456
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## -1.0857 -4.2015
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## -3.7613 -1.7170
##
## Individual P-values (lower = worse):
## sum nonzero
## 3.001081e-02 5.217428e-02
## nodefactor.sum.leader.TRUE absdiff.sum.faction.id.1
## 3.105329e-03 6.558839e-01
## absdiff.sum.faction.id.2 absdiff.sum.faction.id.3
## 2.776086e-01 2.651983e-05
## absdiff.sum.faction.id.4 nodesqrtcovar.centered
## 1.690290e-04 8.597738e-02
## Joint P-value (lower = worse): 0.1839626 .
##
## 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).summary(zach.pois)## Call:
## ergm(formula = zach ~ sum + nonzero + nodefactor("leader") +
## absdiffcat("faction.id") + nodesqrtcovar(TRUE), response = "contexts",
## reference = ~Poisson, control = control.ergm(MCMLE.trustregion = 100,
## MCMLE.maxit = 50), verbose = TRUE)
##
## Iterations: 25 out of 50
##
## Monte Carlo MLE Results:
## Estimate Std. Error MCMC % z value Pr(>|z|)
## sum 0.99487 0.08267 0 12.034 < 1e-04 ***
## nonzero -3.86078 0.24887 0 -15.513 < 1e-04 ***
## nodefactor.sum.leader.TRUE 0.20323 0.06394 0 3.178 0.00148 **
## absdiff.sum.faction.id.1 -0.09620 0.07937 0 -1.212 0.22552
## absdiff.sum.faction.id.2 -0.66352 0.19386 0 -3.423 0.00062 ***
## absdiff.sum.faction.id.3 -0.87451 0.21102 0 -4.144 < 1e-04 ***
## absdiff.sum.faction.id.4 -1.13282 0.23646 0 -4.791 < 1e-04 ***
## nodesqrtcovar.centered 1.79866 0.21340 0 8.428 < 1e-04 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Null Deviance: 0.0 on 561 degrees of freedom
## Residual Deviance: -837.9 on 553 degrees of freedom
##
## Note that the null model likelihood and deviance are defined to be 0.
## This means that all likelihood-based inference (LRT, Analysis of
## Deviance, AIC, BIC, etc.) is only valid between models with the same
## reference distribution and constraints.
##
## AIC: -821.9 BIC: -787.3 (Smaller is better.)What does it tell us? The negative coefficient on nonzero suggests zero-inflation, the faction leaders clearly have more activity than others, and the more ideologically separated two individuals are, the less they interact. Over and above that, there is some additional heterogeneity in how active individuals are: if i has a lot of interaction with j, it is likely that i has more with j’. Could this mean a form of preferential attachment?
We can try seeing whether there is some friend of a friend effect above and beyond that. This can be done by fitting a model with transitivity and seeing whether the coefficient is significant, or we can perform a simulation test. In the following
simulateunpacks thezach.poisERGM fit, extracting the formula, the coefficient, and the rest of the information.nsimsays how many networks to generateoutput="stats"says that we only want to see the simulated statistics, not the networksmonitor=~transitiveweights("geomean", "sum", "geomean")says that in addition to the statistics used in the fit, we wantsimulateto keep track of the transitive weights statistic.
We do not need to worry about degeneracy in this case, because we are not actually using that statistic in the model, only “monitoring” it.
# Simulate from model fit:
zach.sim <- simulate(zach.pois, monitor=~transitiveweights("geomean", "sum", "geomean"), nsim=1000, output="stats")# What have we simulated?
(colnames(zach.sim))## [1] "sum"
## [2] "nonzero"
## [3] "nodefactor.sum.leader.TRUE"
## [4] "absdiff.sum.faction.id.1"
## [5] "absdiff.sum.faction.id.2"
## [6] "absdiff.sum.faction.id.3"
## [7] "absdiff.sum.faction.id.4"
## [8] "nodesqrtcovar.centered"
## [9] "transitiveweights.geomean.sum.geomean"# How high is the transitiveweights statistic in the observed network?
zach.obs <- summary(zach ~ transitiveweights("geomean", "sum", "geomean"), response="contexts")
zach.obs## transitiveweights.geomean.sum.geomean
## 288.9793Let’s plot the density of the simulated values of transitive weights statistic:
par(mar = c(5, 4, 4, 2) + 0.1)
# 9th col. = transitiveweights
plot(density(zach.sim[, 9]))
abline(v = zach.obs)
# Where does the observed value lie in the simualted? This is a p-value for the
# Monte-Carlo test:
min(mean(zach.sim[, 9] > zach.obs), mean(zach.sim[, 9] < zach.obs)) * 2## [1] 0.594Looks like individual heterogeneity and faction alignment account for appearance of triadic effects. (Notably, the factions themselves may be endogenous, if social influence is a factor. Untangling selection from influence is hard enough when dynamic network data are available. We cannot do it here.)