A Code

All source materials used to write this thesis are available in the Reed College Library’s digital thesis archive. This appendix contains the most important code used in the final chapter.

A.1 Network simulation code

if(!require(devtools))
  install.packages("devtools", repos = "http://cran.rstudio.com")
if(!require(thesisdown))
  devtools::install_github("ismayc/thesisdown")
library(thesisdown)

library(statnet)
library(Bergm)
library(tidyverse)
library(broom)
library(ggraph)
library(purrr)
library(Rcpp)
library(knitr)

The following C++ code, compiled into a shared library using the Rcpp package, gives function definitions used in simulating networks (Eddelbuettel & Francois, 2011).

#include <Rcpp.h>

using namespace Rcpp;

//[[Rcpp::export]]
IntegerMatrix make_adj(List vert_df){
  NumericVector vert = vert_df["vert"];
  NumericVector nbhd = vert_df["nbhd"];
  NumericVector grp = vert_df["grp"];
  double n = vert.size();
  
  IntegerMatrix adj_mat(n,n);
  
  for(int i = 0; i < n; i++){
    for(int j = 0; j < n; j++){
      if(nbhd[i] == nbhd[j]){
        if(grp[i] == grp[j]){
          adj_mat(i,j) = as<int>(rbinom(1,1,.3));
          continue;
        } else {
          adj_mat(i,j) = as<int>(rbinom(1, 1, .1));
          continue;
          }
      } else{       
        adj_mat(i,j) = as<int>(rbinom(1, 1, 50/(n*n)));
        continue;
      }
    }
  }

  for(int i = 1; i < n/10; i++){
    for(int j = 0; j < n; j++){
      for(int k = 0; k < n; k++){
        for(int l = 0; l < n; l++){
          if(nbhd[j] == i &&
             nbhd[k] == i &&
             nbhd[l] == i){
            if(adj_mat(j, k) == 1 &&
               adj_mat(k, l) == 1 &&
               adj_mat(j, l) == 0){
              adj_mat(j, l) = as<int>(rbinom(1, 1, .75));
            }
          }
        }
      }
    }
  }

  return adj_mat;
}

//[[Rcpp::export]]
IntegerVector make_attr(List vert_df, IntegerMatrix adj_mat){
  IntegerVector attr = vert_df["attr"];
  IntegerVector nbhd = vert_df["nbhd"];
  int n = attr.size();
  for(int i = 0; i < n; i++){
    for(int j = 0; j < n; j++){
      if(attr[i] == 1 && adj_mat(i,j) == 1 && nbhd[i] == nbhd[j]){
        attr[j] = as<int>(rbinom(1, 1, .75));
      }else if(adj_mat(i,j) == 1 && nbhd[i] == nbhd[j]){
        attr[j] = as<int>(rbinom(1, 1, .1));
      }
    }
  }
  return attr;
}

This R code generates networks and computes the three chosen statistics for each one. The final result is a data set of networks and corresponding statistics used to produce density and quantile-quantile plots.

make_net <- function(n){
  vert_df <- tibble(vert = 1:n)
  nbhds <- 1:floor(n/10)
  vert_df <- vert_df %>% 
    mutate(nbhd = sample(nbhds, n, replace = TRUE),
           grp = rbinom(n, 1, .5),
           attr = rbinom(n, 1, .3))
  
  adj_mat <- make_adj(vert_df)
  
  vert_df$attr <- make_attr(vert_df, adj_mat)
  
  net <- network(adj_mat)
  net %v% "group" <- vert_df$grp
  net %v% "attribute" <- vert_df$attr
  net %v% "nbhd" <- vert_df$nbhd
  
  return(net)
}

node_match <- function(net){
  edge_list <- as.edgelist(net)
  grps <- net %v% "group"
  num <- 0
  for(i in 1:nrow(edge_list)){
    if(grps[edge_list[i, 1]] == grps[edge_list[i,2]]) num <- num + 1
  }
  num/nrow(edge_list)
}
n_reps <- 500
results <- tibble(num = c(rep(20, n_reps),
                        rep(50, n_reps),
                        rep(100, n_reps),
                        rep(200, n_reps)))
results <- results %>% 
  mutate(net = map(num, make_net),
         mean_deg = map_dbl(net, ~ mean(degree(.x))),
         mean_attr = map_dbl(net, ~ mean(.x %v% "attribute")),
         prop_homo = map_dbl(net, node_match))

results <- results %>% 
  gather(stat, value, -(num:net))

A.2 Jackknife standard error code

This code computes the neighborhood level standard error estimates for statistics of the graph.

n_vert <- 200
net <- make_net(n_vert)

nbhds <- net %v% "nbhd"
num_nbhds <- length(unique(nbhds))
boot3 <- numeric(num_nbhds)

jack_df <- tibble(mean_deg = numeric(num_nbhds),
                  mean_attr = numeric(num_nbhds),
                  prop_homo = numeric(num_nbhds))

for(i in 1:num_nbhds){
  temp_net <- rm_vertex(net, which(nbhds == i))
  jack_df$mean_deg[i] <- mean(degree(temp_net))
  jack_df$mean_attr[i] <- mean(temp_net %v% "attribute")
  jack_df$prop_homo[i] <- node_match(temp_net)
}

jack_df <- jack_df %>% 
  gather(stat, value) %>% 
  group_by(stat) %>% 
  summarise(jack_sd = sqrt((num_nbhds - 1)/num_nbhds * 
                             sum((value - mean(value))^2)))

References

Eddelbuettel, D., & Francois, R. (2011). Rcpp: Seamless R and C++ integration. Journal of Statistical Software, (40, 8).