Skip to content

Commit

Permalink
weighed-average vs pick-or-average strategies
Browse files Browse the repository at this point in the history
  • Loading branch information
mjaquiery committed Mar 4, 2019
1 parent fcf13ad commit f9be620
Show file tree
Hide file tree
Showing 5 changed files with 949 additions and 67 deletions.
23 changes: 23 additions & 0 deletions diary.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,3 +49,26 @@ Decisionmaking in humans seems to be more driven by a bimodal choose-self or ave

The model should be run on the more basic manipulations to begin with.

## 25/02/2019

Using the choose-self-or-average strategy under the basic manipulations of advice noise and bad advice, the agents show an extreme sensitivity to the parameter values. Whether they remain stable at any parameter value with an outcome between 0 and 1, or whether only 0 and 1 are stable points at all parameter values. There's potentially a saddle effect where parameter values below some critical number tend eventually to stability at 0, and parameter values above it tend eventually to stability at 1.

Note that 0 and 1 refer to probabilities of choosing self rather than averaging, not to the resultant mean advice-weight. Because simple averaging always occurs where the probability is 0, the own-advice-weight ranges from .5 to 1 for these models.

![Noisy advice](results/2019-02-21_23-50-17_d1-a2_graph.png "Note the remarkable sensitivity to parameter values")

![Bad advice](results/2019-02-23_15-15-17_Uncapped_bad-advice_graph.png "An even closer examination of the potential saddle point.")

## 01/02/2019

Futher exploration in model 19 demonstrates that the relationship is not as clear as previously thought. Over 3000 generations it is plausible that the stability-at-extremes effect emerges, but simulation over more generations is required.

![Weighted-average strategy](results/2019-03-01_01-57-22_Weighted-Average-with bad advice_graph.png "The existence of a fairly stable mid-point solution seems to be detectable at .98 weighting")

![Pick-or-average strategy](results/2019-03-01_05-55-08_Pick-or-Average-with bad advice_graph.png "Plausibly even the most midrange value (.98) is tending towards a stable extreme value.")

For computational reasons this needs to be broken up into multiple rounds. Each generation acts as a save point for the models, so they can be resumed from previous saves without compromise. The core package needs to be updated to handle taking an existing dataframe as input.

## 04/03/2019

A **critical** mistake in Model 20 means previous result all used the pick-or-average strategy. The models must be rerun now that Model 20 has been fixed to use the correct weighted-average function where applicable.
88 changes: 59 additions & 29 deletions evoSim/evoSim/R/evoSim.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -303,14 +303,14 @@ selectParents <- function(modelParams, agents, world, ties) {
# the others get weighted by relative fitness which are transformed to +ve values
tmp$fitness <- tmp$fitness - min(tmp$fitness) + 1
# scale appropriately
while(any(tmp$fitness < 10))
while (any(tmp$fitness < 10))
tmp$fitness <- tmp$fitness * 10
# and round off
tmp$fitness <- round(tmp$fitness)
tickets <- vector(length = sum(tmp$fitness)) # each success buys a ticket in the draw
i <- 0
for(a in 1:nrow(tmp)) {
tickets[(i+1):(i+1+tmp$fitness[a])] <- a
for (a in 1:nrow(tmp)) {
tickets[(i + 1):(i + 1 + tmp$fitness[a])] <- a
i <- i + 1 + tmp$fitness[a]
}
winners <- sample(tickets, modelParams$agentCount, replace = T)
Expand All @@ -326,6 +326,8 @@ selectParents <- function(modelParams, agents, world, ties) {
#' @param decisionCount number of decisions made each generation
#' @param generationCount number of generations (duration of model)
#' @param mutationChance probability of a mutation occuring
#' @param genZero data frame of agents to be used as the starting generation. If
#' this is NULL the first generation is created using the makeAgentsFun
#' @param learnRate scaling factor for updating connection weights
#' @param other list used for passing other arguments to custom functions.
#' Accessed by \code{modelParams$other}
Expand Down Expand Up @@ -376,6 +378,7 @@ evoSim <- function(agentCount,
decisionCount,
generationCount,
mutationChance,
genZero = NULL,
learnRate = 0.00,
other = NULL,
makeAgentsFun = NULL,
Expand Down Expand Up @@ -413,6 +416,7 @@ evoSim <- function(agentCount,
recordDecisions = recordDecisions)

n <- modelParams$agentCount * modelParams$generationCount

agents <- data.frame(id = 1:n,
genId = rep(1:modelParams$agentCount, modelParams$generationCount),
generation = rep(1:modelParams$generationCount, each = modelParams$agentCount),
Expand All @@ -423,8 +427,39 @@ evoSim <- function(agentCount,
advice = rep(NA, n),
finalDecision = rep(NA, n))

tStart <- Sys.time()

# overwrite generation by spawning from genZero if necessary
# N.B. the ties are regenerated, so this will need rewriting if ties are
# important to remember over cycle breaks
if (!is.null(genZero)) {
lastGen <- max(genZero$generation)
parents <- modelParams$selectParents(modelParams, genZero,
list(generation = lastGen,
decision = 0),
NULL)

tmp <- modelParams$makeAgents(modelParams, genZero, parents = parents)

# Update the generation values for compatability
agents$generation <- rep((lastGen + 1):(lastGen + modelParams$generationCount),
each = modelParams$agentCount)
} else {
tmp <- modelParams$makeAgents(modelParams)
}

tmp$agents$fitness <- 0

# ensure all column names appear in both hardcoded and user-supplied data frames
agents[, names(tmp$agents)[which(!(names(tmp$agents) %in% names(agents)))]] <- NA

# fill in the initial generation data in the correct columns
agents[agents$generation == min(agents$generation), names(agents) %in% names(tmp$agents)] <-
tmp$agents[,names(tmp$agents)[order(match(names(tmp$agents), names(agents)))]]
ties <- tmp$ties

n <- n * modelParams$decisionCount
if(recordDecisions)
if (recordDecisions)
decisions <- data.frame(id = rep(agents$id, each = modelParams$decisionCount),
genId = rep(agents$genId, each = modelParams$decisionCount),
generation = rep(agents$generation, each = modelParams$decisionCount),
Expand All @@ -436,22 +471,13 @@ evoSim <- function(agentCount,
advisor = rep(NA, n),
advice = rep(NA, n),
finalDecision = rep(NA, n))

tStart <- Sys.time()
tmp <- modelParams$makeAgents(modelParams)
tmp$agents$fitness <- 0
# ensure all column names appear in both hardcoded and user-supplied data frames
agents[, names(tmp$agents)[which(!(names(tmp$agents) %in% names(agents)))]] <- NA

# fill in the initial generation data in the correct columns
agents[agents$generation==1, names(agents) %in% names(tmp$agents)] <-
tmp$agents[,names(tmp$agents)[order(match(names(tmp$agents), names(agents)))]]
ties <- tmp$ties

for(g in 1:modelParams$generationCount) {

genStart <- min(agents$generation)

for (g in genStart:max(agents$generation)) {
# Make decisions
mask <- which(agents$generation==g)
for(d in 1:modelParams$decisionCount) {
mask <- which(agents$generation == g)
for (d in 1:modelParams$decisionCount) {
# Correct answer (same for all agents)
world <- list(generation = g, decision = d)
world$state <- modelParams$getWorldState(modelParams, world)
Expand All @@ -466,33 +492,37 @@ evoSim <- function(agentCount,
# Update connections
ties <- modelParams$updateConnections(modelParams, agents, world, ties)
# Record decision
if(recordDecisions)
if (recordDecisions & (g %% recordDecisions == 0))
decisions[decisions$generation == g & decisions$decision == d,
c('id', 'genId', 'generation', 'fitness', 'egoBias',
'initialDecision', 'advisor', 'advice', 'finalDecision')] <-
agents[mask, c('id', 'genId', 'generation', 'fitness', 'egoBias',
'initialDecision', 'advisor', 'advice', 'finalDecision')]
}
if(g==modelParams$generationCount)
if (g == modelParams$generationCount)
next()
# Evolve
parents <- modelParams$selectParents(modelParams, agents, world, ties)
tmp <- modelParams$makeAgents(modelParams, agents[agents$generation == g, ], parents)
# combine new generation into full data frame while allowing for mismatched column orders
agents[agents$generation==g+1, names(agents) %in% names(tmp$agents)] <-
agents[agents$generation == g + 1, names(agents) %in% names(tmp$agents)] <-
tmp$agents[,names(tmp$agents)[order(match(names(tmp$agents), names(agents)))]]
ties <- tmp$ties
}

tEnd <- Sys.time()
tElapsed <- as.numeric(format(tEnd,"%s")) - as.numeric(format(tStart,"%s"))
# Return an output containing inputs and the agents data frame
modelData <- list(model = modelParams,
# to avoid duplicating info we reduce agents dataframe to remove columns about decisions
agents = agents[,-which(names(agents) %in%
c('initialDecision', 'advisor', 'advice', 'finalDecsision'))],
duration = tElapsed)
if(recordDecisions)
modelData <- c(modelData, list(decisions = decisions))

if (recordDecisions) {
modelData <- list(model = modelParams,
# to avoid duplicating info we reduce agents dataframe to remove columns about decisions
agents = agents[,-which(names(agents) %in%
c('initialDecision', 'advisor', 'advice', 'finalDecsision'))],
decisions = decisions,
duration = tElapsed)
} else {
modelData <- list(model = modelParams, agents = agents, duration = tElapsed)
}
return(modelData)
}
89 changes: 51 additions & 38 deletions model18.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -47,8 +47,9 @@ runModel <- function(spec) {
data <- evoSim(agentCount = spec$agents,
agentDegree = spec$degree,
decisionCount = spec$decisions,
generationCount = 3000,
generationCount = 1000,
mutationChance = 0.01,
# recordDecisions = 100,
other = list(sensitivity = spec$sensitivity,
sensitivitySD = spec$sensitivitySD,
startingEgoBias = spec$startingEgoBias,
Expand Down Expand Up @@ -79,13 +80,10 @@ runModel <- function(spec) {
adviceNoise = rep(NA, modelParams$agentCount),
confidenceScaling = sample(1:5, 1))

# Overwrite the agents' egobias and sensitivity by
# Overwrite the agents' egobias by
# inheritance from parents where applicable
if (!is.null(parents)) {
if (modelParams$other$manipulation)
evolveCols <- c('egoBias', 'sensitivity')
else
evolveCols <- c('egoBias')
evolveCols <- c('egoBias')
makeAgents.agents[ , evolveCols] <- previousGeneration[parents, evolveCols]
# mutants inherit from a normal distribution
for (x in evolveCols) {
Expand Down Expand Up @@ -138,7 +136,9 @@ runModel <- function(spec) {

# Decision functions - uncapped continuous (default), capped continuous, and discrete ####
uncappedDecisionFun <- function(modelParams, agents, world, ties, initial = F) {
adviceNoise <- ifelse(modelParams$other$manipulation, modelParams$other$adviceNoise, 0)
adviceNoise <- ifelse(modelParams$other$manipulation > 0,
modelParams$other$adviceNoise[modelParams$other$manipulation],
0)
mask <- which(agents$generation == world$generation)
if (initial) {
# initial decision - look and see
Expand Down Expand Up @@ -166,7 +166,9 @@ uncappedDecisionFun <- function(modelParams, agents, world, ties, initial = F) {
}

cappedDecisionFun <- function(modelParams, agents, world, ties, initial = F) {
adviceNoise <- ifelse(modelParams$other$manipulation, modelParams$other$adviceNoise, 0)
adviceNoise <- ifelse(modelParams$other$manipulation > 0,
modelParams$other$adviceNoise[modelParams$other$manipulation],
0)
mask <- which(agents$generation == world$generation)
if (initial) {
# initial decision - look and see
Expand Down Expand Up @@ -245,7 +247,9 @@ categoricalFitnessFun <- function(modelParams, agents, world, ties) {

# Advice functions - there's noisy advice and bad advice ####
noisyAdviceFun <- function(modelParams, agents, world, ties) {
adviceNoise <- ifelse(modelParams$other$manipulation, modelParams$other$adviceNoise, 0)
adviceNoise <- ifelse(modelParams$other$manipulation > 0,
modelParams$other$adviceNoise[modelParams$other$manipulation],
0)
mask <- which(agents$generation == world$generation)
agents$advisor[mask] <- apply(ties, 1, function(x) sample(which(x != 0),1))
# Fetch advice as a vector
Expand All @@ -260,7 +264,9 @@ noisyAdviceFun <- function(modelParams, agents, world, ties) {
}

badAdviceFun <- function(modelParams, agents, world, ties) {
badAdviceProb <- ifelse(modelParams$other$manipulation, modelParams$other$adviceNoise/100, 0)
badAdviceProb <- ifelse(modelParams$other$manipulation > 0,
modelParams$other$adviceNoise[modelParams$other$manipulation]/100,
0)
mask <- which(agents$generation == world$generation)
agents$advisor[mask] <- apply(ties, 1, function(x) sample(which(x != 0),1))
# Fetch advice as a vector
Expand Down Expand Up @@ -291,26 +297,26 @@ confAdviceFun <- function(modelParams, agents, world, ties) {
return(agents)
}

for (decisionType in 3) {
for (adviceType in 4) {
for (decisionType in 1) {
for (adviceType in 2) {
# Storage path for results
resultsPath <- ifelse(ARC,'results/','results/')
time <- format(Sys.time(), "%F_%H-%M-%S")
resultsPath <- paste0(resultsPath,'d',decisionType,'a',adviceType,'_',time)

# Clear the result storage variables
suppressWarnings(rm('rawdata'))

# Parameter space to explore
specs <- list()
for (s in c(1, 10))
for (x in c(F, T))
specs[[length(specs) + 1]] <- list(agents = 1000, degree = 10,decisions = 30,
sensitivity = s,sensitivitySD = s,
startingEgoBias = .45, # .45 is neither optimal nor extreme
adviceNoise = 5,
manipulation = x,
wd = getwd())
for (x in 1:5)
specs[[length(specs) + 1]] <- list(agents = 1000, degree = 10,
decisions = 30,
sensitivity = s, sensitivitySD = s,
startingEgoBias = .45, # .45 is neither optimal nor extreme
adviceNoise = seq(.9, 1, .02),
manipulation = x,
wd = getwd())

if (decisionType == 1) {
for (i in 1:length(specs)) {
Expand Down Expand Up @@ -370,6 +376,9 @@ for (decisionType in 3) {

# Run parallel repetitions of the model with these settings
startTime <- Sys.time()

print(paste('Processing models:', specs[[1]]$shortDesc))

if (!ARC) {
degreeResults <- lapply(specs, runModel)
} else {
Expand All @@ -387,12 +396,14 @@ for (decisionType in 3) {
print(paste0('...complete.'))
print(Sys.time() - startTime)

print('Estimated data size:')
print(object.size(rawdata), units = 'auto')
print(paste0('Estimated data size: ', object.size(rawdata) * 1e-6, 'MB'))
print('Saving data...')
# Save data
save(rawdata, file = paste(resultsPath, 'rawdata.Rdata'))
print('...saved rawdata...')

# Save data with a nice name
resultsPath <- paste0(resultsPath, time, '_',
sub(" ", "-",
sub(" decisions with ", "_", specs[[1]]$shortDesc)))
save(rawdata, file = paste0(resultsPath, '_rawdata.Rdata'))
# Smaller datafile for stopping me running out of memory during analysis
allAgents <- NULL
#allDecisions <- NULL
Expand All @@ -405,42 +416,44 @@ for (decisionType in 3) {
rd$agents$sdSensitivity <- rep(rd$model$other$sensitivitySD,nrow(rd$agents))
rd$agents$startingEgoBias <- rep(rd$model$other$startingEgoBias,nrow(rd$agents))
rd$agents$manipulation <- rep(rd$model$other$manipulation,nrow(rd$agents))
rd$agents$manipulationValue <- rd$model$other$adviceNoise[rd$agents$manipulation]
rd$agents$description <- rep(rd$model$other$shortDesc, nrow(rd$agents))
# only take a subset because of memory limitations
# allAgents <- rbind(allAgents, rd$agents)
allAgents <- rbind(allAgents, rd$agents[rd$agents$generation %% 50 == 1
| (rd$agents$generation %% 25 == 1 & rd$agents$generation < 250), ])
#allDecisions <- rbind(allDecisions, rd$decisions[rd$decisions$generation %in% allAgents$generation, ])
}
#toSave <- list(allAgents, allDecisions)
save(allAgents, file = paste(resultsPath, 'rawdata_subset.Rdata'))
#rm('toSave')
print('...saved subset...')

save(allAgents, file = paste0(resultsPath, '_rawdata-subset.Rdata'))

# Plot
ggplot(allAgents,
aes(x = generation, y = egoBias,
colour = manipulation)) +
colour = factor(manipulationValue))) +
geom_hline(yintercept = 0.5, linetype = 'dashed') +
stat_summary(geom = 'point', fun.y = mean, size = 3, alpha = 0.25) +
stat_summary(fun.data = mean_cl_boot, fun.args = (conf.int = .99), geom = 'errorbar', size = 1) +
scale_y_continuous(limits = c(0,1)) +
facet_grid(~meanSensitivity, labeller = label_both) +
labs(title = allAgents$description[1]) +
style
ggsave(paste(resultsPath, 'graph.png'))
ggsave(paste0(resultsPath, '_graph.png'))

print('...saved graph...')
print('...data saved.')
}
}

ggplot(allAgents[allAgents$meanSensitivity == 1 & allAgents$generation %% 1000 == 1, ],
aes(x = egoBias, y = sensitivity, colour = fitness)) +
geom_point() +
facet_wrap(~generation, labeller = label_both) +
# stat_summary(geom = 'point', fun.y = mean, size = 3) +
style
# For all of the agents, look at how the egocentric bias affects the fitness.
# If there really is a step function, we might expect a non-normal distribution,
# and that distribution's peak should shift as the parameter value changes.
allAgents$egoBiasBin <- cut(allAgents$egoBias, 10)
ggplot(allAgents, aes(x = egoBiasBin, y = fitness, colour = factor(meanSensitivity))) +
stat_summary(geom = 'point', fun.y = mean) +
stat_summary(geom = 'errorbar', fun.data = mean_cl_normal) +
labs(title = "Fitness by egoBias") +
facet_wrap(~manipulationValue, labeller = label_both) +
style

# Cleanup
if (ARC)
Expand Down
Loading

0 comments on commit f9be620

Please sign in to comment.