Appendix A — Permissions and Code for Chapter 2
This chapter contains the code used for the analysis in chapter 2 and the ODD (Overview, Design concepts, Details) protocol for the agent-based model ArchMatNet. The code is organized into three sections: R, Functions, and R Environment. The R section contains the code for the analysis, the Functions section contains the functions used in the analysis, and the R Environment section contains the R environment used to run the analysis. The final section contains the ODD protocol which is an enhanced description of the agent-based model.
R Code
library(tidyverse)
library(here)
library(rio)
library(ggpubr)
library(igraph)
library(ggnetwork)
library(tictoc)
# library(patchwork)
library(magick)
library(BayesFactor)
i_am("Analysis/Chap2Analysis.R")
source("Analysis/functions.R")
vars <- c(
"uniqueTraits",
"nTicks",
"huntingGroupSize",
"pHunting",
"genderedMigration",
"pNewAlliance",
"maximumAlliesN",
"pTrading",
"visibility",
"pLoss",
"pLuckyLeap",
"pVisiting",
"learningMethod",
"gradualEnvChange",
"pEnvChange",
"aggregationDuration",
"transmissionRate",
"pMigration",
"nOfBands",
"campPopulation",
"potteryPrestige",
"pNewVariant",
"pNewObject",
"aggregationFreq"
)
#### load results ####
library(arrow)
runs <- read_feather(here("data/modelOutput.parquet"))
#### look at range of experiments ####
expDF <- runs %>%
select(expID) %>%
group_by(expID) %>%
count()
#### examine network comparisons ####
table_df <- tibble()
plotdf <- runs %>%
filter(uniqueTraits == 1) %>%
# filter(pTrading > 0) %>%
# filter(nOfBands > 3) %>%
select(contains("_corZscore")) %>%
pivot_longer(contains("_corZscore"), names_to = "Network Comparison", values_to = "Correlation Comparison") %>%
mutate(`Network Comparison` = gsub("_corZscore", "", `Network Comparison`)) %>%
filter(!str_detect(`Network Comparison`, "Hunting")) %>%
drop_na()
meanAll <- plotdf %>%
group_by(`Network Comparison`) %>%
summarize(mean = mean(`Correlation Comparison`, na.rm = T)) %>%
arrange(mean)
table_df <- meanAll %>%
mutate(metric = "Graph Correlation")
meanPotPoint <- plotdf %>%
filter(`Network Comparison` != "Point_Pot") %>%
filter(str_detect(`Network Comparison`, "Point|Pot")) %>%
mutate(`Network Comparison` = case_when(`Network Comparison` %>% str_detect("Point") ~ "Point", TRUE ~ "Pot")) %>%
group_by(`Network Comparison`) %>%
summarize(mean = mean(`Correlation Comparison`)) %>%
arrange(`Network Comparison`)
plotSum <- plotdf %>%
group_by(`Network Comparison`) %>%
summarize(mean = mean(`Correlation Comparison`, na.rm = T)) %>%
arrange(desc(mean))
p <- plotdf %>%
mutate(`Network Comparison` = fct_relevel(`Network Comparison`, plotSum$`Network Comparison`)) %>%
ggplot(aes(x = `Network Comparison`, y = `Correlation Comparison`, fill = `Network Comparison`)) +
geom_boxplot(alpha = .25, notch = T) +
guides(fill = "none") +
theme_bw() +
xlab("") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
p
# ggsave(here('figures/corComparisonZScoreBoxplot.png'),p,width = 7,height = 4.5, dpi = 600)
# ggsave(here('figures/corComparisonZScoreBoxplot.tiff'),p,width = 7,height = 4.5, compression = "lzw", dpi = 900)
plotdf <- runs %>%
filter(uniqueTraits == 1) %>%
select(contains("_QAP")) %>%
pivot_longer(contains("_QAP"), names_to = "Network Comparison", values_to = "QAP Test") %>%
mutate(`Network Comparison` = gsub("_QAP", "", `Network Comparison`)) %>%
filter(!str_detect(`Network Comparison`, "Hunting")) %>%
drop_na()
meanAll <- plotdf %>%
group_by(`Network Comparison`) %>%
summarize(mean = mean(`QAP Test`)) %>%
arrange(mean)
table_df <- bind_rows(table_df, meanAll %>%
mutate(metric = "QAP Test"))
meanPotPoint <- plotdf %>%
filter(`Network Comparison` != "Point_Pot") %>%
filter(str_detect(`Network Comparison`, "Point|Pot")) %>%
mutate(`Network Comparison` = case_when(`Network Comparison` %>% str_detect("Point") ~ "Point", TRUE ~ "Pot")) %>%
group_by(`Network Comparison`) %>%
summarize(mean = mean(`QAP Test`)) %>%
arrange(`Network Comparison`)
plotSum <- plotdf %>%
group_by(`Network Comparison`) %>%
summarize(mean = mean(`QAP Test`, na.rm = T)) %>%
arrange(mean)
p <- plotdf %>%
mutate(`Network Comparison` = fct_relevel(`Network Comparison`, plotSum$`Network Comparison`)) %>%
ggplot(aes(x = `Network Comparison`, y = `QAP Test`)) +
geom_point(position = position_jitter(), alpha = .1) +
theme_bw() +
xlab("") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
p
# ggsave(here('figures/QAPTest.tiff'),p,width = 7,height = 4.5, compression = "lzw", dpi = 900)
# browseURL(here('figures/QAPTest.tiff'))
# table_df %>%
# pivot_wider(names_from = metric, values_from = mean) %>%
# write_csv("tables/MetricResults.csv")
#### explore variables ####
experiments <- runs$expID %>%
unique() %>%
sort()
experiments
# unique traits
i <- 24
variables_df <- runs %>%
select(expID, any_of(vars), any_of(contains("_corZscore"))) %>%
pivot_longer(any_of(contains("_corZscore"))) %>%
filter(!str_detect(name, "Hunting")) %>%
filter(str_detect(name, "Pot|Point")) %>%
pivot_longer(any_of(vars), names_to = "Variable", values_to = "parameter") %>%
mutate(parameter = factor(parameter)) %>%
group_by(name, Variable, parameter) %>%
summarize(cor = mean(value, na.rm = T), count = n())
visiting_point_df <- variables_df %>%
filter(name == "Visiting_Point_corZscore") %>%
filter(Variable == vars[i])
visiting_point_df
runs %>%
select(uniqueTraits, nTicks, Visiting_Point_corZscore) %>%
mutate_at(vars(-Visiting_Point_corZscore), factor) %>%
ggplot(aes(x = Visiting_Point_corZscore, group = nTicks, color = nTicks)) +
geom_boxplot() +
theme_minimal() +
facet_wrap(~uniqueTraits, scales = "free_y") +
ggtitle("uniqueTraits")
variables_df <- runs %>%
filter(uniqueTraits == 1) %>%
select(-uniqueTraits) %>%
select(expID, any_of(vars), any_of(contains("_corZscore"))) %>%
pivot_longer(any_of(contains("_corZscore"))) %>%
filter(!str_detect(name, "Hunting")) %>%
filter(str_detect(name, "Pot|Point")) %>%
filter(name == "Visiting_Pot_corZscore") %>%
janitor::remove_constant()
library(mgcv)
dvars <- vars[which(vars %in% names(variables_df))]
fmla <- as.formula("value ~ nTicks + huntingGroupSize + pHunting +
factor(genderedMigration) + pNewAlliance + maximumAlliesN +
pTrading + factor(visibility) + pLoss + pLuckyLeap +
pVisiting + factor(learningMethod) + factor(gradualEnvChange) +
pEnvChange + aggregationDuration + transmissionRate +
pMigration + nOfBands + campPopulation +
factor(potteryPrestige) + pNewVariant + pNewObject +
aggregationFreq")
gam_mdl <- gam(fmla, data = variables_df)
summary(gam_mdl)
gam_summary <- summary(gam_mdl)
# Extract parametric and smooth terms
parametric_df <- as.data.frame(gam_summary$p.table)
smooth_df <- as.data.frame(gam_summary$s.table)
# Add type column and combine tables
parametric_df$type <- "Parametric"
# Format table
names(parametric_df) <- c("Estimate", "Std. Error", "t-value", "p-value", "Type")
parametric_df <- rownames_to_column(parametric_df, "Parameter")
parametric_df <- parametric_df %>%
mutate_if(is.numeric, round, 2)
# write_csv(parametric_df,here('tables/ParametricEstimates.csv'))
# variables that don't make a difference
xvars <- c("pNewAlliance", "maximumAlliesN", "pLoss", "pEnvChange", "aggregationDuration", "potteryPrestige", "pNewVariant", "transmissionRate", "pMigration", "nOfBands", "pNewVariant", "pLuckyLeap", "pHunting")
dvars <- vars[which(!vars %in% xvars)]
dvars <- dvars[which(dvars %in% names(variables_df))]
# fmla = as.formula(paste("value ~ ",paste(dvars,collapse = "+")))
fmla <- as.formula("value ~ nTicks +
huntingGroupSize +
factor(genderedMigration) +
pTrading + factor(visibility) +
pVisiting + factor(learningMethod) + factor(gradualEnvChange) + campPopulation + pNewObject +
aggregationFreq")
gam_mdl2 <- gam(fmla, data = variables_df)
summary(gam_mdl2)
xvars2 <- c("pNewAlliance", "maximumAlliesN", "pLoss", "pEnvChange", "aggregationDuration", "potteryPrestige", "pNewVariant", "transmissionRate", "pMigration", "nOfBands", "pNewVariant", "pLuckyLeap", "genderedMigration", "pHunting", "pTrading", "gradualEnvChange", "huntingGroupSize", "pNewObject", "pVisiting", "campPopulation")
dvars <- vars[which(!vars %in% xvars2)]
dvars <- dvars[which(dvars %in% names(variables_df))]
# fmla = as.formula(paste("value ~ ",paste("factor(",dvars,")",collapse = "+")))
fmla <- as.formula("value ~
nTicks +
factor(visibility) +
factor(learningMethod) +
aggregationFreq")
gam_mdl3 <- gam(fmla, data = variables_df)
summary(gam_mdl3)
gam_nTicks <- gam(value ~ factor(nTicks), data = variables_df)
gam_visibility <- gam(value ~ factor(visibility), data = variables_df)
gam_pVisiting <- gam(value ~ pVisiting, data = variables_df)
gam_learningMethod <- gam(value ~ factor(learningMethod), data = variables_df)
gam_aggregationFreq <- gam(value ~ aggregationFreq, data = variables_df)
gam_campPopulation <- gam(value ~ campPopulation, data = variables_df)
summary(gam_nTicks)
summary(gam_visibility)
summary(gam_pVisiting)
summary(gam_learningMethod)
summary(gam_aggregationFreq)
summary(gam_campPopulation)
variables_df %>%
select(-any_of(xvars2)) %>%
pivot_longer(any_of(vars), names_to = "Variable", values_to = "parameter") %>%
mutate(parameter = factor(parameter)) %>%
rename(cor = value) %>%
ggplot(aes(x = parameter, y = cor, color = Variable)) +
geom_boxplot() +
guides(color = "none") +
theme_minimal() +
facet_wrap(~Variable, scales = "free_x")
# ggsave(here('figures/VariableCorrelations.tiff'),width = 7,height = 4.5, compression = "lzw", dpi = 900)
# ggsave(here('figures/VariableCorrelations.png'),width = 7,height = 4.5, dpi = 150)
# browseURL(here('figures/VariableCorrelations.tiff'))
#### look at individual networks ####
tictoc::tic()
g <- runs %>%
filter(nOfBands == 5) %>%
filter(uniqueTraits == 0) %>%
filter(!is.na(wordMatrixToRowListCampVisitingAdjacency)) %>%
slice_sample(n = 100) %>%
mutate_all(as.character) %>%
select(-runNumber) %>%
rowid_to_column("runNumber") %>%
pivot_longer(-runNumber,
names_to = "property",
values_to = "value"
) %>%
filter(property %in% vars | str_detect(property, "wordMatrixToRowListCamp|runNumber")) %>%
mutate(property = str_remove_all(
property,
"wordMatrixToRowListCamp|Adjacency"
)) %>%
pivot_wider(names_from = property, values_from = value) %>%
pivot_longer(Visiting:Pot, names_to = "property", values_to = "adj") %>%
filter(!is.na(adj)) %>%
mutate(property = tolower(property)) %>%
# rowwise() %>%
mutate(adj = formatNLAdjTables(adj)) %>%
pivot_wider(names_from = property, values_from = adj)
# ungroup()
tictoc::toc()
bands <- g$point[1][[1]] %>%
select(name = from) %>%
distinct_all() %>%
bind_rows(g$point[1][[1]] %>%
select(name = to)) %>%
distinct_all() %>%
mutate(band = sort(rep(0:4, 3)))
getGraphDF <- function(network, i = 1) {
cutoff <- quantile(g[[network]][[i]]$weight, .5)
gdf <- g[[network]][[i]] %>%
filter(weight > cutoff)
n <- graph_from_data_frame(gdf, directed = F)
# summary(n)
plotdf <- ggnetwork(n) %>%
select(-band) %>%
left_join(bands) %>%
mutate(
band = factor(band),
network = !!network
)
return(plotdf)
}
variables <- c("visiting", "learning", "trading", "point", "pot")
combined_df <- map_df(variables, getGraphDF) %>%
mutate(network = factor(network, levels = variables))
p <- ggplot(combined_df, aes(x = x, y = y, xend = xend, yend = yend)) +
geom_edges(aes(color = weight, alpha = weight), linewidth = 1) +
geom_nodelabel(aes(fill = band, label = name), size = 3, color = "white") +
scale_fill_viridis_d() +
theme_void() +
coord_equal(ratio = 1 / 2, clip = "off", expand = T) +
theme(
legend.position = "none",
plot.background = element_rect(colour = "#404040", fill = NA, linewidth = .5),
plot.title = element_text(hjust = 0.5),
plot.margin = unit(c(2, 2, 2, 2), "mm"),
strip.placement = "outside",
strip.text = element_text(margin = margin(b = 10))
) +
facet_wrap(~network)
p
# ggsave("MaterialCultureNetworksPaper\\figures\\SameTraits.png",dpi = 600, width = 8, height = 4.5, bg = "white")
# fn = "MaterialCultureNetworksPaper\\figures\\SameTraits.tiff"
# ggsave(fn,dpi = 900, width = 7, height = 5, bg = "white", compression = "lzw")
# img = image_read(fn)
# img = image_trim(img)
# img = image_write(img,fn)
#
# image1 = image_read("MaterialCultureNetworksPaper\\figures\\SameTraits.tiff")
# image2 = image_read("MaterialCultureNetworksPaper\\figures\\UniqueTraits.tiff")
# image1_annotated <- image_annotate(image1, text = "a",
# gravity = "northwest",
# location = "+20+20",
# size = 200, color = "black")
#
# image2_annotated <- image_annotate(image2, text = "b",
# gravity = "northwest",
# location = "+20+20",
# size = 200, color = "black")
# stacked_images <- image_append(c(image1_annotated, image2_annotated), stack = TRUE)
# image_write(stacked_images, "MaterialCultureNetworksPaper\\figures\\combined_networks.tiff")
# g1 = getGraph("visiting",2)
# g2 = getGraph("learning",2)
# g3 = getGraph("trading",2)
# g4 = getGraph("point",2)
# g5 = getGraph("pot",2)
#
# g1 + g2 + g3 + g4 + g5
# ggsave("MaterialCultureNetworksPaper\\figures\\SameTraits2.png",dpi = 600, width = 8, height = 4.5, bg = "white")
g <- runs %>%
filter(nOfBands > 3) %>%
filter(uniqueTraits == 1) %>%
filter(!is.na(wordMatrixToRowListCampVisitingAdjacency)) %>%
slice_sample(n = 100) %>%
mutate_all(as.character) %>%
select(-runNumber) %>%
rowid_to_column("runNumber") %>%
pivot_longer(-runNumber,
names_to = "property",
values_to = "value"
) %>%
filter(property %in% vars | str_detect(property, "wordMatrixToRowListCamp|runNumber")) %>%
mutate(property = str_remove_all(
property,
"wordMatrixToRowListCamp|Adjacency"
)) %>%
pivot_wider(names_from = property, values_from = value) %>%
pivot_longer(Visiting:Pot, names_to = "property", values_to = "adj") %>%
filter(!is.na(adj)) %>%
mutate(property = tolower(property)) %>%
# rowwise() %>%
mutate(adj = formatNLAdjTables(adj)) %>%
pivot_wider(names_from = property, values_from = adj)
variables <- c("visiting", "learning", "trading", "point", "pot")
combined_df <- map_df(variables, getGraphDF) %>%
mutate(network = factor(network, levels = variables))
p <- ggplot(combined_df, aes(x = x, y = y, xend = xend, yend = yend)) +
geom_edges(aes(color = weight, alpha = weight), linewidth = 1) +
geom_nodelabel(aes(fill = band, label = name), size = 3, color = "white") +
scale_fill_viridis_d() +
theme_void() +
coord_equal(ratio = 1 / 2, clip = "off", expand = T) +
theme(
legend.position = "none",
plot.background = element_rect(colour = "#404040", fill = NA, linewidth = .5),
plot.title = element_text(hjust = 0.5),
plot.margin = unit(c(2, 2, 2, 2), "mm"),
strip.placement = "outside",
strip.text = element_text(margin = margin(b = 10))
) +
facet_wrap(~network)
p
# fn = "MaterialCultureNetworksPaper\\figures\\UniqueTraits.tiff"
# ggsave(fn,dpi = 900, width = 7, height = 5, bg = "white", compression = "lzw")
#
# img = image_read(fn)
# img = image_trim(img)
# img = image_write(img,fn)
#### centrality ####
largeNetworks <- runs %>%
filter(nOfBands > 3) %>%
filter(uniqueTraits == 1)
results_point_e <- list()
results_pot_e <- list()
for (i in 1:nrow(largeNetworks)) {
print(paste(i, "of", nrow(largeNetworks)))
x <- largeNetworks[i, ]
n1 <- matrix2zscore(x[["wordMatrixToRowListCampVisitingAdjacency"]])
diag(n1) <- 0
ec1 <- sna::evcent(n1)
n2 <- matrix2zscore(x[["wordMatrixToRowListCampPointAdjacency"]])
diag(n2) <- 0
ec2 <- sna::evcent(n2)
n3 <- matrix2zscore(x[["wordMatrixToRowListCampPotAdjacency"]])
diag(n3) <- 0
ec3 <- sna::evcent(n3)
r1 <- cor(ec2, ec1) %>% round(., 2)
r2 <- cor(ec3, ec1) %>% round(., 2)
results_point_e[[i]] <- r1
results_pot_e[[i]] <- r2
}
results_point_e <- unlist(results_point_e)
results_pot_e <- unlist(results_pot_e)
hist(results_point_e)
hist(results_pot_e)
results_point_e %>% summary()
results_pot_e %>% summary()
cor(largeNetworks$nOfBands, results_point_e)
plot(largeNetworks$nOfBands, results_point_e)
results_point_d <- list()
results_pot_d <- list()
for (i in 1:nrow(largeNetworks)) {
print(paste(i, "of", nrow(largeNetworks)))
x <- largeNetworks[i, ]
n1 <- matrix2zscore(x[["wordMatrixToRowListCampVisitingAdjacency"]])
diag(n1) <- 0
dc1 <- sna::degree(n1)
n2 <- matrix2zscore(x[["wordMatrixToRowListCampPointAdjacency"]])
diag(n2) <- 0
dc2 <- sna::degree(n2)
n3 <- matrix2zscore(x[["wordMatrixToRowListCampPotAdjacency"]])
diag(n3) <- 0
dc3 <- sna::degree(n3)
r1 <- cor(dc2, dc1) %>% round(., 2)
r2 <- cor(dc3, dc1) %>% round(., 2)
results_point_d[[i]] <- r1
results_pot_d[[i]] <- r2
}
results_point_d <- unlist(results_point_d)
results_pot_d <- unlist(results_pot_d)
hist(results_point_d)
hist(results_pot_d)
results_point_d %>% summary()
results_pot_d %>% summary()
results_point_b <- list()
results_pot_b <- list()
for (i in 1:nrow(largeNetworks)) {
print(paste(i, "of", nrow(largeNetworks)))
x <- largeNetworks[i, ]
n1 <- matrix2zscore(x[["wordMatrixToRowListCampVisitingAdjacency"]])
diag(n1) <- 0
bc1 <- sna::betweenness(n1)
n2 <- matrix2zscore(x[["wordMatrixToRowListCampPointAdjacency"]])
diag(n2) <- 0
bc2 <- sna::betweenness(n2)
n3 <- matrix2zscore(x[["wordMatrixToRowListCampPotAdjacency"]])
diag(n3) <- 0
bc3 <- sna::betweenness(n3)
r1 <- cor(bc2, bc1) %>% round(., 2)
r2 <- cor(bc3, bc1) %>% round(., 2)
results_point_b[[i]] <- r1
results_pot_b[[i]] <- r2
}
results_point_b <- unlist(results_point_b)
results_pot_b <- unlist(results_pot_b)
hist(results_point_b)
hist(results_pot_b)
results_point_b %>% summary()
results_pot_b %>% summary()
results_point_c <- list()
results_pot_c <- list()
for (i in 1:nrow(largeNetworks)) {
print(paste(i, "of", nrow(largeNetworks)))
x <- largeNetworks[i, ]
n1 <- matrix2zscore(x[["wordMatrixToRowListCampVisitingAdjacency"]])
diag(n1) <- 0
bc1 <- sna::closeness(n1)
n2 <- matrix2zscore(x[["wordMatrixToRowListCampPointAdjacency"]])
diag(n2) <- 0
bc2 <- sna::closeness(n2)
n3 <- matrix2zscore(x[["wordMatrixToRowListCampPotAdjacency"]])
diag(n3) <- 0
bc3 <- sna::closeness(n3)
r1 <- cor(bc2, bc1) %>% round(., 2)
r2 <- cor(bc3, bc1) %>% round(., 2)
results_point_c[[i]] <- r1
results_pot_c[[i]] <- r2
# ng = n1
# ng[ng > .75] = 1
# ng[ng <= .75] = 0
# g = ng %>% graph_from_adjacency_matrix(mode = "undirected", weighted = F)
# degree_distribution(g)
}
results_point_c <- unlist(results_point_c)
results_pot_c <- unlist(results_pot_c)
hist(results_point_c)
hist(results_pot_c)
results_point_c %>% summary()
results_pot_c %>% summary()
# examine bad network correlation
library(mgcv)
points_df <- runs %>%
filter(uniqueTraits == 1) %>%
select(expID, any_of(vars), any_of(contains("_corZscore"))) %>%
pivot_longer(any_of(contains("_corZscore"))) %>%
filter(name == "Visiting_Point_corZscore") %>%
pivot_longer(vars, names_to = "variable", values_to = "parameter") %>%
group_by(variable, parameter) %>%
summarize(cor = mean(value, na.rm = T)) %>%
mutate(parameter = as.character(parameter))
run_sample <- runs %>%
filter(uniqueTraits == 1) %>%
filter(pTrading > 0) %>%
filter(learningMethod == 50) %>%
filter(Visiting_Point_corZscore < 0) %>%
slice_sample(n = 1)
run_sample$experiment
run_sample$runNumber
run_sample$nTicks
default_vars <- defaults %>%
mutate(default = TRUE) %>%
select(parameters = Parameter, default_values = Default, default)
#### ERGM ####
library(statnet)
library(Rglpk)
library(ergm)
library(dplyr)
library(igraph)
# compare visiting and other networks to pot and point
# compare nonunique traits to unique traits
# statistics to include mcmc.diagnostics - check simulated distributions have a peak at the estimate and goodness of fit diagnostics
# are the same generative processes responsible for all networks? If we just had the material culture then would we still come to the same conclusions?
set.seed(1010)
run_ergm <- function(i, network, formula_type) {
formula <- switch(formula_type,
"edges" = ns ~ edges + nodematch("band"),
"triangle" = ns ~ edges + triangle + nodematch("band"),
"threetrail" = ns ~ edges + triangle + threetrail + nodematch("band"),
"gwesp" = ns ~ edges + gwesp(0.25, fixed = TRUE) + threetrail + nodematch("band"),
stop("Unknown formula type")
)
library(statnet)
library(igraph)
library(dplyr)
g <- readRDS("data/graphSample.Rds")
tryCatch(
{
cutoff <- quantile(g[[network]][[i]]$weight, .75)
gdf <- g[[network]][[i]] %>%
filter(weight > cutoff)
n <- graph_from_data_frame(gdf, directed = F)
nodeList <- gdf %>%
select(id = from, band) %>%
distinct_all()
ns <- intergraph::asNetwork(n)
ns %v% "band" <- nodeList$band
controlVals <- control.ergm(
parallel = 14,
seed = 1010,
MCMLE.confidence = 0.9,
MPLE.maxit = 1000
)
x <- tryCatch(
ergm(formula, control = controlVals),
error = function(e) {
print(e)
return(NA)
}
)
e <- tryCatch(
summary(x),
error = function(e) {
return(NA)
}
)
gf <- tryCatch(
gof(x),
error = function(e) {
return(NA)
}
)
f <- tryCatch(
x$formula %>% as.character() %>% paste(collapse = " "),
error = function(e) {
return(NA)
}
)
aic <- tryCatch(
e$aic,
error = function(e) {
return(NA)
}
)
gf_pval <- tryCatch(
gf$pval.model[1, 5],
error = function(e) {
return(NA)
}
)
edges <- tryCatch(
e$coefficients["edges", "Pr(>|z|)"],
error = function(e) {
return(NA)
}
)
result <- tibble(
run = i,
network = network,
formula = f,
aic = aic,
gf_pval = gf_pval,
edges = edges
)
saveRDS(
result,
paste0(
"data/ergm_results_",
i,
"_",
network,
"_",
formula_type,
".rds"
)
)
return(result)
},
error = function(e) {
return(NA)
}
)
}
run_and_monitor_processes <- function(bg, timeout = 120) {
start_times <- Sys.time()
# Monitor processes
while (TRUE) {
Sys.sleep(1) # Check every second
# Check the status of each process
for (name in names(bg)) {
if (bg[[name]]$is_alive()) {
run_time <- as.numeric(difftime(Sys.time(), start_times, units = "secs"))
if (shiny::isTruthy(run_time > timeout)) {
bg[[name]]$kill() # Kill process if it exceeds timeout
cat("Process", name, "killed after", run_time, "seconds\n")
}
}
}
# Exit the loop if all processes are done
if (all(!sapply(bg, function(p) p$is_alive()))) {
break
}
}
}
# run_ergm(1, "trading", "edges")
library(callr)
bg <<- list()
for (i in 1:10) {
print(i)
bg <- list()
for (network in c("trading", "learning", "visiting", "point", "pot")) {
print(network)
bg <- list()
for (formula_type in c("edges", "triangle", "threetrail", "gwesp")) {
l <- list(i = i, network = network, formula_type = formula_type)
bg[[paste0(i, "-", network, "-", formula_type)]] <- r_bg(run_ergm, args = l, supervise = FALSE)
}
run_and_monitor_processes(bg, timeout = 120)
}
}
ls <- list.files("data/", pattern = "ergm_results", full.names = T)
ls <- ls[which(!str_detect(ls, ".csv"))]
results <- lapply(ls, readRDS) %>% bind_rows()
results_summary <- results %>%
group_by(network, formula) %>%
summarize(mean_gf_pval = mean(gf_pval, na.rm = T), mean_edges = mean(edges, na.rm = T), mean_aic = mean(aic, na.rm = T), count = n()) %>%
ungroup() %>%
filter(count > 15)
results_table <- results_summary %>%
select(network, formula, gof = mean_gf_pval, aic = mean_aic) %>%
filter(network %in% c("visiting", "point")) %>%
mutate(formula = case_when(
str_detect(formula, "gwesp") ~ "edge + gwesp + threetrail + band",
str_detect(formula, "threetrail") ~ "edge + triangle + threetrail + band",
TRUE ~ "edge + band"
)) %>%
arrange(network, formula)
results_table
results_table %>% write_csv("tables/ergm_results_table.csv")
results_summary2 <- results %>%
group_by(formula) %>%
summarize(mean_gf_pval = mean(gf_pval), mean_edges = mean(edges), mean_aic = mean(aic))
results_summary2
# graphical representation of how an ergm works
library(igraph)
set.seed(1010)
g <- readRDS("data/graphSample.Rds")
i <- 10
network <- "point"
cutoff <- quantile(g[[network]][[i]]$weight, .75)
gdf <- g[[network]][[i]] %>%
filter(weight > cutoff)
n <- graph_from_data_frame(gdf, directed = F)
# nodeList = g[[network]][[i]] %>% select(id = from, band) %>%
# distinct_all()
nodeList <- structure(
list(id = c(
"5", "6", "7", "8", "9", "10", "11", "12",
"13", "14", "15", "16", "17", "18", "19"
), band = c(
"0", "0", "0",
"1", "1", "1", "2", "2", "2", "3", "3", "3", "4", "4", "4"
)),
row.names = c(
NA,
-15L
), class = c("tbl_df", "tbl", "data.frame")
)
g_df <- tibble(id = V(n)$name) %>%
left_join(nodeList)
V(n)$band <- g_df$band
band_colors <- c("0" = "red", "1" = "blue", "2" = "green", "3" = "orange", "4" = "purple")
V(n)$color <- band_colors[V(n)$band]
# png("figures/ergm_example.png", width = 800, height = 800, bg = "black")
# plot(n, main = "", vertex.label.color = "black")
# title(main = "Observed Network", col.main = "white")
# dev.off()
# g_random <- rewire(n, keeping_degseq(niter = vcount(n) * 10))
# png("figures/ergm_example_random.png", width = 800, height = 800, bg = "black")
# plot(g_random, main = "", vertex.label.color = "black")
# title(main = "Random Network", col.main = "white")
# dev.off()
# g <- make_ring(15) %>%
# add_edges(c(1,5, 2,6, 3,7, 4,8, 5,9, 6,10, 7,11, 8,12, 9,13, 10,14, 11,15))
# plot(g)
ns <- intergraph::asNetwork(n)
# plot(ns)
formula_type <- "gwesp"
formula <- switch(formula_type,
"edges" = ns ~ edges + nodematch("band"),
"triangle" = ns ~ edges + triangle + nodematch("band"),
"threetrail" = ns ~ edges + triangle + threetrail + nodematch("band"),
"gwesp" = ns ~ edges + gwesp(0.25, fixed = TRUE) + threetrail + nodematch("band"),
stop("Unknown formula type")
)
random_ergm <- ergm(formula)
summary(random_ergm)
pgof <- gof(random_ergm)
pdf(glue::glue("figures/gof_{network}_{formula_type}.pdf"))
plot(pgof)
dev.off()
# test probability that random networks are correlated
g <- readRDS("data/graphSample.Rds")
i <- sample(1:100, 1)
i
network <- "visiting"
cutoff <- quantile(g[[network]][[i]]$weight, .75)
gdf <- g[[network]][[i]] %>%
filter(weight > cutoff)
n <- graph_from_data_frame(gdf, directed = F)
nodeList <- gdf %>%
select(id = from, band) %>%
distinct_all()
ns <- intergraph::asNetwork(n)
ns %v% "band" <- nodeList$band
controlVals <- control.ergm(
parallel = 14,
seed = 1010,
MCMLE.confidence = 0.9,
MPLE.maxit = 1000
)
ergm_model <- ergm(ns ~ edges + triangle + threetrail + nodematch("band"), control = controlVals)
# ns_visiting = ns
nsim <- simulate(ergm_model, nsim = 10000)
gcor_values <- numeric(length = length(nsim))
# ns_matrix = as.matrix(ns)
ns_matrix <- as.matrix(ns_visiting)
# Loop through each simulated network and calculate gcor
for (i in seq_along(nsim)) {
nsim_matrix <- as.matrix(nsim[[i]]) # Convert simulated network to matrix
ns_row_names <- rownames(ns_matrix)
ns_col_names <- colnames(ns_matrix)
nsim_row_names <- rownames(nsim_matrix)
nsim_col_names <- colnames(nsim_matrix)
# Union of row and column names (all possible nodes in both matrices)
all_row_names <- union(ns_row_names, nsim_row_names)
all_col_names <- union(ns_col_names, nsim_col_names)
# Initialize matrices to the larger set of row and column names, filled with zeros
new_ns_matrix <- matrix(0,
nrow = length(all_row_names), ncol = length(all_col_names),
dimnames = list(all_row_names, all_col_names)
)
new_nsim_matrix <- matrix(0,
nrow = length(all_row_names), ncol = length(all_col_names),
dimnames = list(all_row_names, all_col_names)
)
# Fill the new matrices with the values from the original matrices
new_ns_matrix[ns_row_names, ns_col_names] <- ns_matrix
new_nsim_matrix[nsim_row_names, nsim_col_names] <- nsim_matrix
gcor_values[i] <- gcor(new_ns_matrix, new_nsim_matrix) # Calculate gcor
}
hist(gcor_values)
gcor_values %>% summary()
network <- "visiting"
cutoff <- quantile(g[[network]][[i]]$weight, .75)
gdf <- g[[network]][[i]] %>%
filter(weight > cutoff)
n <- graph_from_data_frame(gdf, directed = F)
nodeList <- gdf %>%
select(id = from, band) %>%
distinct_all()
ns <- intergraph::asNetwork(n)
ns %v% "band" <- nodeList$band
controlVals <- control.ergm(
parallel = 14,
seed = 1010,
MCMLE.confidence = 0.9,
MPLE.maxit = 1000
)
ergm_model <- ergm(ns ~ edges + triangle + threetrail + nodematch("band"), control = controlVals)
nsim <- simulate(ergm_model, nsim = 1000)
gcor_values <- numeric(length = length(nsim))
ns_matrix <- as.matrix(ns)
# Loop through each simulated network and calculate gcor
for (i in seq_along(nsim)) {
nsim_matrix <- as.matrix(nsim[[i]]) # Convert simulated network to matrix
gcor_values[i] <- gcor(ns_matrix, nsim_matrix) # Calculate gcor
}
# hist(gcor_values)
# gcor_values %>% summary()
#### export session info ####
# Get the R version
r_version <- paste0("R version=", paste(R.version$major, R.version$minor, sep = "."))
# Get the loaded packages and their versions
loaded_packages <- sessionInfo()$otherPkgs
loaded_versions <- sapply(loaded_packages, function(x) paste(x$Package, x$Version, sep = "="))
# Write the R version and package versions to the file
writeLines(c(r_version, loaded_versions), con = "analysis/r_environment.qmd")
Functions
# custom functions
library(tictoc)
library(magrittr)
library(ggnetwork)
library(ggthemes)
#' Load all data by specifying one complete name or the partial date file name
#'
#' @param fileName
#'
#' @return
#' @export
#'
#' @examples
loadRun <- function(fileName) {
if (file.exists(here(stringr::str_glue("results/{fileName}")))) {
fileName %<>% str_remove_all("run_[A-z]*_")
}
map(
list.files(here("results"),
pattern = fileName,
full.names = T
) %>% .[str_detect(., fileName)],
~
read_csv(.x, show_col_types = F)
) %>% set_names(c(
"assemblages",
"links",
"nodes",
"variables"
))
}
# x = runResults$assemblages$pointAssemblage[[7]]
#' Convert NetLogo list format into R list
bracket2list <- function(x) {
if (str_detect(x, "[0-9]")) {
x %>%
str_remove_all("\\]|\\[") %>%
str_split_fixed(" ", n = Inf) %>%
matrix(ncol = 2, byrow = T) %>%
as_tibble(.name_repair = "minimal") %>%
setNames(c("tick", "tool"))
} else {
NA
}
}
#' Calculate manhattan similarity
#'
#' @param analysis dataframe formatted through formatAssemblage function
#' @param type pot or point
#' @param group band or site
#' @param cutoff quantile cutoff to keep links
#' @param similarity either 'BR' for Brainerd-Robinson or 'WJ' for weighted Jaccard
#'
#' @return dataframe of nodes and dataframe of links
#' @export
#'
#' @examples
calcDistanceTools <- function(analysis,
toolType = "pot",
group = "band",
cutoff = 0,
similarity = "WJ") {
tmp <- analysis %>%
filter(type == toolType)
m <- table(tmp[[group]], tmp$tool) %>%
prop.table(1) %>%
{
. * 100
}
# calculate similarity
if (similarity == "BR") {
g <- dist(m, "manhattan") %>%
as.matrix()
g <- 1 - (g / 200)
diag(g) <- 0
# g[] = round(g, 2)
} else {
g <- weightedJaccard(m)
# g[] = round(g, 2)
}
links <- g %>%
as.data.frame() %>%
rownames_to_column("from") %>%
pivot_longer(-from,
names_to = "to",
values_to = "weight"
) %>%
# remove duplicate links
dplyr::rowwise() %>%
mutate(combo = paste(min(from, to), max(from, to))) %>%
arrange(desc(to)) %>%
distinct(combo, .keep_all = T) %>%
select(-combo)
links %<>%
filter(weight > 0)
cutoff <- quantile(links$weight, cutoff)
links %<>%
filter(weight >= cutoff) %>%
mutate_at(vars(weight), round, 2)
if (group != "band" & "band" %in% names(analysis)) {
nodes <- links %>%
select(from, to) %>%
pivot_longer(c(from, to), values_to = "name", names_to = "var") %>%
select(name) %>%
arrange(name) %>%
distinct_all() %>%
full_join(
analysis %>%
select(name = site, band) %>%
mutate_at(vars(name), as.character),
by = "name"
) %>%
mutate(camp = name) %>%
distinct_all() %>%
mutate_at(vars(band), factor)
} else {
nodes <- analysis %>%
ungroup() %>%
select(any_of(group)) %>%
mutate(name = !!as.name(group)) %>%
distinct_all()
}
return(list(links = links, nodes = nodes))
}
# function to remove labels from names
rmLabel <- function(x) {
map_chr(
x,
~ str_extract_all(.x, "[0-9]") %>%
unlist() %>%
paste(collapse = "")
)
}
# links = runResults$links
# nodes = runResults$nodes
calcDistanceLinks <- function(links,
nodes,
colS = "all",
group,
cutoff = 0) {
if (colS == "all") {
colS <- c("learning", "hunting", "visiting", "trading")
}
df <- links %>%
select(end1, end2, all_of(colS)) %>%
mutate_at(vars(end1, end2), rmLabel) %>%
left_join(
nodes %>%
select(end1 = who, band, camp) %>%
mutate(end1 = end1 %>% as.character()),
by = "end1"
) %>%
left_join(
nodes %>%
select(end2 = who, band, camp) %>%
mutate(end2 = end2 %>% as.character()),
by = "end2",
suffix = c("", "2")
)
linksMod <- df %>%
select(group, any_of(str_glue("{group}2")), all_of(colS)) %>%
# create unique paths
mutate(path = map2_chr(
!!as.name(group),
!!as.name(str_glue("{group}2")),
~ sort(c(.x, .y)) %>% paste(collapse = "_")
)) %>%
group_by(path) %>%
summarize_at(vars(all_of(colS)), sum) %>%
ungroup() %>%
separate(path, into = c(group, str_glue("{group}2")), sep = "_") %>%
mutate_at(vars(group, str_glue("{group}2")), as.character) %>%
distinct_all() %>%
mutate_at(vars(group, str_glue("{group}2")), factor) %>%
filter(!!as.name(group) != !!as.name(str_glue("{group}2"))) %>%
mutate(total = rowSums(select(., .dots = any_of(colS))))
newLinks <- map(colS %>% c("total"), function(c) {
x <- quantile(linksMod[[c]], cutoff)
linksTmp <- linksMod %>%
select(
from = all_of(group),
to = all_of(str_glue("{group}2")),
weight = !!as.name(c)
) %>%
filter(weight >= x) %>%
arrange(gtools::mixedsort(from))
}) %>% setNames(colS %>% c("total"))
newNodes <- df %>%
select(any_of(contains(c("camp", "band")))) %>%
distinct_all() %>%
rowid_to_column() %>%
pivot_longer(-rowid) %>%
mutate_at(vars(name), list(~ str_remove_all(.x, "2"))) %>%
distinct_all() %>%
pivot_wider(
names_from = name,
values_from = value,
values_fn = list
) %>%
unchop(c(band, camp)) %>%
select(-rowid) %>%
distinct_all() %>%
select(camp, band) %>%
mutate_all(factor) %>%
mutate(name = !!as.name(group))
return(list(links = newLinks, nodes = newNodes))
}
plotGraph <- function(links,
nodes,
title = "",
weighted = T,
color = NA,
shape = NA,
labels = T,
nodeSize = 10,
save = F,
layout = NA,
plotDevice = "svg",
groupByBand = F) {
# layout are coordinates for each node in a data frame format
# old <- .Random.seed
# on.exit({
# .Random.seed <<- old
# })
# set.seed(5010)
if (is.data.frame(layout)) {
if (is.matrix(links)) {
links <- links %>%
as.data.frame() %>%
rownames_to_column("from") %>%
pivot_longer(-from, names_to = "to", values_to = "weight")
}
plotdf <- links %>%
rename(name = from) %>%
full_join(layout, by = "name") %>%
left_join(
layout %>% rename(
to = name,
xend = x,
yend = y
),
by = "to"
)
} else if (is.data.frame(links)) {
n <- igraph::graph_from_data_frame(links, directed = F)
} else if (is.matrix(links)) {
colnames(links) <- rownames(links) <- nodes$name
n <- igraph::graph_from_adjacency_matrix(links, mode = "undirected", weighted = T, diag = F)
} else {
stop("check class of links")
}
if (!is.data.frame(layout)) {
if (groupByBand) {
V(n)$band <- nodes$band
nG <- contract(n, as.integer(V(n)$band) + 1, vertex.attr.comb = "first")
V(nG)$name <- nodes$band %>% unique()
nG <- simplify(nG, edge.attr.comb = "sum")
plotdf <- ggnetwork::ggnetwork(nG) %>%
as_tibble()
} else {
plotdf <- ggnetwork::ggnetwork(n) %>%
as_tibble()
}
}
plotdf %<>%
mutate(weight = weight %>% as.numeric()) %>%
left_join(nodes %>% mutate_at(vars(name), as.character),
by = "name"
)
gg <-
plotdf %>%
ggplot(aes(
x = x,
y = y,
xend = xend,
yend = yend
))
if (weighted == T) {
print("adding weights")
gg <- gg +
geom_edges(aes(color = weight), size = .75) +
scale_color_viridis_c(
direction = -1,
breaks = c(
min(plotdf$weight, na.rm = T),
max(plotdf$weight, na.rm = T)
),
labels = c("low", "high")
)
} else {
gg <- gg +
geom_edges(
color = "black",
size = .5,
alpha = .5
)
}
if (!is.na(shape) && shape == "band" & is.na(color)) {
print("adding shape")
gg <- gg +
geom_point(aes(shape = !!as.name(shape)),
size = nodeSize
) +
scale_shape_manual(values = c(15, 17, 19))
} else if (isTRUE(shape == "band" & !is.na(color))) {
print("adding shape and color")
gg <- gg +
geom_point(
aes(
shape = !!as.name(shape),
fill = !!as.name(color)
),
size = nodeSize
) +
scale_fill_brewer(palette = "Set1") +
scale_shape_manual(values = c(21, 22, 24))
} else if (shape != "band" & !is.na(shape)) {
stop("too many shapes -- maybe try color?")
} else if (is.na(shape) & !is.na(color)) {
print("adding color")
gg <- gg +
geom_point(aes(fill = !!as.name(color)),
size = nodeSize,
shape = 21
) +
scale_fill_brewer(palette = "Set1")
} else {
gg <- gg +
geom_point(
size = nodeSize,
color = "black",
shape = 19
)
}
if (labels == T) {
print("adding labels")
gg <- gg +
geom_nodetext(aes(label = name),
# size = 3,
fontface = "bold",
color = "white"
)
}
# determine if legend shape should be overridden
if (!is.na(shape)) {
f1 <- guide_legend()
} else if (isTRUE(!is.na(color) && shape != color) | !is.na(color)) {
f1 <- guide_legend(override.aes = list(shape = 21))
} else {
f1 <- guide_legend()
}
gg <- gg +
theme_blank() +
theme(
legend.position = "bottom",
legend.background = element_rect(fill = F),
plot.margin = margin(0, 0, 0, 0)
) +
coord_cartesian(clip = "off") +
guides(
fill = f1,
shape = guide_legend(direction = "horizontal"),
color = guide_colorbar(
direction = "horizontal",
title.position = "top"
)
) +
ggtitle(title)
print(gg)
if (save == T) {
if (plotDevice == tolower("svg")) {
ggsave(
filename = glue::glue("plots/{title}.svg"),
plot = gg
)
} else if (plotDevice == tolower("png")) {
ggsave(
filename = glue::glue("plots/{title}.png"),
plot = gg,
dpi = 450
)
} else {
stop("how about svg or png instead?")
}
}
return(gg)
}
#' Calculate weighted Jaccard
weightedJaccard <- function(m) {
#' weighted jaccard similarity matrix setup
#' from:
#' https://rpubs.com/lgadar/weighted-jaccard
r <- matrix(0, nrow = nrow(m), ncol = nrow(m))
rownames(r) <- rownames(m)
colnames(r) <- rownames(m)
# weighted jaccard
pairs <- t(combn(1:nrow(m), 2))
for (i in 1:nrow(pairs)) {
num <-
sum(sapply(1:ncol(m), function(x) {
(min(m[pairs[i, 1], x], m[pairs[i, 2], x]))
}))
den <-
sum(sapply(1:ncol(m), function(x) {
(max(m[pairs[i, 1], x], m[pairs[i, 2], x]))
}))
r[pairs[i, 1], pairs[i, 2]] <- num / den
r[pairs[i, 2], pairs[i, 1]] <- num / den
}
r[which(is.na(r))] <- 0
diag(r) <- 0
return(r)
}
#' Format netlogo adjacency tables as matrices
#'
#' @param x Netlogo matrix
#'
#' @return
#' @export
#'
#' @examples
string2MatrixR <- function(x) {
m <- x %>%
str_remove_all("\\[\\[|\\]\\]") %>%
str_split_1(pattern = fixed("] [")) %>%
as.matrix() %>%
apply(1, str_split_1, pattern = " ")
m[] <- as.double(m[])
return(m)
}
library(Rcpp)
sourceCpp("analysis/string2matrix.cpp")
#' Format Netlogo Adjacency Tables
#'
#' @param p run number
#' @param nB number of bands
#' @param cutoff cutoff value
#' @param x adjacency table formatted as string
#'
#' @return
#' @export
#'
#' @examples
formatNLAdjTables <- function(x, p = NULL, standardize = T) {
result <- purrr::map(x, function(x1) {
m <- string2matrix(x1)
nB <- ncol(m) / 3
nms <- nB:(nB * 3 + nB - 1)
bands <- data.frame(band = as.character(sort(rep(0:(nB - 1), 3))), camp = as.character(nms))
colnames(m) <- rownames(m) <- nms
df <- as.data.frame(m) %>%
rownames_to_column("from") %>%
pivot_longer(all_of(nms %>% as.character()),
names_to = "to",
values_to = "weight"
)
if (!is.null(p)) {
df <- df %>%
mutate(property = p) %>%
mutate(id = map2_chr(from, to, function(f, t) {
sort(c(f, t)) %>%
paste(collapse = "_")
})) %>%
distinct(id, .keep_all = T) %>%
select(-id) %>%
separate(
col = property,
into = c("run_number", "property"),
sep = "_",
extra = "drop"
)
}
df <- df %>%
rowwise() %>%
mutate(id = sort(c(from, to)) %>% paste(collapse = "_")) %>%
ungroup() %>%
distinct(id, .keep_all = T) %>%
filter(from != to) %>%
filter(weight > 0)
if (standardize) {
df <- df %>%
mutate(weight = weight / max(weight))
}
df <- df %>%
# mutate_at(vars(from,to),as.character) %>%
left_join(bands %>% rename(from = camp), by = "from")
return(df)
})
# nodes = tibble(camp = nms,
# band = rep((1:nB - 1), 3) %>% sort) %>%
# mutate_all(as.character) %>%
# mutate(name = camp)
# return(list(links = links, nodes = nodes))
return(result)
}
#' Format assemblages into a data frame structure
formatAssemblages <- function(assemblages) {
# tictoc::tic()
tmp <- assemblages %>%
# remove the aggregation sites
filter(!is.na(band)) %>%
mutate_at(vars(band, site), as.character)
pot <- tmp %>%
select(band, site, potAssemblage) %>%
mutate(potAssemblage = potAssemblage %>%
map(bracket2list)) %>%
unnest(potAssemblage) %>%
mutate(type = "pot") %>%
select(band, site, tick, type, tool)
point <- tmp %>%
select(band, site, pointAssemblage) %>%
mutate(pointAssemblage = pointAssemblage %>%
map(bracket2list)) %>%
unnest(pointAssemblage) %>%
mutate(type = "point") %>%
select(band, site, tick, type, tool)
# print(tictoc::toc()$elapsed)
return(bind_rows(pot, point))
}
#' Create attribute network
attrNetworks <- function(runResults,
attr = "pointAttributes",
group = "who") {
links <- runResults$nodes %>%
select(who, band, camp, !!as.name(attr)) %>%
mutate(attr = !!as.name(attr) %>% str_remove_all("\\[|\\]")) %>%
select(-!!as.name(attr)) %>%
separate_rows(attr, sep = " ")
# table(links$who,links$attr)
if (group == "camp") {
links %<>%
select(-who) %>%
distinct_all()
} else if (group == "band") {
links %<>%
select(-who, -camp) %>%
distinct_all()
}
d <-
table(
links[[group]],
links$attr
)
n <- nrow(d)
m <- matrix(nrow = n, ncol = n)
rownames(m) <- rownames(d)
colnames(m) <- rownames(d)
for (i in 1:n) {
for (j in 1:n) {
r1 <- which(d[i, ] > 0)
r2 <- which(d[j, ] > 0)
m[i, j] <- intersect(
r1,
r2
) %>%
length()
}
}
m
diag(m) <- 0
m[which(m < 2)] <- 0
g <- list()
g$links <- as.data.frame(m) %>%
rownames_to_column("from") %>%
pivot_longer(2:ncol(.), names_to = "to", values_to = "weight") %>%
filter(weight > 0) %>%
mutate(id = map2_chr(
from, to,
~ paste0(min(.x, .y), "_", max(.x, .y))
)) %>%
distinct(id, .keep_all = T) %>%
select(-id)
g$nodes <- links %>%
mutate(name = !!as.name(group)) %>%
mutate_all(as.character) %>%
mutate_at(vars(any_of(c(
"who", "name", "camp", "band"
))), list(function(x) {
factor(x, levels = unique(x) %>% gtools::mixedsort())
}))
return(g)
}
#' Select a proportion of the data either randomly or nonrandomly
sampleData <- function(df,
p, # proportion of each site to keep
method # c("random","nonrandom")
) {
s <- df %>%
group_by(site) %>%
count() %>%
mutate(
c = floor(n * p),
c = case_when(c == 0 ~ 1, TRUE ~ c)
)
if (method == "random") {
rs <- df %>%
group_by(site) %>%
group_split()
r <- map2_dfr(rs, s$c, ~ {
.x %>%
slice_sample(n = .y)
})
} else {
rs <- df %>%
slice(gtools::mixedorder(site)) %>%
arrange(site) %>%
rowid_to_column() %>%
group_by(site) %>%
group_split()
r <- map2_dfr(
rs,
s$c, function(rsi, si) {
print("asdf")
indx <- rsi %>%
group_by(site) %>%
slice_sample(n = 1) %>%
pull(rowid)
indx2 <- indx + si
if (indx2 > max(rsi$rowid)) {
tmp1 <- rsi %>%
filter(
rowid >= indx,
rowid <= max(rsi$rowid)
)
n <- si - nrow(tmp1)
indx <- rsi$rowid[which.min(rsi$rowid)]
indx2 <- indx + n
tmp2 <- rsi %>%
filter(
rowid >= indx,
rowid < indx2
)
rsi <- bind_rows(tmp1, tmp2)
} else {
rsi %<>%
filter(
rowid >= indx,
rowid < indx2
)
}
return(rsi)
}
)
}
return(r)
}
compareMissing <- function(df,
group,
prop, # proportion of each site to keep
method # c("random","nonrandom")){
) {
df %<>% rename(group = !!as.name(group))
original <- table(df$site, df$tool) %>%
prop.table(1) %>%
{
. * 100
} %>%
weightedJaccard()
results <- map_dfr(prop, function(p) {
print(paste("running", method, p))
reduced <- sampleData(df = df, p = p, method = method)
sampled <- table(reduced$site, reduced$tool) %>%
prop.table(1) %>%
{
. * 100
} %>%
weightedJaccard()
l <- list(original, sampled)
distance <- try(round(nd.dsd(l,
out.dist = T,
type = "Lap"
)$D, 3))
if (inherits(distance, "try-error")) {
distance <- NA_real_
}
result <- tibble(
count = nrow(reduced),
proportion = p,
method = method,
distance = unlist(distance)[1]
)
return(result)
})
return(results)
}
formatMatrix <- function(net, q = 0, binarize = F, rescale = T, groupByBand = F) {
if (inherits(net, "list")) {
net <- net[[1]]
}
if (inherits(net[1, 1], "character")) {
net <- apply(net, 1, as.numeric)
}
if (groupByBand) {
nB <- ncol(net) / 3
bNms <- 0:(nB - 1)
nms <- (nB):(ncol(net) + (nB - 1))
colnames(net) <- rownames(net) <- nms
n <- igraph::graph_from_adjacency_matrix(net, mode = "undirected", weighted = T, diag = F)
igraph::V(n)$band <- bNms
nG <- igraph::contract(n, as.integer(V(n)$band) + 1, vertex.attr.comb = "first")
igraph::V(nG)$name <- bNms
nG <- igraph::simplify(nG, edge.attr.comb = "sum")
net <- igraph::as_adjacency_matrix(nG, attr = "weight") %>% as.matrix()
}
if (q > 0) {
cutoff <- quantile(net, q)
net[which(net <= cutoff)] <- 0
if (binarize) {
net[which(net > cutoff)] <- 1
}
}
if (rescale) {
net <- log(net + 1)
net <- scales::rescale(net)
}
rownames(net) <- colnames(net) <- 1:ncol(net)
return(net)
}
adjacency_compare <- function(net1, net2, q = 0, binarize = F, rescale = T) {
net1 <- formatMatrix(net1, q = q, binarize = binarize, rescale = rescale)
net2 <- formatMatrix(net2, q = q, binarize = binarize, rescale = rescale)
net_diff <- abs(net1 - net2)
out <-
1 - (sum(net_diff) / (nrow(net_diff) * (nrow(net_diff) - 1)))
return(out)
}
adjacency_compare2 <- function(net1, net2) {
if (!identical(dim(net1), dim(net2))) {
warning("matrices must be the same size: returning NA")
return(NA)
} else if (any(is.na(net1)) || any(is.na(net2))) {
warning("matrices must not be NA: returning NA")
return(NA)
} else if (any(is.infinite(net1)) || any(is.infinite(net2))) {
warning("matrices must not be Infinite: returning NA")
return(NA)
} else {
diag(net1) <- diag(net2) <- 0
net_diff <- abs(net1 - net2)
out <- 1 - sum(net_diff) / (nrow(net_diff) * nrow(net_diff))
return(out)
}
}
adjacency_compare3 <- function(net1, net2) {
net_diff <- sum(abs(net1 - net2))
return(net_diff)
}
cor_compare <- function(net1, net2) {
tryCatch(
{
if (!identical(dim(net1), dim(net2))) {
warning("matrices must be the same size: returning NA")
return(NA)
} else if (any(is.na(net1)) || any(is.na(net2))) {
warning("matrices must not be NA: returning NA")
return(NA)
} else if (any(is.infinite(net1)) || any(is.infinite(net2))) {
warning("matrices must not be Infinite: returning NA")
return(NA)
} else {
diag(net1) <- diag(net2) <- 0
out <- sna::gcor(net1, net2)
return(out)
}
},
error = function(e) {
return(NA)
}
)
}
qapTestFun <- function(net1, net2) {
tryCatch(
{
if (!identical(dim(net1), dim(net2))) {
warning("matrices must be the same size: returning NA")
return(NA)
} else if (any(is.na(net1)) || any(is.na(net2))) {
warning("matrices must not be NA: returning NA")
return(NA)
} else if (any(is.infinite(net1)) || any(is.infinite(net2))) {
warning("matrices must not be Infinite: returning NA")
return(NA)
} else {
diag(net1) <- diag(net2) <- 0
out <- sna::qaptest(list(net1, net2), sna::gcor, g1 = 1, g2 = 2)$pleeq
return(out)
}
},
error = function(e) {
return(NA)
}
)
}
adjacency_correlation <- function(net1, net2, q = 0, binarize = T, rescale = F) {
net1 <- formatMatrix(net1, q = q, binarize = binarize, rescale = rescale)
net2 <- formatMatrix(net2, q = q, binarize = binarize, rescale = rescale)
library(igraph)
g1 <- graph_from_adjacency_matrix(net1, mode = "undirected", diag = F)
g2 <- graph_from_adjacency_matrix(net2, mode = "undirected", diag = F)
library(multinet)
net <- ml_empty("comparison")
add_igraph_layer_ml(net, g1, "net1")
add_igraph_layer_ml(net, g2, "net2")
m <- layer_comparison_ml(net, method = "pearson.degree", mode = "all")
diag(m) <- NA
return(max(m, na.rm = T))
}
adjacency_laplacian <- function(net1, net2, q = 0, binarize = F, rescale = T) {
net1 <- formatMatrix(net1, q = q, binarize = binarize, rescale = rescale)
net2 <- formatMatrix(net2, q = q, binarize = binarize, rescale = rescale)
library(NetworkDistance)
l <- list(net1, net2)
m <- round(nd.dsd(l,
out.dist = T,
type = "Lap"
)$D, 3)
return(m)
}
adjacency_laplacian2 <- function(net1, net2) {
if (!identical(dim(net1), dim(net2))) {
warning("matrices must be the same size: returning NA")
return(NA)
} else if (any(is.na(net1)) || any(is.na(net2))) {
warning("matrices must not be NA: returning NA")
return(NA)
} else if (any(is.infinite(net1)) || any(is.infinite(net2))) {
warning("matrices must not be Infinite: returning NA")
return(NA)
} else {
library(NetworkDistance)
l <- list(net1, net2)
m <- round(nd.dsd(l,
out.dist = T,
type = "Lap"
)$D, 3)
return(m)
}
}
compareAdjacencyResults <- function(runs) {
results <- runs %>%
mutate_all(as.character) %>%
pivot_longer(-run_number,
names_to = "property",
values_to = "value"
) %>%
mutate(property = str_remove_all(
property,
"matrix_to_row_list_camp_|_adjacency"
)) %>%
pivot_wider(names_from = property, values_from = value) %>%
pivot_longer(visiting:pot, names_to = "property", values_to = "adj") %>%
unite("id", run_number, property, remove = F)
nB <- runs %>%
select(n_of_bands) %>%
dplyr::pull() %>%
unique()
formatted <- map2_df(
results$adj, results$id,
formatNLAdjTables, nB
)
nms <- nB:(nB * 3 + nB - 1)
nodes <- tibble(
camp = nms,
band = rep((1:nB - 1), 3) %>% sort()
) %>%
mutate_all(as.character) %>%
mutate(name = camp)
groups <- formatted$property %>%
unique() %>%
sort()
comparison <- tibble()
for (i in formatted$run_number %>%
unique() %>%
as.integer()) {
print(i)
try({
run <- formatted %>%
filter(run_number == i)
nms <- run$property %>%
unique() %>%
sort()
run %<>%
group_by(property) %>%
group_split() %>%
setNames(nms)
# plotGraph(run$hunting,nodes, "hunting",shape = 'band', color = 'band')
cutoff <- quantile(run$learning$weight, 0.5)
run$learning %<>% filter(weight > cutoff)
# plotGraph(run$learning,nodes, "learning",shape = 'band', color = 'band')
cutoff <- quantile(run$trading$weight, 0.5)
run$trading %<>% filter(weight > cutoff)
# plotGraph(run$trading,nodes, "trading",shape = 'band', color = 'band')
cutoff <- quantile(run$visiting$weight, 0.5)
run$visiting %<>% filter(weight > cutoff)
# plotGraph(run$visiting,nodes, "visiting",shape = 'band', color = 'band')
cutoff <- quantile(total$weight, 0.5)
# plotGraph(total,nodes, "total",shape = 'band', color = 'band')
cutoff <- quantile(run$point$weight, 0.5)
run$point %<>% filter(weight > cutoff)
# plotGraph(run$point,nodes, "point",shape = 'band', color = 'band')
cutoff <- quantile(run$pot$weight, 0.5)
run$pot %<>% filter(weight > cutoff)
# plotGraph(run$pot,nodes, "pot",shape = 'band', color = 'band')
# multilayer network
library(igraph)
library(multinet)
net <- ml_empty("comparison")
ig1 <- graph_from_data_frame(d = run$pot, directed = F, vertices = nodes)
E(ig1)$weight <- run$pot$weight
add_igraph_layer_ml(n = net, g = ig1, name = "pots")
ig2 <- graph_from_data_frame(d = run$point, directed = F, vertices = nodes)
E(ig2)$weight <- run$point$weight
add_igraph_layer_ml(n = net, g = ig2, name = "points")
ig3 <- graph_from_data_frame(d = run$hunting, directed = F, vertices = nodes)
E(ig3)$weight <- run$hunting$weight
add_igraph_layer_ml(n = net, g = ig3, name = "hunting")
ig4 <- graph_from_data_frame(d = run$learning, directed = F, vertices = nodes)
E(ig4)$weight <- run$learning$weight
add_igraph_layer_ml(n = net, g = ig4, name = "learning")
ig5 <- graph_from_data_frame(d = run$trading, directed = F, vertices = nodes)
E(ig5)$weight <- run$trading$weight
add_igraph_layer_ml(n = net, g = ig5, name = "trading")
ig6 <- graph_from_data_frame(d = run$visiting, directed = F, vertices = nodes)
E(ig6)$weight <- run$visiting$weight
add_igraph_layer_ml(n = net, g = ig6, name = "visiting")
m2raw <- layer_comparison_ml(net, method = "pearson.degree", mode = "all")
m2 <- round(m2raw, 2) %>%
as.data.frame() %>%
rownames_to_column("from") %>%
pivot_longer(-from, names_to = "to") %>%
filter(from != to) %>%
rowwise() %>%
mutate(variable = paste0(
paste(sort(c(from, to)), collapse = "->"),
"-pearson.degree"
), .before = 1) %>%
select(-from, -to) %>%
distinct_all()
# m %>% knitr::kable()
ig1m <- igraph::as_adjacency_matrix(ig1)
ig2m <- igraph::as_adjacency_matrix(ig2)
ig3m <- igraph::as_adjacency_matrix(ig3)
ig4m <- igraph::as_adjacency_matrix(ig4)
ig5m <- igraph::as_adjacency_matrix(ig5)
ig6m <- igraph::as_adjacency_matrix(ig6)
result <- bind_rows(
tibble(
variable = "pots->points-AdjComparison",
value = adjacency_compare(ig1m, ig2m)
),
tibble(
variable = "pots->hunting-AdjComparison",
value = adjacency_compare(ig1m, ig3m)
),
tibble(
variable = "pots->learning-AdjComparison",
value = adjacency_compare(ig1m, ig4m)
),
tibble(
variable = "pots->trading-AdjComparison",
value = adjacency_compare(ig1m, ig5m)
),
tibble(
variable = "pots->visiting-AdjComparison",
value = adjacency_compare(ig1m, ig6m)
),
tibble(
variable = "points->hunting-AdjComparison",
value = adjacency_compare(ig2m, ig3m)
),
tibble(
variable = "points->learning-AdjComparison",
value = adjacency_compare(ig2m, ig4m)
),
tibble(
variable = "points->trading-AdjComparison",
value = adjacency_compare(ig2m, ig5m)
),
tibble(
variable = "points->visiting-AdjComparison",
value = adjacency_compare(ig2m, ig6m)
),
tibble(
variable = "hunting->learning-AdjComparison",
value = adjacency_compare(ig3m, ig4m)
),
tibble(
variable = "hunting->trading-AdjComparison",
value = adjacency_compare(ig3m, ig5m)
),
tibble(
variable = "hunting->visiting-AdjComparison",
value = adjacency_compare(ig3m, ig6m)
),
tibble(
variable = "learning->trading-AdjComparison",
value = adjacency_compare(ig4m, ig5m)
),
tibble(
variable = "learning->visiting-AdjComparison",
value = adjacency_compare(ig4m, ig6m)
),
tibble(
variable = "trading->visiting-AdjComparison",
value = adjacency_compare(ig5m, ig6m)
),
m2
)
vars <- results %>%
filter(run_number == i) %>%
select(-id, -property, -adj) %>%
distinct_all() %>%
pivot_longer(-run_number, names_to = "variable")
result <- vars %>%
bind_rows(result %>%
mutate_at(vars(value), round, digits = 2) %>%
mutate(run_number = i) %>%
mutate_all(as.character))
comparison <- bind_rows(comparison, result)
})
}
return(comparison)
}
matrix2rank <- function(adj) {
m <- string2matrix(adj)
v <- sort(unique(as.vector(m)), decreasing = T)
vRank <- 1:length(v)
m[] <- match(m, v)
diag(m) <- 0
return(m)
}
# net1 = matrix2rank(runs$wordMatrixToRowListCampLearningAdjacency[2])
# net2 = matrix2rank(runs$wordMatrixToRowListCampPotAdjacency[2])
# adjacency_compare2(net1,net2)
matrix2rank2 <- function(adj) {
m <- adj
v <- sort(unique(as.vector(m)), decreasing = T)
vRank <- 1:length(v)
m[] <- match(m, v)
diag(m) <- 0
return(m)
}
matrix2zscore <- function(adj) {
m <- string2matrix(adj)
v <- (m - mean(m)) / sd(m)
m[] <- v
diag(m) <- 0
m[] <- round(m, 3)
m <- scales::rescale(m)
return(m)
}
matrix2zscoreBand <- function(adj) {
m <- string2matrix(x)
nc <- ncol(m) / 3
bands <- rep(0:(nc - 1), 3) %>% sort()
camps <- 1:length(bands)
dict <- tibble(bands, camps)
m <- as.data.frame(m) %>%
rowid_to_column("from") %>%
pivot_longer(-from, names_to = "to") %>%
mutate(to = str_remove_all(to, "V") %>% as.integer()) %>%
left_join(dict, by = c("from" = "camps")) %>%
rename(bands1 = bands) %>%
left_join(dict, by = c("to" = "camps")) %>%
rename(bands2 = bands) %>%
group_by(bands1, bands2) %>%
summarize(value = sum(value)) %>%
ungroup() %>%
pivot_wider(names_from = bands2, values_from = value) %>%
column_to_rownames("bands1") %>%
as.matrix()
diag(m) <- 0
v <- (m - mean(m)) / sd(m)
m[] <- v
diag(m) <- 0
m[] <- round(m, 3)
m <- scales::rescale(m)
return(m)
}
matrix2one <- function(adj) {
m <- string2matrix(adj)
m[] <- m / max(m)
diag(m) <- 0
m[] <- round(m, 3)
return(m)
}
# net1 = matrix2rank(runs$wordMatrixToRowListCampLearningAdjacency[2])
# net2 = matrix2rank(runs$wordMatrixToRowListCampPotAdjacency[2])
# adjacency_compare2(net1,net2)
formatNLAdjTablesOld <- function(x, p, nB = 3) {
nms <- nB:(nB * 3 + nB - 1)
links <- x %>%
strsplit(" |\\[|\\]") %>%
unlist() %>%
as.numeric() %>%
.[which(!is.na(.))] %>%
matrix(
nrow = length(nms),
ncol = length(nms),
byrow = T
) %>%
set_colnames(nms) %>%
set_rownames(nms) %>%
as.data.frame() %>%
rownames_to_column("from") %>%
pivot_longer(all_of(nms %>% as.character()),
names_to = "to",
values_to = "weight"
) %>%
mutate(
property = p,
id = map2_chr(from, to, function(f, t) {
sort(c(f, t)) %>%
paste(collapse = "_")
})
) %>%
distinct(id, .keep_all = T) %>%
select(-id) %>%
separate(
col = property,
into = c("run_number", "property"),
sep = "_",
extra = "drop"
) %>%
filter(from != to) %>%
filter(weight > 0)
return(links)
}
#' Format assemblages into a data frame structure
eigenDiff <- function(adj1, adj2) {
m1 <- matrix2rank(adj1)
m2 <- matrix2rank(adj2)
e1 <- sna::evcent(m1)
e2 <- sna::evcent(m2)
return(sum(abs(e1 - e2)))
}
# tmp = runs %>% slice(1:100)
# subset common vertices from 2 networks
network_subset_common <-
function(net1, net2, namevar = "vertex.names") {
net1_names <-
network::get.vertex.attribute(net1, attrname = namevar)
net2_names <-
network::get.vertex.attribute(net2, attrname = namevar)
common_names <- net1_names[which(net1_names %in% net2_names)]
out1 <- net1 %s% which(net1 %v% namevar %in% common_names)
out2 <- net2 %s% which(net2 %v% namevar %in% common_names)
return(list(out1, out2))
}R Environment
R version=4.4.3 callr=3.7.6 ergm=4.8.1 network=1.19.0 mgcv=1.9-1 nlme=3.1-167 Rcpp=1.0.14 ggthemes=5.1.0 magrittr=2.0.3 BayesFactor=0.9.12-4.7 Matrix=1.7-2 coda=0.19-4.1 magick=2.8.5 tictoc=1.2.1 ggnetwork=0.5.13 igraph=2.1.4 ggpubr=0.6.0 arrow=19.0.1 lubridate=1.9.4 forcats=1.0.0 stringr=1.5.1 purrr=1.0.4 readr=2.1.5 tidyr=1.3.1 tibble=3.2.1 ggplot2=3.5.1 tidyverse=2.0.0 here=1.0.1 flextable=0.9.7 gt=0.11.1 dplyr=1.1.4 rio=1.2.3
ODD
“ArchMatNet: Archaeological Material Networks”
ODD: Overview, Design concepts, and Details
Cecilia Padilla-Iglesias, Claudine Gravel-Miguel, and Robert Bischoff1
23 December 2022
Purpose
The purpose of the model is to investigate how different factors affect the ability of researchers to reconstruct prehistoric social networks from artifact stylistic similarities, as well as the overall diversity of cultural traits observed in archaeological assemblages. Given that cultural transmission and evolution are affected by multiple interacting phenomena, our model allows to simultaneously explore six sets of factors that may condition how social networks relate to shared culture between individuals and groups:
- Factors relating to the structure of social groups
- Factors relating to the cultural traits in question
- Factors relating to individual learning strategies
- Factors relating to the environment
- Factors relating to the context in which different types of cultural traits are learned
- Factors relating to the method used to reconstruct ancient social networks
Process overview and scheduling
At the beginning of each time-step, the model checks what type of environment is set. If it is only one large biome, it rolls the die to see if the environment will change. If it does, the new environment is a random value between 0 and 99. Patches get that new environment.
If the user wants highResolutionData, the model checks if it is time to export data (every 1000 time-steps). If it is the case, it calls procedures that compile data and the procedure that exports it to the correct CSV file.
People remove any superfluous alliances they may have accumulated in previous tick, using the cleanUpAlliances procedure.
People run the possibility to invent new types and new variants and create one pot and point each. A person may create a new pot and/or point each tick. Each new object is created based on a randomly selected trait from the person’s trait list. Each object also records the tick when it was created for use in subsequent analyses.
The model then checks if it is time to aggregate (based on the aggregation-freq slider). If it is, it sets aggregation to T, and calls the aggregate procedure, which moves people to their respective aggregation camps.
While aggregating, the model keeps track of the aggregation’s length. If it has reached the aggregationDuration limit, people go back to their own camp. When they are aggregating, people “visit” with all the people they are aggregated with, which is recorded in the interactions’ “visiting” counter, as well as the inter-camps campInteractions “visiting” counters. Aggregated people can go hunting with a certain probability.
When an aggregation is not occurring, people can do different activities, with a given probability. They can create a new alliance with someone from another camp. They can go visit an ally. They can go hunting.
When hunting, people gather other people from their band and go to a different location together. When the hunt occurs during an aggregation, the hunting site is within 15 patches of the aggregation site. When the hunt occurs at other times, the hunting site is simply within the band territory. On a hunt, the hunters calculate their pointPrestige, and have different probabilities to drop a point, trade points with one another, and transmit trait knowledge to one another. Hunters keep track of how many times they hunted together. This interaction is also recorded by the campInteractions counter between camps. When the hunt is over, everyone returns to their own camp.
When visiting, one person chooses one ally to go visit at random. They move to that ally’s camp and visit with all the allies present there. The people who are involved in the visit keep track of how many times they visited one another, which is also recorded by the campInteractions visiting counter. The people who participate in the visit have different probabilities to drop one pot, drop one point, trade objects with one another, and transmit trait knowledge to one another. When those activities are done, the visitor returns to their own camp.
If recordPrestige is T, camps call the recordCampPrestige procedure. This procedure helps to follow the progression of mean and max prestige of the people of each camp, to see how those values change over time.
All people have a certain probability to call the loseTrait procedure, in which they forget one pot or point trait.
All people call the transmitCulture procedure to share knowledge of traits, and tradeItems to trade objects from their inventories.
Then, all people call the dropItem procedure, in which they have a certain probability to drop one pot or one point.
After the first 1000 ticks, people check if they want to migrate or not at every tick, based on a probability set by the user. When genderedMigration is T, only people with gender 1 roll the die and they have twice the probability to migrate as everyone does when genderedMigration is F. When migration occurs, the person who is migrating chooses another person to switch place with. When genderedMigration is T, only people of gender 1 can be selected. The two people switch home camps.
The model checks to see if it has reached the maximum number of ticks. If it has and the user wants to see the metrics, it computes and prints all the metrics to the Command Center, exports results as a CSV file if saveResults is T, and stops.
Design concepts
Basic principles
We built this model to allow testing assumptions that archaeologists have made about the possibility to reconstruct social networks from material culture, as well as improve our methods to reconstruct said networks.
Emergence
We expect styles to emerge over time within well-connected people who share the same traits and innovate on the shared known traits.
Adaptation
The people do not adapt to their environment or their context. Their behavior is constant throughout the run.
Objectives
One objective people have is to learn fitter traits (when learning via Prestige).
Learning
People regularly learn from the people around them. This is the base of the cultural transmission subprocedure. While the learning is not considered adaptation, people will learn the traits that are most popular at any given time, traits which change as innovations gain popularity or are lost.
Prediction
There is no prediction.
Sensing
People sense the environment onto which they hunt, and it affects their prestige. When meeting one another, people can sense which traits the others know. People also know who their allies are, and how often they have visited them. Camps and bands sense the knowledge of their inhabitants. None of the knowledge involve errors. Interaction Interactions make the most important part of the model. People hunt together and visit one another. When interacting, people can trade items or share knowledge of some traits.
Stochasticity
To keep the model general, we use stochasticity often. A lot of the activities are done with a certain probability, controlled by Interface sliders. Traded items are chosen at random amongst traders’ inventories. When creating a new variant, the type to innovate upon is chosen randomly amongst the traits known by the innovator. Forgotten traits are also selected at random. This use of stochasticity is also because we do not know how often and how those activities were done in the past. Therefore, it is simpler to use probabilities to get a wider range of behaviors.
Collectives
People are grouped by camps, which are grouped by bands. There are three camps per band. People aggregate and hunt with the people from their band, therefore, their group affect who they will interact with the most, which in turn affects the similarity of their known traits and created objects.
Observation
The model can produce different types of outputs. It can calculate and print the following metrics (in the following names, X can be point or pot): 1) AssemblageSizeX: the sum of pots or points discarded at camps and bands; 2) interactionFreqY (where Y can be hunting, visiting, learning, and trading): the sum of all Y interactions between people; 3) brDistBandX: mean Brainerd Robinson distance between known traits, grouped by bands; 4) brDistCampX: mean Brainerd Robinson distance between known traits, grouped by camps; 5) meanPeoplePrestige: mean hunting prestige of all people; 6) modBandX: modularity metric which evaluates if there is more between-band links than within-band links, based on the Jaccard similarity of traits known by people; 7) modCampX: modularity metric which evaluates if there is more between-camp links than within-camp links, based on the Jaccard similarity of traits known by people; 8) networkMetricX: difference between the clustering coefficient of the network created by people who interacted with one another (hunting and learning) and the clustering coefficient of the network created by similarities in objects (point and pot); 9) PhiSTXBand: PhiST measure to see if variation in known traits is structured around bands; 10) PhiSTXCamp: PhiST measure to see if variation in known traits is structured around camps; 11) nTraitX: sum of the known traits at the end of the simulation (includes duplicates); and 12) elapsedTime: time elapsed since the start of the run.
In addition, the model can produce CSV outputs recording the parameter values, all the interactions between people, as well as the characteristics of all objects discarded at camp/band sites.
Input data
This model does not import external data.
Submodels
Social activities
People meet one another through hunting and visiting. Hunting is done with bandmates only, but visiting is done between allies, which can be of other bands, as hunter-gatherers in the ethnographic literature tend to travel further distances to visit acquaintances than to carry out subsistence activities (Hewlett et al. 1982; Weissner, 1977; Hill et al. 2014)
When someone initiates a hunt, they select a certain number of bandmates to bring with them. They all move to the hunting spot, which is either within the band territory, or within 15 patches of the aggregation site (when the hunt happens during an aggregation). All people taking part in one hunt update their and their camps’ interaction “hunting” counter. Hunting people calculate the fitness of their known traits in the current environment, by calculating how many of their point traits were invented in the current environment. They then have a certain probability of dropping a point, trading items and/or transmit trait knowledge with a fellow hunter. After that, all hunters return to their own camp.
To visit an ally, a person has first to have allies. Allies are formed randomly and with a given probability at each tick. When forming an alliance, a person chooses another person at random. That ally can be of any camp other than ego’s. When someone has more than the limit of allies allowed, they remove one. To select which alliance will go, the person uses weighted sampling, with the number of visits serving as the inverted weight (alliances of people who visit one another often will have less chances of being removed).
When visiting an ally, a person moves to the camp of one of their allies, choosing at random. If they have other allies in that camp, they all participate in the visit and they all update their and their camps’ interaction “visiting” counter. They then have a certain probability of dropping a point and a pot, trading items and/or transmit trait knowledge with a fellow ally. After that, the visitor returns to their own camp. Item exchange
The tradeItems procedure controls exchange between people. It is only run if the person has a point or pot to give and has found a trading partner on the same patch who also has an object to trade in its inventory. This is done first for points, and then for pots. The trading partners for both objects do not need to be the same, but they can be. Both partners choose one item at random from their inventory, and swap them. The interaction link between the agents and the campInteraction links between their camps record that interaction.
Cultural transmission
Cultural transmission happens when the transmitCulture is called (in GO, hunt, and visit). This procedure identifies if the transmission will take place (based on probabilities) and the type of transmission. For pots, and then points, people roll the die to see if they will transmit information. When successful, they will call transmitPrestige with a probability set by the learningMethod slider. Every other time, they will call transmitConformism. These procedures return a set of new traits, following the rules detailed below.
After learning the new traits, people craft a new item (calling the craftItem procedure) and remove any duplicate traits they might have just learned. The transmitConformism [ toolType ] reporter takes as toolType either “pot” or “point” and is used to determine which attributes are copied.
- Ego asks all people on the same patch (including itself) to select their traits of the given toolType.
- If visibility changes how often traits are transmitted, everyone uses the getMostlyVisibleTraits reporter to get mostly visible traits from their list. High visibility traits are chosen for transmission 98% of the time. This value is arbitrary but generally fits expectations from a broad reading of the ethnographic literature.
- Everyone pools their list together and the reporter reports the most common trait. Then, the interaction links between all people update their learning counter (+1).
The transmitPrestige [ toolType ] takes as toolType either “pot” or “point” and is used to determine which attributes are copied.
- Ego identifies who on the same patch they will learn from. For points, they select the person with the highest pointPrestige. If potPrestige is T, they choose the person with the highest potPrestige, but if potPrestige is F, they choose someone at random. This person is now the “teacher.”
- The teacher adds its list of traits to a temporary global list.
- If visibility affects transmission, the list is filtered through the getMostlyVisibleTraits reporter.
- The reporter chooses and reports one of the traits in the teacher’s list. The interaction links between the teacher and Ego updates its learning counter (+1).
The getMostlyVisibleTraits [ traitList ] reporter selects and reports mostly visible tools from the traitList given. There is some error as the visibility of a tool becomes its probability of being chosen.
- This first separates the traitList into two lists. One that focuses on the main types, and one that focuses on types with variants. For the main types, it simply pairs this new list with a list of 2s of the same length, as all main tools have high visibility. For variants, it pairs the list with each variants’ visibility obtained from the visibilityTable.
- It then merges both paired lists and samples traits from that list, using their score as the probability to be chosen. The list of chosen traits is as long as the original list and includes replacements.
- The reporter then removes any duplicates and reports this list.
Cultural innovation
Innovation can occur through the procedures luckyLeap (where a person creates a new main type) and createVariant (where a person creates a variant of a main type they know) (see Kolodny et al. 2015; Creanza et al. 2017 for a similar implementation of cultural dynamics in an ABM. This allows creating innovations that are dependent on other innovations, which is a key feature of cultural innovations. For example, the invention of a new snare hunting method may enable new kinds of game hunting and in turn the further invention of specialized tools for processing these novel foods. Such dependence in culture is what makes it cumulative over time.These procedures add the newly invented trait to a person’s repertoire and to the master pot and point traits tables. Main types are coded as integers, and their variants are coded as letters.
The luckyLeap [ nextTrait toolType ] procedure is called by a person (ego) in Go if a random value is below the threshold set in pLuckyLeap (Interface). The value is multiplied by 2 for pots so that pots are twice less likely to see random innovations.
The nextTrait is the number of the next trait in the list of main tools to invent, calling getNextNumber to calculate it. The getNextNumber reporter takes the list of pot or point attributes ego has and reports the number following the number ego already knows. That number is converted as a string to be used by luckyLeap. The toolType can be either “pot” or “point” and is used to determine which traits are innovated.
- Ego adds the nextTrait to its relevant trait list (point or pot).
- It then removes any duplicates.
- If this type is not in its relevant master table, it adds an entry for it.
- If toolType is “point” and the traitEnvironment table does not have it already, ego adds an entry to that table with the key = trait, and the value = the environmentZone of the patch where ego stands.
The createVariant [ toolType ] procedure is called by a person in Go if a roll die is successful. It takes one main type from the agent’s known traits and adds a new variant to it.
- It first populates a temporary variable with the list of traits the agent knows for the given toolType.
- It then chooses one trait at random from that list. It uses the getTraitNumber reporter to get the number of the trait (main type). It also identifies the last variant the agent owns of that type, using the getLastVariant reporter.
- The procedure then checks if ego has reached the limit of possible variants (“o”). If they have not, they identify the next variant in the alphabet. This letter becomes the newVariation.
- Ego merges the newVariation to the main type to create a new type with variant.
- If toolType is “point” and the traitEnvironment table does not have it already, ego adds an entry to that table with the key = trait, and the value = the environmentZone of the patch where ego stands.
Creating the archaeological record
The dropItem procedure removes objects from a person’s inventory based on the pointsBreakageRate or potsBreakageRate parameters, which are now set to 0.3. If the person dropping an object is hunting, the probability of discard is unchanged. If not hunting, the probability is divided by twenty to represent less likelihood of losing a point while engaged in other activities. This value was chosen to represent the probabilities of a point being discarded at a camp based on Claudine Gravel-Miguel et al. 2021. The potsBreakageRate does not vary. When dropping an item, a person removes the first item in its inventory (the oldest item made). If the person is on a camp or aggregation site, that site records the discarded object in its pot or point assemblage.
Global Extensions Description {-} These are the extensions used in the model. Note that stats is not bundled and must be downloaded and placed in the extensions folder of the Netlogo app or in the same directory as the model. sr is commented out as it requires additional setup, but can be enabled by uncommenting.
extensions [
profiler
rnd
nw
table
matrix
stats
; sr
]
Global Variables Description {-} These are the global variables defined for the model.
globals [
inventorySizePoints ; size of point inventory allowed per person
inventorySizePots ; size of pot inventory allowed per person
groupColors ; list of colors possible for each group
currentView ; determines whether environmental zones or territories are displayed
aggregation ; Monitor whether people are currently aggregating or not
aggregationCounter ; A counter to monitor the aggregation period
hideLinksInteractions ; whether links should start out hidden or not
randomWho ; random person to record activity for
traitLetters ; for trait variations (new variants on main types)
pointTraitsTable ; table that records all point traits created
potTraitsTable ; table that records all pot traits created
traitEnvironments ; table recording environment in which each point trait got created (the trait, not the object itself)
visibilityTable ; table that records which letter traits are visible (2) or not (1)
potBreakageRate; rate that pots break -- do not break while hunting
pointsBreakageRate; rate that points break -- more likely to break while hunting
campVisitingAdjacency ; record visits between their agent and other people
campTradingAdjacency ; record trades between their agent and other people
campHuntingAdjacency ; record hunts between their agent and other people
campLearningAdjacency ; record learning between their agent and other people
campPointAdjacency ; record the number of similar points between camps
campPotAdjacency ; record the number of similar pots between camps
fileHighRes ; file to record the high resolution dataset
]
Agent Types Description {-} Each agent type forms a nested hierarchy–bands contain camps which contain people.
breed [ bands band ] ; bands represent a higher-level grouping of people
breed [ camps camp ] ; camps represent a small group of people
breed [ people person ] ; people belong to a camp and its band
undirected-link-breed [ interactions interaction ] ; keep track of interactions between people
undirected-link-breed [ campInteractions campInteraction ] ; to record interactions between people, compiled by camps
undirected-link-breed [ tempLinks tempLink ] ; links created temporarily between individuals, camps, or bands to calculate different metrics
Agent Properties Description {-} These are the variables specific to each type of agent.
camps-own [
potAssemblage ; pots deposited at camp
pointAssemblage ; points deposited at camp
ownedBy ; which band owns the camp
otherCamps ; record camps that are not this one
otherCampsMyband ; records camps that are in this camp's band
meanPrestige ; record mean prestige per camp
maxPrestige ; record prestige of most prestigious individual in a camp
allPoints ; point traits known by the camp's people (with duplicates)
allPots ; point traits known by the camp's people (with duplicates)
]
bands-own [
potAssemblage ; pots deposited at band
pointAssemblage ; points deposited at band
allPoints ; point traits known by the band's people (with duplicates)
allPots ; pot traits known by the band's people (with duplicates)
]
people-own [
activity ; variable to hold the activity the agent is engaged in
originalHome ; location of originalHome camp
homeBand ; origin band for each person
homeCamp ; camp the person lives at
pointTraits ; point traits known (can be used to create a point or be transfered to someone else)
potTraits ; pot traits known (can be used to create a pot or be transfered to someone else)
pointInventory ; points the person possesses
potInventory ; pots the person possesses
pointPrestige ; value influenced by the traits possessed
potPrestige ; value given at random
huntResults ; list of last 10 hunting trip success rates
myAllies ; people with whom person has an alliance
gender ; 0 or 1 to represent two different genders (for gendered migration only)
]
Patch Variables Description {-} These variables are specific to patches.
patches-own [
environmentZone ; type of environment for patch
pOwnedBy ; which band owns the patch
]
Link Variables Description {-} These variables are specific to links.
links-own [
ally? ; switch to determine whether links will be "ON" or "OFF"
lastActive ; counter to record the last tick that the link was active if switched off
jaccard ; record the Jaccard distance
; keep track of every time people interact for the following categories
learning
visiting
hunting
trading
]
Procedures setup Description {-} Setup the model.
to setup
; sr:setup ; setup R extension -- comment out if not using adjcorr reporter
; reset everything
reset-timer ; track time
stop-inspecting-dead-agents ; remove any open turtle windows for now dead turtles
clear-all ; reset everything but time
reset-ticks ; reset time
; set a random seed for each run if option is true - this allows the model to record the seed for reproducibility
if randomizeSeed [ set randomSeed random 10000 ]
random-seed randomSeed ; set random seed
; track one person's activity
if recordPerson [
if file-exists? "results/outputActivity.csv" [ file-delete "results/outputActivity.csv" ] ; overwrites the file if it already exists
file-open "results/outputActivity.csv"
file-print "activity, tick"
file-close
set randomWho nofBands * 3 + nofBands + campPopulation
show ( word "tracked person is " randomWho ) ; prints who is followed in the command center
]
; create a new file to log into data at every 1000 tick if wanted
if HighResolutionData [
set fileHighRes (word "results/outputHighRes-" filenamebase ".csv")
file-open fileHighRes
]
; set some defaults values or parameter settings
set potBreakageRate 0.3
set pointsBreakageRate 0.3
set hideLinksInteractions true ; by default, links are not visible (for World clarity)
set aggregation false ; the simulation starts with people in their respective camps
set inventorySizePoints 10 ; people can only carry up to 10 points
set inventorySizePots 10 ; people can only carry up to 10 pots
set traitLetters ( list "a" "b" "c" "d" "e" "f" "g" "h" "i" "j" "k" "l" "m" "n" "o" ) ; creates the possible set of new variant letters
; procedures to create master tables, patches with environments, bands, camps, and people
setupMasterTables
if visibility [ setupVisibility ] ; assign visibility to specific variants if it affects transmission
setupPatches
setupBandsAndTerritories
setupCampsAndPeople
; create empty matrices (or adjacency tables) that log:
set campVisitingAdjacency matrix:make-constant count camps count camps 0 ; visiting interactions between camps at each wanted iteration
set campTradingAdjacency matrix:make-constant count camps count camps 0 ; trading interactions between camps at each wanted iteration
set campHuntingAdjacency matrix:make-constant count camps count camps 0 ; hunting interactions between camps at each wanted iteration
set campLearningAdjacency matrix:make-constant count camps count camps 0 ; learning interactions between camps at each wanted iteration
set campPointAdjacency matrix:make-constant count camps count camps 0 ; similarities in points between camps at each wanted iteration
set campPotAdjacency matrix:make-constant count camps count camps 0 ; similarities in pots between camps at each wanted iteration
end
setupPatches Description {-} Procedure to create patch environments (1 or 3).
to setupPatches
; define initial viewpoint as environmental zones
set currentView "zone"
; if chooser gradualEnvChange on - homogenous environmentZone that changes over time
ifelse gradualEnvChange [
ask patches [ set environmentZone 67 ] ; all patches get the same environmentZone
][
let pZones [ 47 67 37 ] ; colors for patches
; each pZone finds an unassigned patch, make it their own and create a buffer around it
foreach pZones [
a -> ask one-of patches with [ environmentZone = 0 ][
set environmentZone a
set pOwnedBy a ; borrowing this variable to speed up the code as it is not set yet
]
]
; sample patches to speed up code
let assignedPatches patches with [ pOwnedBy != 0 ] ; only the 3 patches that have an environment (will be the center of each environment)
; Each patch without an environmentZone takes the environmentZone of the closest patch with one
ask patches with [ environmentZone = 0 ][
set environmentZone [ environmentZone ] of min-one-of assignedPatches [ distance myself ]
]
]
ask patches [ set pcolor environmentZone ]
end
setupBandsAndTerritories Description {-} Procedure to create bands and their territories. Each band is separated by a minimum distance.
to setupBandsAndTerritories
; create invisible bands that will define the extent of the band's territory
; important - this code will separate each band by a specific radius, if the
; radius is too high then the code will print an error
let radiusBand 15 ; radius for min distance between band centers
let nBands ( count bands ) ; number of bands (update through the while loop)
set groupColors [ 14 24 34 44 54 64 74 84 94 104 114 124 134 ] ; define colors for bands
let tmpColor groupColors ; list of colors the bands take from (and remove the color they took so that each band has a unique color)
let counter 0 ; to break the while loop in certain circumstances
while [( nBands < nOfBands ) and ( counter < 10 )][
let potentialPatches 0
create-bands 1 [
set hidden? false
set shape "circle"
set size 3
set label who
set label-color white
; set color from list but make sure it is unique
set color ( item 0 tmpColor )
set tmpColor ( remove-item 0 tmpColor )
move-to one-of patches
; if there are other bands too close, the band will move somewhere else
if any? other bands in-radius radiusBand [
set potentialPatches patches with [ not any? bands in-radius radiusBand ]
ifelse any? potentialPatches [
move-to one-of potentialPatches
][
set counter 10 ; to kill the while loop
die
]
]
]
set nBands ( count bands )
set counter counter + 1
]
; warning message if radius is too high and not all bands can be placed
if nBands < nOfBands [ print "did not succeed in placing all bands" ]
; define default variable for pOwnedBy
ask patches [ set pOwnedBy 9999 ]
; setup bands as aggregation sites (and their inventories)
ask bands [
set pointAssemblage []
set potAssemblage []
set size 10
set label ""
set hidden? false
set shape "campsite"
ask patches in-radius radiusBand [ set pOwnedBy [ who ] of myself ] ; set band territory
]
let tmpPatches up-to-n-of 100 patches with [ pOwnedBy != 9999 ]
ask patches with [ pOwnedBy = 9999 ][
set pOwnedBy [ pOwnedBy ] of min-one-of tmpPatches [ distance myself ]
]
end
setupCampsAndPeople Description {-} Procedure to set up camps and their people.
to setupCampsAndPeople
; create camp sites within band territory
; follows the same rules above by not letting camps be too close
; but in this case all the camps are created first and then moved
let nOfCamps 3 ; number of camps per band
let radiusCamp 10 ; radius for camps to be spread out
ask bands [
hatch-camps nOfCamps [
; initialize lists to hold deposited items
set pointAssemblage []
set potAssemblage []
set size 6
set label ""
set hidden? false
set shape "campsite"
set ownedBy [ who ] of myself
set meanPrestige [] ; list of mean prestige of the camp per tick
set maxPrestige [] ; list of max prestige of the camp per tick
]
]
; set additional variables after all camps are created
ask camps [
set otherCamps other camps
set otherCampsMyband other camps with [ ownedBy = [ ownedBy ] of myself ]
create-campInteractions-with otherCamps [
set hidden? true
]
]
; move camps to their location
ask bands [
let tmpPatches patches with [ pOwnedBy = [ who ] of myself ] ; identifies the band's territory
ask camps with [ ownedBy = [ who ] of myself ][ move-to one-of tmpPatches ] ; asks their camp to move to their territory
]
; if the camp is too close to someone else, they try to move
ask camps [
if any? other camps in-radius radiusCamp [
let potentialPatches patches with [ pOwnedBy = [ OwnedBy ] of myself and not any? camps in-radius radiusCamp ]
ifelse any? potentialPatches [
move-to one-of potentialPatches
][
show ( word [ who ] of one-of other camps in-radius radiusCamp " is too close" )
]
]
]
; create user specified number of people per camp
ask camps [
hatch-people campPopulation [
set shape "person" ; make it look like a person
set size 2 ; change size
set homeBand [ pOwnedBy ] of patch-here ; set home band
set homeCamp myself ; set home camp
set originalhome [ who ] of myself; set original home camp
set pointInventory [] ; create empty list to hold point inventory
set potInventory [] ; create empty list to hold pot inventory
set pointTraits [] ; create empty list to hold point traits
set potTraits [] ; create empty list to hold pot traits
set huntResults [] ; list of hunting results
set myAllies ( turtle-set ) ; set this as turtle-set that will be filled later
set potPrestige precision ( random-normal 10 5 ) 2 ; each person gets a value from random distribution with mean 10 and standard deviation 5 (arbitrary)
; get its first known trait
let minTrait 1 ; set the min trait to create objects from (by default is 1)
if uniqueTraits [ set minTrait [ who ] of homeCamp * 10 ] ; but if we want different values per camp, this changes it trait set for each camp
luckyLeap ( word minTrait ) "point" ; this also creates one point at the beginning
luckyLeap ( word minTrait ) "pot" ; this does NOT create a pot at the beginning
]
; each camp asks 1/2 of its population to set their gender to 1
setupMasterTables Description {-} Observer procedure. Creates 3 tables, one for point traits, one for pot traits, and one for points’ environment (in which environment aare they fit?).
to setupMasterTables
set pointTraitsTable table:make ; Table that will take on the learned values of main traits
table:put pointTraitsTable "mainType" "variations" ; For easier visualization of which values are what
set potTraitsTable table:make ; Table that will take on the learned values of main traits
table:put potTraitsTable "mainType" "variations" ; For easier visualization of which values are what
set traitEnvironments table:make ; Table that will take on the learned values of main traits
table:put traitEnvironments "trait" "environment" ; For easier visualization of which values are what
end
setupVisibility Description {-} Observer procedure. Creates visibility table.
to setupVisibility
set visibilityTable table:make
end
go Description {-} This is the main procedure that calls the others. In this order: the visibility table is updated; the environment is changed (if necessary), the adjacency tables are updated, alliances are cleaned, innovations are made, people craft new object, aggregation occurs (if applicable), people do activities, prestige gets calculated, Traits are lost, traits and items are transmitted, objects are dropped, migration occurs (if applicable). Then, if the run has reached its end, metrics are outputted (if applicable).
to go
; update visibility table
if visibility [ updateVisibility ]
; change the environmentZone
if gradualEnvChange and random-float 100 < pEnvChange [ ; if the environment is set to change over time and the roll die is successfull...
print ( word "Climate change at tick: " ticks )
let newZone random 100 ; random number (integer)
ask patches [ set environmentZone newZone ]
]
; Update adjacency tables when in need of high resolution data (need to figure out timing, right now, it updates every 1000 ticks)
if HighResolutionData and ticks > 0 and remainder ticks 1000 = 0 [ ; this is active only if we need it
compileArtifactSimilarities
compileCampInteractions
exportFineTunedData
; clear camp interaction links and camp assemblages
ask campInteractions [
set learning 0
set visiting 0
set hunting 0
set trading 0
]
ask camps [
set potAssemblage []
set pointAssemblage []
]
]
; people remove alliancess if they are over the limit, before doing anything else
ask people [ cleanUpAlliances ]
; Innovate at very low probability
ask people [
if random-float 1000 < pNewVariant [ createVariant "point" ] ; if die is successful, create a new point variant
if random-float 1000 < pNewVariant [ createVariant "pot"] ; if die is successful, create a new pot variant
if random-float 1000 < pLuckyLeap [ luckyLeap ( getNextNumber "point" ) "point" ] ; if die is successful, create a new point type
if random-float 1000 < pLuckyLeap * 2 [ luckyLeap ( getNextNumber "pot" ) "pot" ] ; if die is successful, create a new pot type (not as often as points)
]
; produce new items if needed
ask people [
craftItem "point" true ; the "true" calls a probability within the craftItem procedure, so that they create new objects at a certain rate
craftItem "pot" true ; the "true" calls a probability within the craftItem procedure, so that they create new objects at a certain rate
]
; control aggregation
if ticks > 0 and remainder ticks aggregationFreq = 0 [ ; If it is time to aggregate then aggregate
set aggregation true
aggregate
]
ifelse aggregation = true [
ifelse aggregationCounter = aggregationDuration [ ; If the aggregation duration period has passed then disperse
set aggregationCounter 0
disperse
][
set aggregationCounter ( aggregationCounter + 1 ) ; if they are still aggregating, update the counter and push people to participate in some activities (hunting only, as they are already visiting)
ask people [
changeActivity "aggregating"
if random-float 100 < pHunting [ ; random decision to hunt based on slider specified odds
changeActivity "hunting"
hunt ; procedure to create travel groups and destinations
]
]
]
][
; trading, hunting and visiting only if not aggregated
ask people [
if random-float 100 < pNewAlliance [ createAlliance ] ; procedure to create alliances who can be outside one's band
]
; do activities (the same person can both go hunting and go visiting in the same tick
ask people [
if random-float 100 < pHunting [ ; random decision to hunt based on slider specified odds
changeActivity "hunting"
hunt ; procedure to create travel groups and destinations
]
if random-float 100 < pVisiting and any? myAllies [ ; you cannot go visit anyone if you dont have any allies
changeActivity "visiting"
visit ; procedure to create travel groups and destinations
]
]
]
if recordPrestige [ ; if recordPrestige switch is on
ask camps [ recordCampPrestige ] ; camps record the prestige variables from their people (mean and max)
]
; give each person chance to lose a trait (forget what they know) - equal chance to lose pot or point trait
ask people [
if random-float 100 < pLoss [
ifelse random 2 = 0 [ ; roll the dice to see if they will lose a point or pot trait that time
set pointTraits ( loseTrait pointTraits )
][
set potTraits ( loseTrait potTraits )
]
]
]
; cultural transmission of known traits and exchange of objects
ask people [
transmitCulture
tradeItems
]
; drop object from inventory at random
ask people [
dropItem "point" true
dropItem "pot" true
]
if ticks >= 1000 [
ask people [
; Migration (or exogamy) procedure
ifelse genderedMigration = TRUE [ ; when genderedMigration occurs, only people with gender 1 will move
changeActivity Parameters {-} Description needs to be a string identifying the new activity. It can be “hunting”, “learning”, “trading”, or “visiting”.
Description {-} Change the activity of the person and if recordPerson is true then records the activity.
to changeActivity [ description ]
set activity description
if recordPerson and who = randomWho [ ; if tracking someone, open the file, write the new activity with timestamp, and close the file
file-open "results/outputActivity.csv"
file-print ( word ( activity ) "," ( ticks ))
file-close
if debug [ show ( word activity " on tick " ticks )] ; if debugging, show the new activity and timestamp
]
end
aggregate Description {-} Observer procedure. Tells people to go to their band’s spot, which serves as the aggregation camp. All people in the same camp add 1 visiting to their counter. People can then create alliance and exchange items with people in the same camp.
to aggregate
if debug [ print "aggregate" ] ; for debugging
; people move to their band
ask people [ move-to band homeBand ]
ask interactions with [ link-length = 0 ][ ; this focuses only on links between people who are in the same aggregation camp.
set visiting visiting + 1 ; this assumes that everyone aggregating with one another has some sort of social contact (counted only once, not every month of aggregation)
updateCampInteractions "visiting" end1 end2 ; updates the camp counter
disperse Description {-} Observer procedure. Called when the aggregation is over. Tells people to return to their home camp.
to disperse
if debug [ print "disperse" ] ; for debugging
set aggregation false ; reset the variable
set aggregationCounter 0 ; restart the counter
ask people [ returnToCamp ] ; people go back to their camp
end
hunt Description {-} People procedure. People go hunting with band mates.
to hunt
let party ( turtle-set self ) ; create turtle-set of hunting or visiting party
let activityLocation 0
let availablePeople 0
; sets the activity location (when aggregating, it is nearby, when not aggregating, it is simply within the band territory)
set activityLocation ifelse-value aggregation [ one-of patches in-radius 15 ][ one-of patches with [ pOwnedBy = [ homeBand ] of myself ]]
move-to activityLocation ; move to the set location
; identify the hunting party (always band mates)
set availablePeople ( turtle-set other people with [ homeBand = [ homeBand ] of myself ])
ask up-to-n-of huntingGroupSize availablePeople [ ; tell a specific number of those people to join the hunt
changeActivity "hunting"
set party ( turtle-set party self )
move-to activityLocation
; add to counters for activity (each member of the party does it once with the person starting the hunt, so everythin counts only once).
updateCampInteractions "hunting" self myself ; this updates the counter of the camps
ask interaction [ who ] of self [ who ] of myself [ set hunting hunting + 1 ] ; update the people counter
]
ask party [
calculatePrestige ; this assumes that the hunt occurred and everyone evaluated how well their knowledge fit the current environment
dropItem "point" true ; this calls the dropItem procedure to see if the people lose one of their points ("true" means that it will happen with a certain probability)
tradeItems ; people hunting together have the opportunity to trade points with one another
transmitCulture ; people hunting together also share knowledge with one another
]
; when all members of the party have had the chance to do all activities, they all return home
ask party [ returnToCamp ] ; this needs to be separated so that people don't leave before others have had the time to share with them
end
createAlliance Description {-} People procedure. People create a new alliance. When aggregating, they choose someone new (who is not already their ally) in the same camp. At other times, they choose someone new from another camp, favoring people from nearesst camps.
to createAlliance
if debug [ print "createAlliance" ] ; for debugging
let newAlly 0 ; temp
changeActivity "making allies"
ifelse aggregation = true [
if any? other people-here [ ; When aggregating, they choose someone in the same band who is not already their ally, who can get a new ally, and who is not from the same camp
set newAlly ( one-of other people-here with [ not member? self myAllies and count myAllies < maximumAlliesN and homeCamp != [ homeCamp ] of myself ])
]
][ ; when not aggregating, they can form an alliance with someone outside their band (they favor people from nearby camps, though)
let allOtherCamps []
let distAllOtherCamps []
ask [ otherCamps ] of homeCamp [ ; this fills in the two temporary lists above (the camps are in the same order)
set allOtherCamps lput self allOtherCamps ; list of camps other than homeCamp
set distAllotherCamps lput precision ( distance [ homeCamp ] of myself ) 0 distAllOtherCamps ; list of their distance from homeCamp (without decimals)
]
let maxDist max distAllOtherCamps ; the furthest camp
set distAllOtherCamps ( map [ x -> abs( x - maxDist )] distAllOtherCamps ) ; inverts the distance so that the closest camp has the highest probability
let pairs ( map list allOtherCamps distAllOtherCamps ) ; Pairs the camps with their probabilities
let selectedCamp first rnd:weighted-one-of-list pairs [[ p ] -> last p ] ; chooses one camp at random, but favoring closest ones
set newAlly one-of people with [ homeCamp = selectedCamp and not member? self myAllies and count myAllies < maximumAlliesN ] ; select one person from that camp who is not yet its ally
]
if is-turtle? newAlly [ ; in case no newAlly was found
ask interaction [ who ] of self [ who ] of newAlly [ set ally? true ] ; activate the alliance
; reset the agentset to only allies from active links and tell the new ally to reset theirs as well
set myAllies ( turtle-set [ other-end ] of my-interactions with [ ally? ])
cleanUpAlliances Description {-} People procedure. People count the number of allies they have and drop one or more if they have more than the limit. Allies with strong ties have less chance to get dropped than allies with weak ties.
to cleanUpAlliances ; to make sure no one ever has more allies than allowed
if debug [ print "cleanUpAlliances" ] ; for debugging
; remove alliances with probability inversely proportional to the number of times you've visited that ally
if count myAllies > maximumAlliesN [ ; if you have more allies than maximum
let allyVisitsList ([ visiting ] of my-interactions with [ ally? ]) ; create a temporary variable with the list of all links' visiting counters
set allyVisitsList ( map [ x -> x + 1 ] allyVisitsList ) ; add one to avoid dividing by zero when no visit have yet occurred
let sumVisits sum allyVisitsList ; calculate the sum of the number of visits to other people
let allyProbList allyVisitsList ; create a copy list that will take in the calculated probabilities of each visit value
foreach range length allyProbList [ itemN -> ; itemN takes the value of each item in the allyProbList list
let currentVisitN item itemN allyProbList ; Reports the value of each link (number if visits)
set allyProbList replace-item itemN allyProbList ( 1 / ( currentVisitN / sumVisits )) ; the probability that the link will be lost will be inversely related to network counter
]
let pairs ( map list allyVisitsList allyProbList ) ; Pairs the link counters with their 0-1 probabilities
let chosenVisitN rnd:weighted-one-of-list pairs [[ p ] -> last p ] ; rnd:weighted-one-of-list reports a random agent from the 'pairs' list, taking their probability (p) into account
let chosenLink one-of interactions with [ visiting = ( first chosenVisitN ) - 1 ] ; The link with the chosen cost (see above) is used as an alliance.
ask chosenLink [
set ally? false ; link becomes inactive (alliance lost)
set lastActive ticks
]
set myAllies ( turtle-set [ other-end ] of my-interactions with [ ally? ]) ; reset agentset to only active links after making some inactive
visit Description {-} People procedure. People go visit one of their allies.
to visit
let party ( turtle-set self ) ; will take on all people participating in the activity
move-to one-of myAllies ; moving to the ally should be the same as moving to its camp.
; as we are choosing amongst allies, a camp with multiple allies will have more chances to be chosen (same as when choosing amongst camps)
let alliesHere people-here with [ member? myself myAllies ] ; identifies all the allies within one camp
ask alliesHere [ ; visit only allies in destination camp
changeActivity "visiting"
set party ( turtle-set party self ) ; add themselves to the party (to do activities later on)
; add to counters for activity (each person visited does it once with the person visiting, so everythin counts only once).
updateCampInteractions "visiting" self myself ; this updates the counter of the camps
ask interaction [ who ] of self [ who ] of myself [ set visiting visiting + 1 ] ; this is done only by the allies, so it counts only once.
]
ask party [ ; then everyone does a few activities
dropItem "point" true ; they will drop a point with a certain probability
dropItem "pot" true ; they will drop a pot with a certain probability
tradeItems ; they will drop a pot with a certain probability
transmitCulture
]
returnToCamp ; this is separate from the "ask alliesHere" because allies stay in their camp. Only ego goes home.
end
calculatePrestige Description {-} People procedure. People take their inventory of points and calculate how many of them are fit in the current environment. This gets added to their huntResults variable. They take the average of the last 10 of those values as their prestige.
to calculatePrestige
if debug [ print "calculatePrestige" ] ; for debugging
if length pointInventory > 0 [ ; prestige is only calculated on the points in the inventory
let onlyPoints ( map [ x -> last x ] pointInventory ) ; points in pointInventory have the tick and then the object. This filters out the tick to keep only the point
let toolList ( map [ x -> table:get traitEnvironments x ] onlyPoints ) ; this extracts the environment in which the trait of each object in inventory was invented
let myZone [ environmentZone ] of patch-here ; variable to determine environmental zone of location
let usefulToolList filter [ i -> i = myZone ] toolList ; keeps only the points' environments that fit the current environment
let fitness length usefulToolList / length toolList ; the ratio of points the person has that is good for the current environment
set huntResults lput fitness huntResults ; add the latest fitness to the list of huntResults
if length huntResults > 10 [ set huntResults remove-item 0 huntResults ] ; this list can only have the latest 10 hunts, so the oldest 11st gets dropped, when appropriate
set pointPrestige mean huntResults ; pointPrestige is the average of the last 10 hunts' fitness
]
end
tradeItems Description {-} People procedure. Two people in the same place exchange one of their point and/or pot with one another. When they do so, they add 1 to their “trading” counter.
to tradeItems
if debug [ print "tradeItems" ] ; for debugging
let egoCurrentActivity activity ; so they can go back to their activity after the exchange
let partnerCurrentActivity 0 ; for the partners to remember what they were doing too
let theirGift 0 ; hold their present
;; POT ;;
if random-float 100 < pTrading and activity != "hunting" and length potInventory > 0 [ ; pot trading cannot occur during a hunt or if you have nothing to trade
let tradingPartners other people-here with [ length potInventory > 0 ] ; find other people who also have pots to trade
if any? tradingPartners [ ; the trade occurs only if they found someone who had something to trade
changeActivity "trading pots"
let myGift one-of potInventory ; choose gift
let posGift position myGift potInventory ; identify the position of the item to be given
set potInventory remove-item posGift potInventory ; remove that specific object from the inventory
ask one-of tradingPartners [ ; exchange with one person
set partnerCurrentActivity activity ; identify the activity the partner was doing
changeActivity "trading pots" ; change activity (for debugging, mostly)
set theirGift one-of potInventory ; they choose one pot to trade
set posGift position theirGift potInventory ; identify the position of the item to be given
set potInventory remove-item posGift potInventory ; remove that specific object from the inventory
set potInventory lput myGift potInventory ; then add the new item
changeActivity partnerCurrentActivity ; reset to previous activity
; update the interaction counters for camps and people links
updateCampInteractions "trading" self myself ; update the counter of the camps
ask interaction [ who ] of self [ who ] of myself [ set trading trading + 1 ] ; update the people counter
]
set potInventory lput theirGift potInventory ; add the partner's gift to inventory
changeActivity egoCurrentActivity ; resets to previous activity
]
]
;; POINT ;;
if random-float 100 < pTrading and length pointInventory > 0 [ ; point trading can only happen if there are points to trade (which should always be the case)
let tradingPartners other people-here with [ length pointInventory > 0 ] ; find other people who also have points to trade
if any? tradingPartners [ ; the trade occurs only if they found someone on the same spot who had something to trade
changeActivity "trading points"
let myGift one-of pointInventory ; choose gift
let posGift position myGift pointInventory ; identify the position of the item to give
set pointInventory remove-item posGift pointInventory ; remove that specific object from the inventory
ask one-of tradingPartners [ ; exchange with one person
set partnerCurrentActivity activity ; identify the activity the partner was doing
changeActivity "trading points" ; change activity (for debugging, mostly)
set theirGift one-of pointInventory ; they choose one point to trade
set posGift position theirGift pointInventory ; identify the position of the item to give
set pointInventory remove-item posGift pointInventory ; remove that specific object from the inventory
set pointInventory lput myGift pointInventory ; then add the new item
changeActivity partnerCurrentActivity ; reset to previous activity
; update the interaction counters for camps and people links
updateCampInteractions "trading" self myself ; update the counter of the camps
ask interaction [ who ] of self [ who ] of myself [ set trading trading + 1 ] ; update the people counter
]
set pointInventory lput theirGift pointInventory
changeActivity egoCurrentActivity ; resets to previous activity
]
]
end
dropItem Parameters {-} itemType can be “point” or “pot” so it can be called individually. prob? is a boolean used to increase flexibility of the procedure (so it can be called by breaking the object as well as if they are over their inventory limit).
Description {-} People procedure. People drop one itemType with a probability set by the object’s breakageRate (if prob? is FALSE) or automatically if they are over their hand limit (prob? TRUE). Hunters are much more likely to drop a point when hunting. If it was the last point they had, they create a new one. If the drop occurs at a camp, that camp records the object’s characteristics.
to dropItem [ itemType prob? ]
if debug [ print "drop-item" ] ; for debugging
ifelse itemType = "point" [
;; POINT;;
; the threshold for breakage is lower if the agent is hunting
let threshold ifelse-value activity = "hunting" [ pointsBreakageRate ][ pointsBreakageRate / 20 ]
; the first part of the line below is if dropItem is called during a hunt or GO, the second is if it is called because the person has reached its hand limit
if ( prob? and random-float 100 < threshold ) or ( not prob? ) [
if ( length pointInventory > 0 )[ ; only drop if an item exists
changeActivity "dropping point"
let site ( turtle-set camps-here bands-here ) ; identifies if the person is at a camp or band
; if there is a site here, it records the dropped point (the oldest of the inventory)
ask up-to-n-of 1 site [ set pointAssemblage lput ( item 0 [ pointInventory ] of myself ) pointAssemblage ]
set pointInventory remove-item 0 pointInventory ; remove the dropped item from inventory
; if the inventory is empty, create a new point (hunters always need to have one point in their inventory)
if ( length pointInventory = 0 )[ set pointInventory lput ( list ticks one-of pointTraits ) pointInventory ]
]
]
][
;;POT ;;
; same procedure as for points, except activity does not affect breakage rate
if ( prob? and random-float 100 < potBreakageRate ) or ( not prob? )[
if ( length potInventory > 0 )[ ; only drop if an item exists
changeActivity "dropping pot"
let site ( turtle-set camps-here bands-here ) ; identifies if the person is at a camp or band
; if there is a site here, it records the dropped pot (the oldest of the inventory)
ask up-to-n-of 1 site [ set potAssemblage lput ( item 0 [ potInventory ] of myself ) potAssemblage ]
set potInventory remove-item 0 potInventory ; remove that pot from the inventory
]
]
]
end
returnToCamp Description {-} People procedure. Agents move back to home camp or to the aggregation location (can be called by disperse or doActivity at the end of a hunt).
to returnToCamp
if debug [ print "returnToCamp" ] ; for debugging
; choose where to move to depending on what's happening right now
updateCampInteractions Parameters {-} interactionType is a string that can be “learning”, “hunting”, “trading”, or “visiting”. agentA and agentB are two different people doing the interaction.
Description {-} Observer procedure. Update the camp interaction link every time agents interact
to updateCampInteractions [ interactionType agentA agentB ]
let campA [ who ] of [ homeCamp ] of agentA ; gets the number of agentA's camp
let campB [ who ] of [ homeCamp ] of agentB ; gets the number of agentB's camp
if campA != campB [ ; Because people interacting within one camp will not be counted here
ask campInteraction campA campB [ ; the remainding lines update the appropriate counter
ifelse interactionType = "learning" [
set learning learning + 1
][
ifelse interactionType = "hunting" [
set hunting hunting + 1
][
ifelse interactionType = "trading" [
set trading trading + 1
][
set visiting visiting + 1
]
]
]
]
]
end
craftItem Parameters {-} itemType can be “point” or “pot” so it can be called individually. prob? is a boolean used to increase flexibility of the procedure (T: it happens with a certain probability, F: it happens automatically).
Description {-} People procedure. People have a certain probability to create an object (pot or point). When created, the object records the current tick as well as its type (one of the types its maker knows).
to craftItem [ itemType prob? ]
if debug [ print "craftItem" ] ; for debugging
; Every tick, each individual has a small chance to produce a new object (when prob? = true).
; Every time people acquire a new knowledge (when prob? = false), the create a new object
let diceRoll random-float 100 ; for easier testing and troubleshooting
;; POINTS ;;
ifelse itemType = "point" [
if ( prob? and diceRoll < pNewObject ) or ( not prob? )[
changeActivity "making a point"
let fitTrait getFitTrait ; chooses the trait among the fitter ones (based on environment)
set pointInventory lput ( list ticks fitTrait ) pointInventory ; add the point to its inventory (includes the tick number it was made)
]
; if that new creation pushed them above their toolkit limit, they drop the oldest point
if length pointInventory > inventorySizePoints [ dropItem "point" false ]
][
;; POTS ;;
if ( prob? and diceRoll < pNewObject ) or ( not prob? ) [
changeActivity "making a pot"
set potInventory lput ( list ticks one-of potTraits ) potInventory ; add the pot to its inventory (includes the tick number it was made)
]
; if that new creation pushed them above their toolkit limit, they drop the oldest pot
if length potInventory > inventorySizePots [ dropItem "pot" false ]
]
end
migrate Description {-} People procedure. One person chooses another camp and trade place with someone from that camp. Migration can be gendered, in which case, only people with gender 1 trade places.
to migrate
if debug [ print "migrate" ] ; for debugging
changeActivity "migrating"
let oldCamp homeCamp ; store old camp for swapping
let possibleCamps 0 ; will take on the ID of possible camps
ask homeCamp [ set possibleCamps other camps ]
set homeCamp one-of possibleCamps ; choose one of the other camps at random as its new home
set homeBand [ ownedBy ] of homeCamp
if aggregation = false [ move-to homeCamp ] ; when not aggregating, move to the new home (when aggregating, stay within the aggregation camp for now)
; find someone who can take its old place
let myID homeCamp ; the new camp
let myNewCampmates people with [ homeCamp = myID ] ; all the people in the new home
; select the migration partner (if gendered, only gender 1 can be selected)
recordCampPrestige Description {-} Camp procedure. Camps calculate the mean and max prestige of their people.
to recordCampPrestige
let myPeople people with [ homeCamp = myself ] ; identify the people who belong to this camp
let meanPrestigeHere precision (( sum [ pointPrestige ] of myPeople ) / campPopulation ) 2 ; Mean prestige of its people (rounded to 2 decimals)
set meanPrestige lput meanPrestigeHere meanPrestige ; adds that value to a list
let maxPrestigeHere ( max [ pointPrestige ] of myPeople ) ; Max prestige of its people
set maxPrestige lput maxPrestigeHere maxPrestige ; add that value to a list
end
luckyLeap Parameters {-} objectType can be “point” or “pot”. nextTrait is the trait that will be innovated (added to the list).
Description {-} People procedure. People create a new object (point or pot) with the trait set to the given nextTrait. They then create a new object.
to luckyLeap [ nextTrait objectType ]
if debug [ print "luckyLeap" ] ; for debugging
ifelse objectType = "point" [
changeActivity "inventing point"
; create a new main POINT type
set pointTraits lput nextTrait pointTraits ; add the next number to the list
set pointTraits remove-duplicates pointTraits ; make sure there are no duplicates
if not table:has-key? pointTraitsTable nextTrait [
table:put pointTraitsTable nextTrait [] ; if this is a brand new type, create a new entry in the table
]
; assign the new type to its environment if not already invented
; but only for points as pots' fitness does not change with environment
if not table:has-key? traitEnvironments nextTrait [
let numberHere [ environmentZone ] of patch-here ; record the environment where the type is created
table:put traitEnvironments nextTrait numberHere ; list recording environment in which each main axis type was created
]
craftItem "point" false ; then create a new point
][
changeActivity "inventing pot"
; create a new main POT type
set potTraits lput nextTrait potTraits ; add the next number to the list
set potTraits remove-duplicates potTraits ; make sure there are no duplicates
if not table:has-key? potTraitsTable nextTrait [
table:put potTraitsTable nextTrait [] ; if this is a brand new type, create a new entry in the table
]
craftItem "pot" false ; then create a new pot
]
end
createVariant Parameters {-} objectType can be “point” or “pot”.
Description {-} People procedure. Add a variant to a mainType already known by ego. Then create an object.
to createVariant [ objectType ]
if debug [ print "addVariant" ]
changeActivity word "inventing variant" objectType
; set up temporary variables
let reachedLimit? false ; will help stop the addition of variation when it reached the limit ("o" for sequential, and all letters used for non-sequential
let newTypeVariation "" ; temp variable that will be filled later
let newVariation "a" ; by default, the newVariation will be "a" so that it can be used if there is not already a variant
; select one known trait of the correct type
let traitKit ifelse-value objectType = "point" [ pointTraits ][ potTraits ] ; define a temporary set of the relevant traits based on objectType
let traitA one-of traitKit ; choose one known trait at random.
let traitNumber getTraitNumber traitA ; focus only on the main number (the main type)
let lastVariant getLastVariant traitA objectType ; gets the last variant that person has of the main type
if lastVariant != "" [ ; if ego knew some variants...
ifelse lastVariant != "o" [ ; ... and that variant is not the last possible one (there can only be 15 variations on each type)...
set newVariation nextLetter lastVariant ; ...get the next possible letter alphabetically
][
set reachedLimit? true ; if the type has reached its variation limit ("o"), the following block of code will not get called
]
]
; After the new type variant has been chosen (if not passed the limit), this updates the traitkit and the tables.
if reachedLimit? = false [ ; only if the type has variations remaining (including "a") does it create a new one
set newTypeVariation ( word traitNumber "-" newVariation ) ; add a letter to the type's variation
; add the new trait to the relevant trait list
ifelse objectType = "point" [
; POINT ;
set pointTraits lput newTypeVariation pointTraits ; add the new trait to the person's known traits list
; determine if the point type already has an environment attached to it, and adds it if not
if not table:has-key? traitEnvironments newTypeVariation [
let numberHere [ environmentZone ] of patch-here ; record the environment number
table:put traitEnvironments newTypeVariation numberHere ; list recording environment in which each type was created
]
; Add this new variant to the master table
let listTraitKits table:get pointTraitsTable ( word traitNumber ) ; focus on all variations already done for this type
let newTraitKit remove-duplicates ( sentence listTraitKits newVariation ) ; add the new type variation to that type's list only if it's not already there
table:put pointTraitsTable ( word traitNumber ) newTraitKit ; replace the list in the table by that updated one
craftItem "point" false ; then create a new point
][ ; POT ;
set potTraits lput newTypeVariation potTraits ; ads the new trait to the person's trait list
; Add this new type variant to the master table
let listTraitKits table:get potTraitsTable ( word traitNumber ) ; focus on all variations already done for this type
let newTraitKit remove-duplicates ( sentence listTraitKits newVariation ) ; add the new type variation to that type's list only if it's not already there
table:put potTraitsTable ( word traitNumber ) newTraitKit ; replace the list in the table by that updated one
]
]
end
loseTrait Parameters {-} traitsList represents the list of known trait by ego (either potTraits or pointTraits).
Return {-} Returns the trait or traits that ego remembers after forgetting one trait (and its variant if it was a main type).
Description {-} People procedure. People forget one of the traits they know, at random. If a main trait is forgotten, ego forgets its variants as well.
to-report loseTrait [ traitsList ]
let lossTrait 0
let rememberedTraitsList traitsList ; by default, if they do have enough traits to lose, they keep what they had
if length traitsList > 1 [ ; do not remove trait if only one trait exists
changeActivity "losing trait"
set lossTrait one-of traitsList ; choose one of the known traits at random
let listToFilter map [ i -> ifelse-value member? lossTrait i [ "loss" ][ i ]] traitsList ; creates a list with "loss" if the trait needs to be lost, and the trait if it is kept
set rememberedTraitsList filter [ i -> i != "loss" ] listToFilter ; filters the list to keep only the traits remembered (remove all "loss")
; if they forgot the only main type they had and its variants (would happen a lot at first), they do not forget anything
if length rememberedTraitsList = 0 [ set rememberedTraitsList traitsList ]
]
report rememberedTraitsList
end
getNextNumber Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the max number already known + 1.
Description {-} People procedure. Called with LuckyLeap only. It finds the next number (type) to innovate.
to-report getNextNumber [ objectType ] ; People reporter, reports the next number for the objectType (can be "pot" or "point")(for luckyLeap)
let tempListTraits [] ; temp variable to allow using this reporter for both points and pots
ifelse objectType = "point" [
set tempListTraits pointTraits ; populate it with point traits
][
set tempListTraits potTraits ; populate it with point traits
]
let listOfNumbers ( map [ trait -> getTraitNumber trait ] tempListTraits ) ; transform the traits into their main trait number
let nextNumber ( max listOfNumbers ) + 1
report ( word nextNumber ) ; get the max number already known and add 1 for the next (as a string)
end
getFitTrait Return {-} Returns one point type with a relatively high fitness, given the current environment.
Description {-} People procedure. Looks at the person’s pointTraits’s fitnesses in the current environment, and then chooses one of the traits from that list (fit traits have more chance to be selected than non-fit traits).
to-report getFitTrait
if debug [ print "getFitTrait" ] ; for debugging
let myZone [ environmentZone ] of patch-here ; environmental zone of the current location
let traitFitnessList ( map [ x -> table:get traitEnvironments x ] pointTraits ) ; this extracts the environment of each point trait known
; for each point trait known, it replaces the trait with its fitness (5 if the trait was invented in the same environment as myZone, 2 if not)
set traitFitnessList map [ x -> ifelse-value ( x = myZone ) [ 5 ] [ 2 ]] traitFitnessList
let pairs ( map list pointTraits traitFitnessList ) ; Pairs the pointTraits counters with their current fitness
let fitTrait first rnd:weighted-one-of-list pairs [[ p ] -> last p ] ; chooses one trait, favoring the fitter ones (based on environment)
report fitTrait ; report the chosen trait
end
deconstructingTrait Parameters {-} trait is a trait known by the person calling this reporter. It can be either a point or a pot trait (not the object made with it).
Return {-} Returns the last part of the trait in its correct format (either a letter as a string or a number as a number).
Description {-} People procedure. Used to define if a trait is a main type or a type with a variant, so the person innovating onto it can add the correct variant (called in createVariant).
to-report deconstructingTrait [ trait ] ; people procedure called in createVariant to set the last item to its type (string or number) without erroring-out
let traitWithType 0 ; temporary variable
carefully [ ; this allows separating the trait without creating errors (if the trait is a main type, it should not be read-from-string).
set traitWithType read-from-string last trait ; that's when the trait is a main type, which is transformed into a number (was a string of the number)
][
set traitWithType last trait ; that's when the trait has a variant, which is a string
]
if traitWithType = e [ set traitWithType "e" ] ; e is read as a number, so this fixes the error it would cause.
report traitWithType
end
getTraitNumber Parameters {-} trait is a trait known by the person calling this reporter. It can be either a point or a pot trait (not the object made with it).
Return {-} Returns the number of the trait as a number (main type).
Description {-} People procedure. Gets the main type of the object (without variants).
to-report getTraitNumber [ trait ]
; first identify the position of the "-" (when there is a variant) or the length of the trait (when it is a main type)
let p position "-" trait
if p = false [ set p length trait ] ; if there was no "-", p became "false," so this gives it the length of the whole trait
; then extract the part of the trait up to p (a.k.a., up to the "-" or the whole trait)
let traitN substring trait 0 p ; removes the variation letter and keeps only the number
set traitN read-from-string traitN ; transform the string into its number form
report traitN
end
nextLetter Parameters {-} currentLetter is the letter of the variant upon which the person will innovate.
Return {-} returns the next letter alphabetically.
Description {-} People procedure. Helps find the next letter when innovating on a trait with variant.
to-report nextLetter [ currentLetter ]
let thisPosition position currentLetter traitLetters ; finds the position of the currentLetter in the alphabetical list set earlier
let nextPosition thisPosition + 1 ; get the next position
report item nextPosition traitLetters ; returns the letter following currentLetter alphabetically
end
traitVisibility Parameters {-} trait is one trait known by ego (can be for points or pots). The trait needs to be one with variant.
Return {-} Returns the visibility assigned to the trait provided.
Description {-} People procedure. Identifies the visibility of a certain trait. Called by checkVisibility only for traits with variants.
to-report traitVisibility [ trait ]
let result 0
carefully [set result table:get visibilityTable trait][; get the assigned visibility of that variant
set result 0
]
report result
end
getMostlyVisibleTraits Parameters {-} traitList can be a pointTraits list or the potTraits list.
Return {-} Returns a list of some of the most visible traits from the traitList provided.
Description {-} People procedure. Samples traits from the traitList with repeat. Visible traits have higher chances to get selected.
to-report getMostlyVisibleTraits [ traitList ]
; ; all main types are highly visible, so we can automatically separate them to assign visibility
; let mainTypes filter [ x -> is-number? deconstructingTrait x ] traitList
;
; let variantTraits filter [ x -> is-string? deconstructingTrait x ] traitList ; filters only the traits with variants
; for variants
let traitsVisibility ( map [ x -> traitVisibility x ] traitList ) ; computes each variants' visibility
let traitsPairedVisibility ( map list traitList traitsVisibility ) ; pairs each variant with its visibility
; filter out low visibility traits
let traitsVisibile []
foreach traitsPairedVisibility [
x -> if item 1 x = 3 [
set traitsVisibile insert-item 0 traitsVisibile item 0 x
]
]
; combine mainTypes and high-visibility variants
; let visibleTraits sentence mainTypes traitsVisibile
report traitsVisibile
end
getVariantsOfMainType Parameters {-} trait is the trait considered, objectType can be “point” or “pot”.
Return {-} Returns a list of all the variants of the trait’s mainType already known by ego (NOT all variants created in the simulation).
Description {-} People procedure. Identifies all the variants of trait’s mainType known by ego.
to-report getVariantsOfMainType [ trait objectType ]
; first get the mainType number
let mainTypeNumber getTraitNumber trait
; create a temp variable that takes on the traits considered here, based on what the objectType is
let traits ifelse-value objectType = "point" [ pointTraits ][ potTraits ]
; then produce a list of that type's variants (from the relevant list)
let mainTypeVariants map [ i -> ifelse-value getTraitNumber i = mainTypeNumber [ deconstructingTrait i ][ "" ]] traits
; This list now has the variants of the main tool or the main tool itself, and "" for all other tools
report mainTypeVariants
end
getLastVariant Parameters {-} trait is the trait considered, objectType can be “point” or “pot”.
Return {-} Returns the last variant ego knows of the type of the trait provided (variant based on alphabetical order).
Description {-} People procedure. Identifies the last variants ego knows of the type of the given trait.
to-report getLastVariant [ trait objectType ]
let variants getVariantsOfMainType trait objectType ; call the reporter to get a list of ego's traits that contains the mainType, its variants, and "" for anything else
let positionVariants map [ i -> position i traitLetters ] variants ; get the position of each variant in the alphabet (returns false if the variant is not part of the alphabet)
let maxVariant ""
carefully [ ; if the type of the trait did not have a variant, this prevents an error and goes to get any other variants that would be for that type
let maxVariantPos max positionVariants ; get the max position (furthest along the alphabet)
set maxVariant item maxVariantPos traitLetters ; get the variant associated with that position in the alphabet
][ ] ; there is nothing in the second bracket, because if there is an error in the first, nothing changes (maxVariant remains "")
report maxVariant
end
updateVisibility Description {-} Observer procedure. Assigns visibility of 1 (low visibility) or 2 (high visibility) to a trait or variant.
to updateVisibility
let a remove-duplicates reduce sentence [potTraits] of people
let b remove-duplicates reduce sentence [pointTraits] of people
let allTraits remove-duplicates sentence a b
foreach allTraits [
k ->
if not table:has-key? visibilityTable k [
table:put visibilityTable k (random 2) + 1 ; if this is a brand new type, create a new entry in the table
]
]
end
transmitCulture Description {-} People procedure. Calls the appropriate procedure to transmit traits based on the transmissionRate, learning method, and object selected.
to transmitCulture
if debug [ print "transmitCulture" ] ; for debugging
let newTrait 0
;; POT ;;
if random-float 100 < transmissionRate [ ; If they are learning from someone
changeActivity "learning pots"
; add one learned trait (the source varies based on the learningMethod)
set newTrait ifelse-value random-float 100 < learningMethod [ transmitPrestige "pot" ][ transmitConformism "pot" ]
if newTrait != 0 [
set potTraits lput newTrait potTraits
craftItem "pot" false ; then create a new pot
set potTraits remove-duplicates potTraits ; make sure that no duplicates are added through transmission
]
]
;; POINT ;;
if random-float 100 < transmissionRate [ ; If they are learning from someone
changeActivity "learning points"
; add one learned trait (the source varies based on the learningMethod)
set newTrait ifelse-value random-float 100 < learningMethod [ transmitPrestige "point" ][ transmitConformism "point" ]
if newTrait != 0 [
set pointTraits lput newTrait pointTraits
craftItem "point" false ; then create a new point
set pointTraits remove-duplicates pointTraits ; make sure that no duplicates are added through transmission
]
]
end
transmitConformism Parameters {-} objectType can be “point” or “pot”
Return {-} Returns the most common trait
Description {-} People procedure. Ego identifies the most common trait known by most people at the same spot. All people who contributed add 1 to their learning counter with ego.
to-report transmitConformism [ objectType ]
if debug [ print "transmitConformism" ] ; for debugging
let learnedTrait 0
let traitList []
let currentTraits []
ask people-here [
; everyone (including the agent asking) adds their traits to the list (if visibility is taken into consideration, the list has more visible traits
set currentTraits ifelse-value objectType = "point" [
ifelse-value visibility and random-float 100 < 98 [ getMostlyVisibleTraits pointTraits ][ pointTraits ]
][
ifelse-value visibility and random-float 100 < 98 [ getMostlyVisibleTraits potTraits ][ potTraits ]
]
set traitList sentence traitList currentTraits ; adds the traits to the list
]
if not empty? currentTraits [
set learnedTrait one-of modes traitList ; choose the trait known the most (mode) by everyone
; update the interaction links counter
ask other people-here [
updateCampInteractions "learning" self myself ; update the counter of the camps
ask interaction-with myself [ set learning learning + 1 ] ; identify all links from the people here to the person asking
]
]
report learnedTrait ; reports that trait
end
transmitPrestige Parameters {-} objectType can be “point” or “pot”.
Return {-} Returns the trait chosen from the person with higher prestige.
Description {-} People procedure. Ego identifies the person with highest prestige around and takes one of its traits. Ego adds 1 to its learning counter with the person with prestige.
to-report transmitPrestige [ objectType ]
if debug [ print "transmitPrestige" ] ; for debugging
let learnedTrait 0 ; variable for the new trait
let currentTraits []
let teacher 0
ifelse objectType = "point" [
set teacher one-of people-here with-max [ pointPrestige ] ; prestige is updated at each hunt
set currentTraits [ pointTraits ] of teacher
][
; set the teacher based on if potteryPrestige is ON or not
set teacher ifelse-value potteryPrestige [ one-of people-here with-max [ potPrestige ]][ one-of people-here ]
set currentTraits [ potTraits ] of teacher
]
if visibility and random-float 100 < 98 [ set currentTraits getMostlyVisibleTraits currentTraits ] ; if visibility is in play, this reduces the traits to the most visible
if not empty? currentTraits [
; choose amongst the teacher's traits one that has higher fitness based on the current environment (points) or at random (pots)
set learnedTrait ifelse-value objectType = "point" [ getFitTrait ][ one-of currentTraits ]
; update the interaction link counter
if teacher != self [
ask teacher [
updateCampInteractions "learning" self myself ; update the counter of the camps
ask interaction-with myself [ set learning learning + 1 ] ; identify the link from the teacher to the person asking
]
]
]
report learnedTrait
end
totalTraits Parameters {-} trait can be “potTraits” or “pointTraits”.
Return {-} Returns the number of traits currently known by the population.
Description {-} Observer procedure. Called by the Trait counts plot.
to-report totalTraits [ trait ]
let traits length remove-duplicates reduce sentence trait ; calculates the number of DIFFERENT traits known by the population
report traits
end
printMetrics Description {-} Observer procedure. Print multiple metrics in the Command Center.
to printMetrics
print word "assemblageSizePoints: " assemblageSize "point"
print word "assemblageSizePots: " assemblageSize "pot"
print word "elapsedTime: " timer
print word "interactionFreqHunting: " interactionFreq "hunting"
print word "interactionFreqLearning: " interactionFreq "learning"
print word "interactionFreqTrading: " interactionFreq "trading"
print word "interactionFreqVisiting: " interactionFreq "visiting"
print word "brDistBandPoints: " getBRBands "point"
print word "brDistBandPots: " getBRBands "pot"
print word "brDistCampPoints: " getBRCamps "point"
print word "brDistCampPots: " getBRCamps "pot"
print word "meanPeoplePrestige: " precision (mean [ pointPrestige ] of people) 3
print word "modBandPoint: " modularityBand "point"
print word "modBandPot: " modularityBand "pot"
print word "modCampPoint: " modularityCamp "point"
print word "modCampPot: " modularityCamp "pot"
print word "networkMetricPoint: " networkMetric "point"
print word "networkMetricPot: " networkMetric "pot"
print word "PhiPointCamp: " PhiSTCamps "point"
print word "PhiPotCamp: " PhiSTCamps "pot"
print word "PhiPointBand: " PhiSTBands "point"
print word "PhiPotBand: " PhiSTBands "pot"
print word "nTraitsPoints: " nTraits "point"
print word "nTraitsPots: " nTraits "pot"
end
assemblageSize Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the number of objects discarded at camps and bands.
Description {-} Observer procedure. Calculates the number of objects (pot or point) discarded at camps and bands.
to-report assemblageSize [ objectType ] ; objectType can be "pot" or "point"
let assemblageLength 0 ; will take on all the camps and bands' assemblage lengths
ask camps [ ; camps add their assemblage length first
set assemblageLength ifelse-value objectType = "pot" [ assemblageLength + length potAssemblage ][ assemblageLength + length pointAssemblage ]
]
ask bands [ ; then bands
set assemblageLength ifelse-value objectType = "pot" [ assemblageLength + length potAssemblage ][ assemblageLength + length pointAssemblage ]
]
report assemblageLength
end
interactionFreq Parameters {-} activityName is a string describing the activity for which we want the interactions (hunting, learning, trading, or visiting).
Return {-} Returns the rounded sum of all interactions of the given activity.
Description {-} Observer procedure. Calculates the sum of the interactions of the given activity.
to-report interactionFreq [ activityName ]
let sumInteractions (ifelse-value
activityName = "hunting" [ sum [ hunting ] of interactions ] ; take on the sum of all interactions of the given activity, based on the name given
activityName = "learning" [ sum [ learning ] of interactions ]
activityName = "trading" [ sum [ trading ] of interactions ]
activityName = "visiting" [ sum [ visiting ] of interactions ])
report precision sumInteractions 3 ; return that value, rounded to the nearest 3 decimals.
end
getBRBands Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the mean of all inter-band distances (Brainerd-Robinson). Scaled where 1 is most similar, and 0 is no similiarity.
Description {-} Observer procedure. Computes the Brainerd-Robinson similiarity between bands.
to-report getBRBands [ objectType ]
; temporary variable
let allDist []
; create links between bands and calculate distances between each linked band
ask bands [
create-tempLinks-with other bands [ set hidden? false ] ; this makes sure distances are not calculated twice
set allDist lput calculateBR objectType allDist ; call the procedure that actually calculates the distance and add it to the master list
]
set allDist ( reduce sentence allDist ) ; remove all square brackets in the master list (that were added as each band reports 2 distances as one list)
report precision mean allDist 3 ; report the mean of all inter-band distances
end
getBRCamps Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the mean of all inter-CAMP distances (Brainerd-Robinson). Scaled where 1 is most similar, and 0 is no similiarity.
Description {-} Observer procedure. Computes the Brainerd-Robinson similiarity between camps of the same band.
to-report getBRCamps [ objectType ]
; temporary variable
let allDist []
; calculate distances between each camp within the same band and calculate distances between each linked camp
ask camps [
create-tempLinks-with other camps with [ ownedBy = [ ownedBy ] of myself ] [ set hidden? false ] ; this makes sure distances are not calculated twice
set allDist lput calculateBR objectType allDist ; call the procedure that actually calculates the distance and add it to the master list
]
set allDist ( reduce sentence allDist ) ; remove all square brackets in the master list (that were added as each camp reports multiple distances as one list)
report precision mean allDist 3 ; report the mean of all inter-band distances
end
calculateBR Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns a list of distances (Brainerd-Robinson) between the agent asking and its linked agents. Scaled where 1 is most similar, and 0 is no similiarity.
Description {-} Observer procedure. Computes the Brainerd-Robinson similiarity of the objectType between linked agents.
to-report calculateBR [ objectType ]
; populate variables with all distinct types of traits from people in that camp or band
listTraits
; temporary variables used below
let allTraitsA 0
let allTraitsB 0
let allDist []
let dist 0
; calculate Brainerd-Robinson distances between linked bands or camps
set allTraitsA ifelse-value objectType = "point" [ allPoints ][ allPots ] ; take on the traits of the band/camp's people
set allTraitsA getTraitFreqTable allTraitsA ; get the frequency of each trait in the list
ask tempLink-neighbors [ ; asks the linked agents (bands or camps)
set allTraitsB ifelse-value objectType = "point" [ allPoints ][ allPots ] ; take on the traits of the band/camp's people
set allTraitsB getTraitFreqTable allTraitsB ; get the frequency of each trait in the list
set dist getFreqDistList allTraitsA allTraitsB ; get a list of the differences in frequency for each trait in both lists
set dist ifelse-value ( sum dist = 0 )[ 0 ][ 1 - ( sum dist / 200 )] ; compute the Brainerd-Robinson dist with the list provided
set allDist lput dist allDist ; add the dist calculated to the allDist list
]
ask tempLinks [ die ] ; always remove the links after we used them
report allDist
end
listTraits Description {-} Observer procedure. Tells bands and camps to create a list with all the traits known by their people (one list for point and one list for pot). Then collate this list into a table that has the frequency of each trait type.
to listTraits
let traitList []
ask bands [
set traitList [ pointTraits ] of people with [ homeBand = [ who ] of myself ] ; compile all the pointTraits of people from this band
set allPoints reduce sentence traitList ; remove all the square brackets around individual traitlists
set allPoints table:counts allPoints ; count the number of each trait and enter it as a table
set traitList [ potTraits ] of people with [ homeBand = [ who ] of myself ] ; compile all the potTraits of people from this band
set allPots reduce sentence traitList ; remove all the square brackets around individual traitlists
set allPots table:counts allPots ; count the number of each trait and enter it as a table
]
ask camps [ ; comments are the same as above, but for camps
set traitList [ pointTraits ] of people with [ homeCamp = myself ]
set allPoints reduce sentence traitList
set allPoints table:counts allPoints
set traitList [ potTraits ] of people with [ homeCamp = myself ]
set allPots reduce sentence traitList
set allPots table:counts allPots
]
end
getTraitFreqTable Parameters {-} allTraits is a table of all traits collected by either the camp or band. That table has the traits known by their people (key) and the count of those traits (value).
Return {-} Returns a table with the frequency of each trait given as a percentage rather than a count.
Description {-} Camp AND Band procedure. Transforms the count of each trait into a percentage of known traits (frequency).
to-report getTraitFreqTable [ allTraits ]
let sumCountTraits sum table:values allTraits ; the sum count of all traits
foreach table:keys allTraits [ key -> ; for each key in the table
let valCount table:get allTraits key ; get the count associated with the key
table:put allTraits key (( valCount / sumCountTraits ) * 100 ) ; replace the count associated with each key with its frequency (percentage of all known traits)
]
report allTraits
end
getFreqDistList Parameters {-} traitsTableA and traitsTableB are two tables of traits with key = traits, and values = frequency of that trait known by the people of the band/camp asking.
Return {-} Returns a list of the frequency differences for all traits found in the two tables.
Description {-} Band AND Camp procedure. Calculate the difference in popularity between all items in two tables.
to-report getFreqDistList [ traitsTableA traitsTableB ]
; temporary variables
let dist []
let tmpA 0
let tmpB 0
let allKeys remove-duplicates sentence table:keys traitsTableA table:keys traitsTableB ; get all keys in both tables of proportions
foreach allKeys [ key -> ; for each key, this calculates the difference between tables A and B
set tmpA ifelse-value table:has-key? traitsTableA key [ table:get traitsTableA key ][ 0 ] ; take the value associated with the key in table A, if present
set tmpB ifelse-value table:has-key? traitsTableB key [ table:get traitsTableB key ][ 0 ] ; take the value associated with the key in table B, if present
set dist lput abs ( tmpA - tmpB ) dist ; calculate their difference and add to the list
]
report dist
end
modularityBand Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the modularity of bands. Modularity is “defined based on the number of in-community links versus the number of between-community links” (manual).
Description {-} Observer procedure. Calculates the modularity of bands, based on the traits known. Creates links between people based on the similarity of their known traits and uses those links to calculate modularity.
to-report modularityBand [ objectType ] ; objectType can be "pot" or "point"
let modul -999 ; because 0 is an actual possible modularity measure, this is used to troubleshooting (if this returns -999, there is a problem)
; reset the tempLinks, which only keeps the links with 25th percentile Jaccard distance (the most similar)
resetTempLinksForModularity objectType
nw:set-context people tempLinks ; set the context to compute network metrics on
let peoplePerBands []
ask Bands [ set peoplePerBands lput people with [ homeband = [ who ] of myself ] peoplePerBands ] ; group people per band
set modul nw:modularity peoplePerBands ; calculate modularity by considering the number of WITHIN-BAND links versus the number of BETWEEN-BAND links
ask tempLinks [ die ] ; always remove the links after we used them
report modul
end
modularityCamp Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the modularity of camps. Modularity is “defined based on the number of in-community links versus the number of between-community links” (manual).
Description {-} Observer procedure. Calculates the modularity of camps, based on the traits known. Creates links between people based on the similarity of their known traits and uses those links to calculate modularity.
to-report modularityCamp [ objectType ]
let modul -999 ; because 0 is an actual possible modularity measure, this is used to troubleshooting (if this returns -999, there is a problem)
; reset the tempLinks, which only keeps the links with 25th percentile Jaccard distance (the most similar)
resetTempLinksForModularity objectType
nw:set-context people tempLinks ; set the context to compute network metrics on
let peoplePerCamps []
ask camps [ set peoplePerCamps lput people with [ homecamp = myself ] peoplePerCamps ] ; groups people per camp
set modul nw:modularity peoplePerCamps ; calculate modularity by considering the number of WITHIN-CAMP links versus the number of BETWEEN-CAMP links
ask tempLinks [ die ] ; always remove the links after we used them
report modul
end
resetTempLinksForModularity Parameters {-} objectType can be “pot” or “point”.
Description {-} Observer procedure. Calculates the Jaccard distance between traits of people. The most similar create a new set of links used to compute modularity.
to resetTempLinksForModularity [ objectType ] ; objectType can be "pot" or "point"
ask interactions [
setupJaccardDistBetweenTraits objectType
]
let interactionsJaccardList [( list jaccard )] of interactions ; get a list of all interaction Jaccards
let tbl stats:newtable-from-row-list interactionsJaccardList ; transform the list into a stat table, necessary for the next step (quantile)
let threshold stats:quantile tbl 0 25 ; compute the Jaccard threshold (the 25th quantile)
ask interactions with [ jaccard <= threshold ][ ask end1 [ create-tempLink-with [ end2 ] of myself ]] ; create tempLinks between people with jaccard distance below threshold
networkMetric Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns a value between -1 and 1. If it is close to 0, it means the network from artifacts represents relatively well the network from interactions, if it is positive, the network from interactions is more clustered than the network from artifacts, if it is negative, the network from artifacts is more clustered than the network from interactions.
Description {-} Observer reporter. Calculates the difference between the clustering coefficient of the network of people and the clustering-coefficient of the network reconstructed from similar objects.
to-report networkMetric [ objectType ]
; if it is close to 0, it means the network from artifacts represents relatively well the network from interactions
; if it is positive, the network from interactions is more clustered than the network from artifacts
; if it is negative, the network from artifacts is more clustered than the network from interactions
ifelse objectType = "point" [
nw:set-context people ( interactions with [ hunting > 0 ]) ; sets the context through which the metrics will be computed (hunting interactions when looking at points)
][
nw:set-context people ( interactions with [ learning > 0 ]) ; sets the context through which the metrics will be computed (learning interactions when looking at pots)
]
let globalClustering mean [ nw:clustering-coefficient ] of people ; calculate the clustering-coefficient of the whole network
; compute the clustering coefficient between camps on the basis on their shared items
ask campInteractions [
let A ifelse-value objectType = "point" [[ pointAssemblage ] of end1 ] [[ potAssemblage ] of end1 ]
setupJaccardDistBetweenTraits Parameters {-} objectType can be “pot” or “point”.
Description {-} Link procedure. Calculates the Jaccard distance between traits of agents linked by the link asking.
to setupJaccardDistBetweenTraits [ objectType ] ; objectType can be "pot" or "point"
let A ifelse-value objectType = "point" [[ pointTraits ] of end1 ] [[ potTraits ] of end1 ] ; get the correct traits to compare
setupJaccardDistBetweenObjects Parameters {-} objectType can be “pot” or “point”.
Description {-} Link procedure. Calculates the Jaccard distance between assemblages (objects) of agents linked by the link asking.
to setupJaccardDistBetweenObjects [ objectType ] ; objectType can be "pot" or "point"
let A ifelse-value objectType = "point" [[ pointAssemblage ] of end1 ] [[ potAssemblage ] of end1 ]
jaccardDist Parameters {-} A and B are two lists (of traits or objects discarded) to be compared.
Return {-} Returns the Jaccard distance between the two lists.
Description {-} Observer procedure. This calculates the Jaccard distance between two lists (A and B).
to-report jaccardDist [ A B ] ; A and B represent two lists
; The Jaccard distance is 1- Jaccard Index
; The Jaccard index for two sets of numbers is: (the Union of both sets - the Intersection of both sets)/The Union of both sets
let AIB [] ; this will be the intersection of both sets
let AUB 0 ; this will be the union of both sets
foreach A [ x ->
if member? x B [ set AIB lput x AIB ] ; if the item is shared, it is added to the list
]
set AIB length AIB ; the number of shared items
set AUB length ( remove-duplicates ( sentence A B )) ; the total number of items
; calculate the Jaccard index if there are traits already
let JDist ifelse-value AUB = 0 [ 2 ][ ( AUB - AIB ) / AUB ] ; to avoid errors if there are no traits (e.g., at beginning of a simulation)
report JDist
end
jaccardDistWeighted Parameters {-} A and B are two lists (of traits or objects discarded) to be compared.
Return {-} Returns the weighted Jaccard distance between the two lists.
Description {-} Observer procedure. This calculates the weighted Jaccard distance between two lists (A and B).
;; SHOULD NOT BE USED FOR PEOPLE INTERACTIONS AS IT IS VERY SLOW
to-report jaccardDistWeighted [ A B ] ; A and B represent two lists
; The weighted Jaccard index is the sum of the minimum of each set divided by the sum of the max of each set
; 1 - Jaccard index is the similarity value
; faster?
let tblA getTraitFreqTable table:counts A ; get a table of proportions of each trait in list
let tblB getTraitFreqTable table:counts B ; get a table of proportions of each trait in list
let allKeys remove-duplicates ( sentence table:keys tblA table:keys tblB ) ; get all keys in both tables of proportions
; will be used to calculate min and max of each set
let minVals 0
let maxVals 0
foreach allKeys [ key -> ; for each key, this calculates the min and max value between tables A and B
let tmpA ifelse-value table:has-key? tblA key [ table:get tblA key ][ 0 ] ; take the value associated with the key in table A, if present
let tmpB ifelse-value table:has-key? tblB key [ table:get tblB key ][ 0 ] ; take the value associated with the key in table B, if present
set minVals minVals + ( min ( list tmpA tmpB )) ; calculate the min of the two values and add to the previous value (cumulative sum)
set maxVals maxVals + ( max ( list tmpA tmpB )) ; calculate the max of the two values and add to the previous value (cumulative sum)
]
let JDist minVals / maxVals ; calculate the weighted Jaccard distance
report precision JDist 3
end
PhiSTCamps Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the PhiST statistic. If it is close to 0, variation is not significantly structured around camps.
Description {-} Observer reporter. Calculate PhiST statistic between people’s inventories according to the camp they belong to.
to-report PhiSTCamps [ objectType ]
let phi -999 ; for troubleshooting to see errors
carefully [
; first compute jaccard among people between camps
ask people [ ; create tempLinks between people who are NOT from the same camp
let peopleOtherCamps other people with [ homeCamp != [ homeCamp ] of myself ]
create-tempLinks-with peopleOtherCamps
]
; calculate jaccard between traits of those linked people
ask tempLinks [ setupJaccardDistBetweenTraits objectType ]
; calculate variance in jaccard distances
let j [ jaccard ] of tempLinks
let varBetween variance j
ask tempLinks [ die ] ; reset the tempLinks
; now variance within camps
ask people [ ; create tempLinks between people who ARE from the same camp
let peopleMyCamp other people with [ homeCamp = [ homeCamp ] of myself ]
create-tempLinks-with peopleMyCamp
]
; calculate jaccard between traits of those linked people
ask tempLinks [ setupJaccardDistBetweenTraits objectType ]
; calculate variance in jaccard distances
let k [ jaccard ] of tempLinks
let varWithin variance k
ask tempLinks [ die ] ; always remove the links after we used them
; calculate PhiST
set phi ( varBetween - varWithin ) / ( varBetween )
][ show "error in PhiST" ]
report phi
end
PhiSTBands Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the PhiST statistic between bands. If it is close to 0, variation is not significantly structured around camps.
Description {-} Observer reporter. Calculate PhiST statistic between people’s inventories according to the camp they belong to.
to-report PhiSTBands [ objectType ] ; if it is close to 0, variation not significantly structured around camps
let phi -999 ; for troubleshooting to see errors
carefully [
; first compute jaccard among people between bands
ask people [ ; create tempLinks between people who are NOT from the same band
let peopleOtherBands other people with [ homeBand != [ homeBand ] of myself ]
create-tempLinks-with peopleOtherBands
]
; calculate jaccard between traits of those linked people
ask tempLinks [ setupJaccardDistBetweenTraits objectType ]
; calculate variance in jaccard distances
let j [ jaccard ] of tempLinks
let varBetween variance j
ask tempLinks [ die ] ; reset the tempLinks
; now variance within bands
ask people [ ; create tempLinks between people who ARE from the same band
let peopleMyBand other people with [ homeBand = [ homeBand ] of myself ]
create-tempLinks-with peopleMyBand
]
; calculate jaccard between traits of those linked people
ask tempLinks [ setupJaccardDistBetweenTraits objectType ]
; calculate variance in jaccard distances
let k [ jaccard ] of tempLinks
let varWithin variance k
ask tempLinks [ die ] ; always remove the links after we used them
; calculate PhiST
set phi ( varBetween - varWithin ) / ( varBetween )
][ show "error in PhiST" ]
report phi
end
nTraits Parameters {-} objectType can be “pot” or “point”.
Return {-} Returns the number of traits known of the given objectType. It considers repetitions (e.g., if 15 agents all know trait “1” only, the result will be 15).
Description {-} Observer procedure. Calculates the number of traits known by the whole population (with repetitions).
to-report nTraits [ objectType ] ; objectType can be "pot" or "point"
let nList []
ask people [ ; each person adds their known trait to the list
set nList ifelse-value objectType = "pot" [ sentence nList potTraits ][ sentence nList pointTraits ]
]
report length nList
end
compileArtifactSimilarities Description {-} Observer procedure. Calculate the similarity of objects between camps, keep only the most similar and update the relevant adjacency tables with the similarity data.
to compileArtifactSimilarities
; For points
ask campInteractions [
setupJaccardDistBetweenObjects "point" ; compute the similar objects
updateAdjacencyTables campPointAdjacency jaccard ; and that that value to the table
]
; For pots
ask campInteractions [
setupJaccardDistBetweenObjects "pot" ; compute the similar objects
updateAdjacencyTables campPotAdjacency jaccard ; and that that value to the table
]
end
compileCampInteractions Description {-} Observer procedure. Asks the links between camps to update the adjacency tables with their counters.
to compileCampInteractions
ask campInteractions [ ; each link between camps add their counter value to the adjacency table
updateAdjacencyTables campVisitingAdjacency visiting
updateAdjacencyTables campTradingAdjacency trading
updateAdjacencyTables campHuntingAdjacency hunting
updateAdjacencyTables campLearningAdjacency learning
]
end
updateAdjacencyTables Parameters {-} matrixName is the name of the adjacency table getting updated, and addedValue is the value that will be added to the table.
Description {-} Link procedure. Updates the values in the adjacency table given (matrixName) by adding addedValue to the values already there.
to updateAdjacencyTables [ matrixName addedValue ]
let rowN [ who ] of end1 - nOfBands ; get the row number associated with the camp at one end (as camps are created after bands, their numbers follow the band count)
fileNameBase Return {-} Returns a filename with the following structure: “day-month-year-time”.
Description {-} Observer procedure. Takes the current time and produces a name with it.
to-report fileNameBase
let dt (remove ":" remove " " date-and-time)
let t substring dt 0 12
let d substring dt 12 length dt
report (word d "-" t)
end
exportResults Description {-} Observer procedure. Exports variable and values to a CSV file with filename set by fileNameBase.
to exportResults
let fileName fileNameBase ; create a new file name based on date and time
print ( word "file name is " fileName )
let fileNameCSV ( word "results/run_outputAssemblages_" fileName ".csv" ) ; the complete CSV file name
if not file-exists? fileNameCSV [ ; so it does not overwrite it, if it already exists.
; save assemblages
file-open fileNameCSV
file-print ( "band,site,pointAssemblage,potAssemblage,MeanPrestige,MaxPrestige" )
ask camps [ file-print ( word ownedBy "," who "," pointAssemblage "," potAssemblage "," MeanPrestige "," MaxPrestige) ]
ask bands [ file-print ( word "" "," who "," pointAssemblage "," potAssemblage )]
file-close
; save variables to a different file
set fileNameCSV ( word "results/run_outputVariables_" fileName ".csv" )
file-open fileNameCSV
file-print ( "variable,value" )
file-print ( word "nOfBands," nOfBands )
file-print ( word "campPopulation," campPopulation )
file-print ( word "pVisiting," pVisiting )
file-print ( word "pHunting," pHunting )
file-print ( word "huntingGroupSize," huntingGroupSize )
file-print ( word "potBreakageRate," potBreakageRate )
file-print ( word "pointsBreakageRate," pointsBreakageRate )
file-print ( word "pNewObject," pNewObject )
file-print ( word "aggregationFreq," aggregationFreq )
file-print ( word "aggregationDuration," aggregationDuration )
file-print ( word "pMigration," pMigration )
file-print ( word "pNewAlliance," pNewAlliance )
file-print ( word "maximumAlliesN," maximumAlliesN )
file-print ( word "learningMethod," learningMethod )
file-print ( word "transmissionRate," transmissionRate )
file-print ( word "gradualEnvChange," gradualEnvChange )
file-print ( word "pEnvChange," pEnvChange )
file-print ( word "nTicks," nTicks )
file-print ( word "randomSeed," randomSeed )
file-print ( word "uniqueTraits," uniqueTraits )
file-close
; save network links to a different file
set fileNameCSV ( word "results/run_outputLinks_" fileName ".csv" )
file-open fileNameCSV
file-print ( "end1,end2,learning,hunting,visiting,trading" )
exportFineTunedData Description {-} Observer procedure. Writes adjacency tables data to a separate file at 1000 ticks interval if wanted.
to exportFineTunedData
file-open fileHighRes ; open the file already created and write the adjacency tables info with ticks (one per row)
file-print ( word ticks ",campVisiting," matrix:to-row-list campVisitingAdjacency )
file-print ( word ticks ",campTrading," matrix:to-row-list campTradingAdjacency )
file-print ( word ticks ",campHunting," matrix:to-row-list campHuntingAdjacency )
file-print ( word ticks ",campLearning," matrix:to-row-list campLearningAdjacency )
file-print ( word ticks ",campPots," matrix:to-row-list campPotAdjacency )
file-print ( word ticks ",campPoints," matrix:to-row-list campPointAdjacency )
file-close
end
getJaccardBand Parameters {-} objectType is either “points” or “pots”.
Description {-} Observer procedure. Calculates the Jaccard distance between the traits found in the camps belonging to each band. Essentially this is a measure where 0 is most similar and 1 is most distant for sharing the same types and variants between bands.
to-report getJaccardBand [ objectType ]
; temporary variable
let allDist []
; create links between bands and calculate distances between each linked band
ask bands [
create-tempLinks-with other bands [ set hidden? true] ; this makes sure distances are not calculated twice
set allDist lput calculateJaccard objectType allDist ; call the procedure that actually calculates the distance and add it to the master list
]
ifelse not empty? allDist [
set allDist ( reduce sentence allDist ) ; remove all square brackets in the master list (that were added as each band reports 2 distances as one list)
] [
set allDist [0]
]
report precision mean allDist 3 ; report the mean of all inter-band distances
end
getJaccardCamp Parameters {-} objectType is either “points” or “pots”.
Description {-} Observer procedure. This is the same as the getJaccardBand procedure except the values are calculated between camps.
to-report getJaccardCamp [ objectType ]
; temporary variable
let allDist []
; create links between bands and calculate distances between each linked band
ask camps [
create-tempLinks-with other camps [ set hidden? true] ; this makes sure distances are not calculated twice
set allDist lput calculateJaccard objectType allDist ; call the procedure that actually calculates the distance and add it to the master list
]
ifelse not empty? allDist [
set allDist ( reduce sentence allDist ) ; remove all square brackets in the master list (that were added as each band reports 2 distances as one list)
] [
set allDist [0]
]
report precision mean allDist 3 ; report the mean of all inter-band distances
end
calculateJaccard Parameters {-} objectType is either “points” or “pots”.
Description {-} Observer procedure. This function calculates the Jaccard distance between two lists.
to-report calculateJaccard [ objectType ]
; populate variables with all distinct types of traits from people in that camp or band
listTraits
; temporary variables used below
let allTraitsA table:make
let allTraitsB table:make
let allDist []
let dist 0
; calculate Jaccard distances between linked bands or camps
set allTraitsA ifelse-value objectType = "points" [ allPoints ][ allPots ] ; take on the traits of the band/camp's people
set allTraitsA remove-duplicates (table:keys allTraitsA)
ask tempLink-neighbors [ ; asks the linked agents (bands or camps)
set allTraitsB ifelse-value objectType = "points" [ allPoints ][ allPots ] ; take on the traits of the band/camp's people
set allTraitsB remove-duplicates table:keys allTraitsB
set dist jaccardDist allTraitsA allTraitsB
set allDist lput dist allDist ; add the dist calculated to the allDist list
]
ask tempLinks [ die ] ; always remove the links after we used them
report allDist
end
adjCompare Parameters {-} m1 and m2 are matrices.
Return {-} Returns a number between zero and one where one represents identical matrices and zero is no similarity.
Description {-} Observer procedure. This function calculates the absolute differences between two matrices. It is a method to compare networks.
to-report adjCompare [m1 m2]
let result 0
carefully[
let nr item 0 matrix:dimensions m1
let md median reduce sentence reduce sentence matrix:to-row-list m1
set m1 matrix:map [x -> ifelse-value (x > md) [1][0]] m1
set md median reduce sentence reduce sentence matrix:to-row-list m2
set m2 matrix:map [x -> ifelse-value (x > md) [1][0]] m2
let netDiff matrix:minus m1 m2
set netDiff matrix:map abs netDiff
let out 1 - (sum reduce sentence reduce sentence matrix:to-row-list netDiff) / (nr * (nr - 1))
set result precision out 3
][set result "NA"]
report result
end
adjCorr Parameters {-} m1 and m2 are matrices
Return {-} Returns a number between -1 and 1 representing the Pearson’s R correlation.
Description {-} Observer procedure. This procedure uses the sr extension which is disabled by default and this code is commented out. To enable. Uncomment the sr extension, follow the sr setup instructions, install R on your computer with the igraph and multinet packages, and uncomment this code. This procedure takes two matrices and converts them to networks. It then calculates the Pearson’s R correlation between the network edges. It is a method to compare networks.
;to-report adjCorr [m1 m2]
; let result 0
; carefully[
; let md median reduce sentence reduce sentence matrix:to-row-list m1
; set m1 matrix:map [x -> ifelse-value (x > md) [1][0]] m1
; set md median reduce sentence reduce sentence matrix:to-row-list m2
; set m2 matrix:map [x -> ifelse-value (x > md) [1][0]] m2
; sr:set "m1" matrix:to-row-list m1
; sr:set "m2" matrix:to-row-list m2
; (sr:run
; "m1 = as.matrix(do.call(rbind,m1))"
; "m2 = as.matrix(do.call(rbind,m2))"
; "rownames(m1) = colnames(m1) = 1:nrow(m1)"
; "rownames(m2) = colnames(m2) = 1:nrow(m2)"
; "print(m1)"
; "g1 = igraph::graph_from_adjacency_matrix(m1, mode = 'undirected', diag = F)"
; "g2 = igraph::graph_from_adjacency_matrix(m2, mode = 'undirected', diag = F)"
; "net = multinet::ml_empty('comparison')"
; "multinet::add_igraph_layer_ml(n = net, g = g1,name = 'm1')"
; "multinet::add_igraph_layer_ml(n = net, g = g2,name = 'm2')"
; "comp = multinet::layer_comparison_ml(net, method = 'pearson.degree', mode = 'all')"
; )
; set result precision sr:runresult "comp[1,2]" 3
; ][set result "NA"]
; report result
;end
GUI elements | Name | Type | Value | Properties | | —- | —- | —– | ———- | | nTicks | INPUTBOX | 2000.0 | entrytype: Number | | randomSeed | INPUTBOX | 9433.0 | entrytype: Number | | randomizeSeed | SWITCH | TRUE | | | saveResults | SWITCH | TRUE | | | pLoss | SLIDER | 33 | min: 0; max: 100; incr: 1 | | aggregationFreq | SLIDER | 12 | min: 6; max: 36; incr: 1 | | pNewObject | SLIDER | 10 | min: 0.05; max: 100; incr: 0.05 | | aggregationDuration | SLIDER | 1 | min: 0; max: 12; incr: 1 | | pMigration | SLIDER | 0.05 | min: 0; max: 10; incr: 0.1 | | transmissionRate | SLIDER | 5 | min: 0; max: 100; incr: 1 | | nOfBands | SLIDER | 3 | min: 3; max: 6; incr: 1 | | campPopulation | SLIDER | 15 | min: 1; max: 100; incr: 1 | | pHunting | SLIDER | 1 | min: 0; max: 100; incr: 1 | | pVisiting | SLIDER | 5 | min: 0; max: 100; incr: 1 | | huntingGroupSize | SLIDER | 5 | min: 1; max: 100; incr: 1 | | pNewAlliance | SLIDER | 0.5 | min: 0; max: 10; incr: 0.05 | | maximumAlliesN | SLIDER | 50 | min: 10; max: 150; incr: 1 | | pLuckyLeap | SLIDER | 1 | min: 0; max: 10; incr: 0.1 | | learningMethod | SLIDER | 50 | min: 0; max: 100; incr: 1 | | recordPrestige | SWITCH | TRUE | | | gradualEnvChange | SWITCH | TRUE | | | pEnvChange | SLIDER | 0.1 | min: 0; max: 20; incr: 0.1 | | uniqueTraits | SWITCH | TRUE | | | visibility | SWITCH | TRUE | | | genderedMigration | SWITCH | TRUE | | | pNewVariant | SLIDER | 6 | min: 0; max: 100; incr: 1 | | potteryPrestige | SWITCH | TRUE | | | pTrading | SLIDER | 5 | min: 0; max: 100; incr: 1 | | recordPerson | SWITCH | TRUE | | | debug | SWITCH | TRUE | | | showMetrics | SWITCH | TRUE | | | HighResolutionData | SWITCH | TRUE | |
References cited
Creanza, N., Kolodny, O., Feldman, M.W., 2017. Greater than the sum of its parts? Modelling population contact and interaction of cultural repertoires. Journal of The Royal Society Interface. 14, 20170171.
Derex, M., Mesoudi, A., 2020. Cumulative Cultural Evolution within Evolving Population Structures. Trends in Cognitive Sciences. 24, 654–667.
Hamilton, M.J., 2021. Collective computation and the emergence of hunter-gatherer small-worlds (preprint). Ecology.
Kolodny, O., Creanza, N., Feldman, M.W., 2015. Evolution in leaps: The punctuated accumulation and loss of cultural innovations. Proceedings of the National Academy of Sciences. 112, E6762–E6769.
Mesoudi, A., 2011. Cultural Evolution: How Darwinian Theory Can Explain Human Culture and Synthesize the Social Sciences. University of Chicago Press.
Migliano, A.B., Battiston, F., Viguier, S., Page, A.E., Dyble, M., Schlaepfer, R., Smith, D., Astete, L., Ngales, M., Gomez-Gardenes, J., Latora, V., Vinicius, L., 2020. Hunter-gatherer multilevel sociality accelerates cumulative cultural evolution. Science Advances. 6, eaax5913.
Powell, A., Shennan, S., Thomas, M.G., 2009. Late Pleistocene Demography and the Appearance of Modern Human Behavior. Science. 324, 1298–1301.
Premo, L.S., Tostevin, G.B., 2016. Cultural Transmission on the Taskscape: Exploring the Effects of Taskscape Visibility on Cultural Diversity. PLOS ONE. 11, e0161766.
Rorabaugh, A., 2015. Modeling Pre-European Contact Coast Salish Seasonal Social Networks and Their Impacts on Unbiased Cultural Transmission. Journal of Artificial Societies and Social Simulation. 18.
Salali, G.D., Chaudhary, N., Bouer, J., Thompson, J., Vinicius, L., Migliano, A.B., 2019. Development of social learning and play in BaYaka hunter-gatherers of Congo. Scientific Reports. 9.
Salali, G.D., Chaudhary, N., Thompson, J., Grace, O.M., van der Burgt, X.M., Dyble, M., Page, A.E., Smith, D., Lewis, J., Mace, R., Vinicius, L., Migliano, A.B., 2016. Knowledge-Sharing Networks in Hunter-Gatherers and the Evolution of Cumulative Culture. Current Biology. 26, 2516–2521.
Yasuoka, H., 2006. Long-Term Foraging Expeditions (Molongo) among the Baka Hunter-Gatherers in the Northwestern Congo Basin, with Special Reference to the “Wild Yam Question.” Human Ecology. 34, 275–296.
authors listed alphabetically by first name–each contributed equally to the agent based model.↩︎