Skip to content

Instantly share code, notes, and snippets.

@unchitta

unchitta/DW1.jl

Last active Oct 3, 2020
Embed
What would you like to do?
mutable struct Agent
id::Int
current_opinion::Float64
previous_opinion::Float64
new_opinion::Float64
neighbors::Array{Int64}
end
N = 100; p=0.05; c=0.3; alpha=0.1; q=1
model = initialize_model(N,p,c,alpha,q)
data = run!(model,100000)
steps = size(data)[1]
using Plots
plot(1:steps, data, label=false)
struct Model
agents::Dict{Int,Agent}
c::Float64
alpha::Float64
q::Int
end
function Model(c::Float64, alpha::Float64, q::Int)
agents = Dict{Int, Agent}()
return Model(agents,c,alpha,q)
end
# implement function to get and set an agent in the model by index
# we can do this by extending Base.getindex and Base.setindex!
function Base.getindex(model::Model, id::Int)
model.agents[id]
# it is implied in julia that this is the thing that will be returned
# without us explicitly having to use the keyword `return`
end
#=
Note: a more concise way to write the function above is
function Base.getindex(model::Model, agent::Agent) = model.agents[id]
=#
function Base.setindex!(model::Model, agent::Agent, id::Int)
# set the value of the agents dict in the model to be this given agent
model.agents[id] = agent
# TODO: need to implement a check to make sure that the id is valid
end
allagents(model::Model) = values(model.agents) # this returns the values in the dict containing ids as keys and agents as values
function add_agent!(model::Model, agent::Agent)
model[agent.id] = agent
return nothing
end
# a function to put agents in the model in a G(N,p) network
function initialize_network!(model::Model, N::Int, p::Float64)
# compute the number of edges according to the ER model:
# number of possible edges in the graph with N nodes
M_max = div(N*(N-1), 2)
# number of edges in a random G(N,p) network (binomially distributed)
M = rand(Binomial(M_max,p))
m=0
while m<M
edge = Tuple(rand(1:N, 2))
edge[1] != edge[2] || continue # if we have a self-loop, generate a new one
a1 = edge[1]; a2 = edge[2]
@inbounds n1 = model[a1].neighbors
# index gives us the index of the first value in n1 greater than or equal to x;
# and length(n1)+1 if a2 is greater than all values in n1
index = searchsortedfirst(n1, a2)
# check if agent 2 is in the list of agent 1's neighbors
# (continue to generate a new edge if edge is already in graph)
@inbounds (index <= length(n1) && n1[index] == a2) && continue
insert!(n1, index, a2)
# using insert! ensures our lists of neighbors will remain sorted and searchsortedfirst can work efficiently
# do the same for agent 2
@inbounds n2 = model[a2].neighbors
index = searchsortedfirst(n2, a1)
insert!(n2, index, a1)
m+=1
end
end
# a function to initialize model with N agents with uniformly distributed initial opinion
function initialize_model(N::Int, p::Float64, c::Float64, alpha::Float64, q::Int)
model = Model(Agent,c,alpha,q)
for i in 1:N
x = rand() # Uniform([0,1]) RV
add_agent!(model, Agent(i,x,x,x,Int[]))
end
initialize_network!(model,N,p)
return model
end
function update_opinions!(model::Model, N::Int, q::Int)
# TODO: check if q <= N
# choose q random nodes
nodes_to_update = rand(1:N, q)
for i in nodes_to_update
agent1 = model[i]
agent2 = model[rand(agent1.neighbors)] # choose a random neighbor
if abs(agent1.current_opinion - agent2.current_opinion) < model.c
agent1.new_opinion = agent1.new_opinion + model.alpha*(agent2.current_opinion - agent1.current_opinion)
agent2.new_opinion = agent2.new_opinion + model.alpha*(agent1.current_opinion - agent2.current_opinion)
end
end
end
# this function will help us determine when to stop our trial
# Reference: Agents.jl (Hegselmann-Kraus model)
function opinions_have_changed(model::Model)
if any(
!isapprox(a.previous_opinion, a.new_opinion; rtol = 1e-12) for a in allagents(model)
)
return false
else
return true
end
end
function run!(model::Model, max_steps::Int)
N = length(allagents(model))
q = model.q
# initialize a [max_steps, N]-array with UnionType{missing,Float64}
data = Array{Union{Missing, Float64}}(missing, max_steps, N)
# run until max_steps or termination criterion is met
tstep = 1
stationary = 0
while (tstep <= max_steps) & (stationary < 100)
update_opinions!(model,N,q)
for agent in allagents(model)
agent.previous_opinion = agent.current_opinion
agent.current_opinion = agent.new_opinion
# collect data here
data[tstep,agent.id] = agent.current_opinion
end
tstep += 1
opinions_have_changed(model) ? stationary+=1 : stationary=0
end
return data[1:tstep-1,:]
end
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment