PK «©ML•ˆp
codemeta.json{"@context": ["https://doi.org/doi:10.5063/schema/codemeta-2.0", "http://schema.org"], "@type": "SoftwareSourceCode", "provider": {"url": "https://www.comses.net", "@type": "Organization", "name": "CoMSES Network (CoMSES)", "@id": "https://www.comses.net"}}PK «©ML;ˆ(î î docs/documentation.txtKnowledge Space Opinion Dynamics Model
Shane T. Mueller
Model appearing in Mueller & Tan (under review) Cognitive perspectives on opinion dynamics: The role of knowledge in consensus formation, opinion divergence, and group polarization
The model is implemented in the R statistical computing language.
It uses the sna library for visualizing knowledge spaces
Simulation uses global-scoped variables for control.
Model functions:
setupKnowledgeLattice(numstates,numfeatures,type)
Sets up an entire lattice of states
getp(type)
When creating a node, gets the 'p' probability of its states
precedes(l1,l2)
Determines whether one node is a subset of another
editprec(matrix)
Increments precedence matrix
newpos(pos,adj,scale,byrow)
Figures out how to position a the graph left-to-right to have things centered
setup(knowledge,n)
Initializes knowledge state
run(n,keep)
runs a specified number of iterations. The given n indicates the number of 10x reps that will
be run, so 1000 is actually 10K interactions.
update(kl)
Updates knowledge space
isvalid(belief,states)
Determines whether a feature-state is a coherent state according to the graph.
Global Control variables
numstates (number of states to use)
numfeatures (number of features representing each state)
numagents (number of agents in the simulation)
mu (probability that a feature will be selected for exchange)
Basic operation:
A minimal example:
numstates <<- 15
numfeatures <<- 20
numagents <<- 10
mu <<- .3 ##convergence parameter
##set up new knowledge space:
kl <<- setupKnowledgeLattice(numstates,numfeatures,type="spread")
setup(kl,numagents) ##create agents within the knowledge space
dat <- run(100,F) ##run for 100x10 interactions.
PK «©ML;ˆ(î î code/documentation.txtKnowledge Space Opinion Dynamics Model
Shane T. Mueller
Model appearing in Mueller & Tan (under review) Cognitive perspectives on opinion dynamics: The role of knowledge in consensus formation, opinion divergence, and group polarization
The model is implemented in the R statistical computing language.
It uses the sna library for visualizing knowledge spaces
Simulation uses global-scoped variables for control.
Model functions:
setupKnowledgeLattice(numstates,numfeatures,type)
Sets up an entire lattice of states
getp(type)
When creating a node, gets the 'p' probability of its states
precedes(l1,l2)
Determines whether one node is a subset of another
editprec(matrix)
Increments precedence matrix
newpos(pos,adj,scale,byrow)
Figures out how to position a the graph left-to-right to have things centered
setup(knowledge,n)
Initializes knowledge state
run(n,keep)
runs a specified number of iterations. The given n indicates the number of 10x reps that will
be run, so 1000 is actually 10K interactions.
update(kl)
Updates knowledge space
isvalid(belief,states)
Determines whether a feature-state is a coherent state according to the graph.
Global Control variables
numstates (number of states to use)
numfeatures (number of features representing each state)
numagents (number of agents in the simulation)
mu (probability that a feature will be selected for exchange)
Basic operation:
A minimal example:
numstates <<- 15
numfeatures <<- 20
numagents <<- 10
mu <<- .3 ##convergence parameter
##set up new knowledge space:
kl <<- setupKnowledgeLattice(numstates,numfeatures,type="spread")
setup(kl,numagents) ##create agents within the knowledge space
dat <- run(100,F) ##run for 100x10 interactions.
PK «©ML–HJÿ` ` ! code/MuellerKnowledgeSpaceModel.R## Knowledge Spaces in Opinion Dynamics
## (c) 2017 Shane T. Mueller (shanem@mtu.edu)
## Software code in the R statistical computing language
## for simulating opinion dynamics using a knowledge-space
## model.
##
##use sna library for visualization
library(sna)
setupKnowledgeLattice <- function (numstates,numfeatures,type="spread")
{
##generate feature states along a spectrum.
kl <- matrix(0,ncol=numfeatures,nrow=(numstates-2))
for(i in 1:(numstates-2))
{
p <- getp(type)
kl[i,] <- runif(numfeatures)
2] <- 0
mat
}
newpos <- function(pos,adj,scale=10,byrow=T)
{
##make it the mean of nodes connected to it.
newpos <- pos
weights <- rep(0,nrow(adj))
for(i in 1:nrow(adj))
{
##get a list of all nodes that connect TO i
to <- adj[,i]
from <- adj[i,]
tofrom <- to|from
##adjust x val to be the average of to and from positions.
print(tofrom)
print(sum(tofrom))
if(sum(tofrom)>1)
{
newpos[i,1] <- mean(pos[tofrom,][,1])
}
weights[i] <- sum(tofrom)+1
}
##
if(byrow)
{
x <- data.frame(newpos,weights=weights)
np.min <- aggregate(newpos[,1],list(newpos[,2]),min)$x
np.max <- aggregate(newpos[,1],list(newpos[,2]),max)$x
# np.mean <- aggregate(newpos[,1],list(newpos[,2]),mean)$x
np.mean <- by(x,x$X2,function(g){weighted.mean(g$X1,g$weights)})
ids <- 1:nrow(adj)
##
np.range <- np.max-np.min
np.range[np.range==0] <- 1 ##rescore 0 range to be 1;
##we divide by this and it shouldn't be undefined.
##this scales for each row
ii <- 1
for(i in levels(as.factor(pos[,2])))
{
tmp <- ids[as.character(pos[,2])==i]
num <- length(tmp)
##scale to fit in -scale..scale
range <- np.range[ii]
newpos[,1][tmp] <- (newpos[,1][tmp] - np.mean[ii])/range * scale+np.mean[ii]
ii <- ii+1
}
}else {
if(1){
min <- min(newpos[,1])
max <- max(newpos[,1])
mean <- mean(newpos[,1])
range <- max-min
newpos[,1] <- (newpos[,1]-mean)/range*scale
}else{
##do nuttin.
}
}
newpos
}
##
setup <- function(knowledge,n)
{
##we won't allow ineligible states, so just
##keep track of knowledge states
#beliefs <<- sample(1:numstates,n,replace=T)
##set up two generally extreme groups
num <- nrow(knowledge)
# prob <- ((1:num)^2 + (num:1)^2)
# prob <- prob/sum(prob)
# print(prob)
##Sample uniformly across the states
beliefs<<-sample(1:num,n,replace=T)
# beliefs<<-sample(c(1,num),n,replace=T,prob=prob)
# beliefs<<-c(1:num,sample(c(1,num),n-num,replace=T))
# beliefs<<-c(rep(1,50),rep(15,50))
}
##note that reps is 10 for each n.
##keep =T restarts simulation.
run <- function(n=1000,keep=F)
{
if(keep)
{
dat <- matrix(0,numagents,(n+1))
dat[,1] <- beliefs
}
reps<- 10
for(i in 1:n+1)
{
for(j in 1:reps)
update(kl)
if(keep)
{
dat[,i] <- beliefs
}
}
# matplot(0:n,t(dat),type="p",pch=".",col=grey(thresh))
if(!keep)
{dat <- beliefs}
dat
}
update <- function(kl)
{
##choose two
pair <- sample.int(numagents,2)
##pick out their states in the knowledge lattice:
kID <- beliefs[pair]
## get their belief vectors
kpair <- kl[kID,]
## create candidate new belief vectors by
## changing each belief.
changeab <- runif(numfeatures)0)
{
beliefs[pair[1]] <<- xa
}
xb <- isvalid(newb,kl)
if(xb>0)
{
beliefs[pair[2]] <<- xb
}
}
isvalid <- function(belief, states)
{
as.numeric(apply(states,1,function(x){(sum(x==belief))==numfeatures})%*% 1:numstates)
}
#################################################
## This is a smiple visualization set-up
#################################################
set.seed(105)
numstates <<- 15
numfeatures <<- 20
numagents <<- 10
mu <<- .3 ##convergence parameter
##set up knowledge lattice so it is distributed
kl <<- setupKnowledgeLattice(numstates,numfeatures)
##set it up so it is polarized
#kl <<- setupKnowledgeLattice(numstates,numfeatures,type="polar")
precedence <- matrix(0,numstates,numstates)
for(i in 1:numstates)
for(j in 1:numstates)
{
precedence[i,j] <- precedes(kl[i,],kl[j,])
}
pedit <- editprec(precedence)
par(mfrow=c(1,2))
pts <-gplot(precedence)
#gplot(pedit,coord=pts)
#par(mfrow=c(1,1))
pos <- cbind(pts[,2],rowSums(kl));gplot(pedit,coord=pos,jitter=F)
#pos <- newpos(pos,pedit,10,byrow=T);gplot(pedit,coord=pos,jitter=F,edge.lty=1,edge.lwd=.1,vertex.cex=2.5,vertex.col="white")
vcol <- rgb(49/255,130/255,189/255) #"#3182BD"
pdf("lattice1.pdf",width=9,height=8)
par(mfrow=c(1,2))
par(mar=c(8,4,12,0))
image(1:ncol(kl),1:nrow(kl),t(kl),col=c("white",vcol),xaxt="n",yaxt="n",
xlab="Feature",ylab="Knowledge state",main="Features in belief states")
segments(0,0:15+.5,25.5,0:15+.5,col="grey")
segments(0:25+.5,0,0:25+.5,15.5,col="grey")
axis(2,1:15,las=1,cex.axis=.75,col="grey80")
par(mar=c(0,5,3,0))
pos <- newpos(pos,pedit,10,byrow=F)
gplot(pedit,coord=pos,jitter=F,edge.lty=1,edge.lwd=.1,vertex.cex=3,vertex.col=vcol,
suppress.axes=T,pad=-2,
main="",ylab="")
title(main="Lattice of Knowledge States",line=0)
title(ylab="Number of positive beliefs",line=-.5)
axis(2,0:20,las=1,line=-2.5)
text(pos[,1],pos[,2],1:numstates)
dev.off()
#dev.copy2eps()
#####################################################################
######## Simulation 1: Looks at consensus formation in distributed
######## knowledge spaces.
########
########
########
numstates <<- 15
numfeatures <<- 20
numagents <<- 100
mu <<- .3 ##convergence parameter
##Boost these numbers to have product of 1000+ for final simulation:
numspaces <- 5 #25 ## number of distinct spaces to simulate
numtries <- 5 #20 ## number of runs per space
runlog <- matrix(0,125,(numfeatures+1))
histlog <- matrix(0,numspaces*numtries,(numfeatures+1))
spacelog <- matrix(0,numspaces*numtries,(numstates))
iterations <- rep(0,numspaces*numtries)
set.seed(1001) ##set seed for convergence to one group
set.seed(1003) ##set a seed for convergence to two; see examples
##testing. Comment this out to see consistency across
##simulations
logi <- 1
for(space in 1:numspaces)
{
##set it up so it is polarized
kl <<- setupKnowledgeLattice(numstates,numfeatures,type="spread")
##compute precedence
precedence <- matrix(0,numstates,numstates)
for(i in 1:numstates)
for(j in 1:numstates)
{
precedence[i,j] <- precedes(kl[i,],kl[j,])
}
pedit <- editprec(precedence)
# pos <- cbind(pts[,2],rowSums(kl))
# pos <- newpos(pos,pedit,10,byrow=T) ##create the graph
posstr <- rowSums(kl)
lasty <- rep(F,numstates)
lasty2 <- rep(F,numstates)
lastyrun <- 0
sameas <- 0
sameas1 <- 0
sameas2 <- 0
ii <- 1
for(try in 1:numtries)
{
setup(kl,numagents)
cont <- 1
i <- 1
while(cont)
{
dat <- run(100,F)
#strength plots the belief strength
strength <- posstr[beliefs]
par(mfrow=c(1,2))
y <- hist(beliefs,breaks=1:(numstates+1)-.5,ylim=c(0,numagents),main=paste(lastyrun,sameas1*1,sameas2*1,sameas*1))
x <- hist(strength,breaks=0:(numfeatures+1)-.5,ylim=c(0,numagents),
main=paste(space,try,cont,i,"\n",logi,"of",numspaces*numtries))
runlog[i,] <- x$counts
ydist <- y$counts>5
# ydist2 <- y$counts>0
sameas1 <-sum(lasty==ydist)==length(ydist)
# sameas2 <- sum(lasty2==ydist2)==length(ydist2)
sameas <- sameas1#|sameas2
if(sameas)
{
lastyrun <- lastyrun+1
}else{
lastyrun <- 0
}
num <- sum(ydist)
i <- i + 1
if(num==1 |lastyrun >25 | i>125)
{
cont<-0
}
lasty <- ydist
# lasty2 <- ydist2
}
histlog[logi,] <- x$counts
spacelog[logi,] <- posstr
iterations[logi] <- i
logi <- logi +1
}
}
rl <- runlog[1:(nrow(runlog)-1),]
rlc <- apply(rl,1,cumsum)
matplot(t(rlc),type="l")
##this let's me save different particular runs for a demonstration figure below
if(0)
{
rcl1 <- rlc
rlc2 <- rlc
}
if(1)
{
write.csv(histlog,"sim1-histlogf.csv")
write.csv(spacelog,"sim1-spacelogf.csv")
write.csv(iterations,"sim1-iterationsf.csv")
}
histlog <- read.csv("sim1-histlogf.csv")[,2:22]
spacelog <- read.csv("sim1-spacelogf.csv")[,2:16]
iterations <- read.csv("sim1-iterationsf.csv")$x
histlog.trim <- histlog[iterations<=125,]
#postscript("sim1-dist.eps",width=6,height=6)
pdf("sim1-dist.pdf",width=6,height=6)
par(las=1,mfrow=c(1,1))
labs <- c(sum(iterations>125),table(rowSums(histlog.trim>=5)))
xvals <- barplot(labs/sum(labs)*100,col="#3182BD", xaxt="n",
main="",ylab="Percent of simulations ",
ylim=c(0,125),las=1,xlim=c(0,6), yaxt="n",
xlab="Number of groups at convergence")
text(xvals,(labs/nrow(histlog.trim)*100)+10,paste("n =",labs))
it.tab <- table(iterations)
it.x <- as.numeric(rownames(it.tab))
points(c(0,it.x)*6/125,100+cumsum(c(0,it.tab))/sum(it.tab)*100/4,type="l",lwd=3)
axis(3, 0:5*6/5,labels=0:5*25)
mtext("Cycles to convergence (x1000)",3,line=2)
mtext("Percentile",2,line=3,at=110,las=0)
axis(1,1:5*1.2-.5,labels=0:4)
axis(2,0:9*10)
axis(2,100+0:5*5,0:5*20)
abline(100,0)
dev.off()
hl1 <- histlog[rowSums(histlog>0)==1,]
hl2 <- histlog[rowSums(histlog>0)==2,]
hl3 <- histlog[rowSums(histlog>0)==3,]
hl4 <- histlog[rowSums(histlog>0)==4,]
pdf("sim1-dist2.pdf", width=9,height=6)
par(mfrow=c(1,1))
plot(colMeans(histlog>5),bty="n",type="l",ylab="Proportion of trials",las=1,xaxt="n",
ylim=c(0,1),lwd=3,
col="#3182BD",xlab="Issues in agreement with extreme view",pch="2")
tmp <- choose(20,0:20)
tmpsum <- sum(tmp)
points(tmp/tmpsum,lwd=3,col="grey",type="l")
axis(1,1:21,labels=0:20)
axis(2,0:5*20,las=1)
points(table(as.numeric(as.matrix(spacelog)))/nrow(spacelog),type="o",lwd=3,pch=16)
legend(9,1,c("Initial","Converged states","Binomial(20,n)"),
lwd=c(3,3,2),pch=c(16,NA,NA),
col=c("black","#3182BD","grey"),lty=1,bty="n")
dev.off()
#####################################################################
######## Simulation 2: Looks at polarized knowledge spaces, and
######## which state the things converge to.
########
######## Note: histogram is plotted by belief strength, not belief state;
######## this does not matter because convergence is typically to the most
######## extreme state.
########
########
numstates <<- 15
numfeatures <<- 20
numagents <<- 100
mu <<- .3 ##convergence parameter
numspaces <- 50
numtries <- 20
histlog <- matrix(0,numspaces*numtries,(numfeatures+1))
spacelog <- matrix(0,numspaces*numtries,(numstates))
iterations <- rep(0,numspaces*numtries)
logi <- 1
for(space in 1:numspaces)
{
##set it up so it is polarized
kl <<- setupKnowledgeLattice(numstates,numfeatures,type="polar")
##compute precedence
precedence <- matrix(0,numstates,numstates)
for(i in 1:numstates)
for(j in 1:numstates)
{
precedence[i,j] <- precedes(kl[i,],kl[j,])
}
pedit <- editprec(precedence)
posstr <- rowSums(kl)
for(try in 1:numtries)
{
setup(kl,numagents)
cont <- 1
i <- 1
while(cont)
{
dat <- run(100,T)
strength <- posstr[beliefs]
par(mfrow=c(1,2))
y <- hist(beliefs,breaks=0:(numstates+1)-.5)
x <- hist(strength,breaks=0:(numfeatures+1)-.5,ylim=c(0,numagents),
main=paste(space,try,cont,i,"\n",logi,"of",numspaces*numtries))
num <- sum(y$counts>0)
i <- i + 1
if(num==2 | i>50)
{
cont<-0
}
}
histlog[logi,] <- x$counts
spacelog[logi,] <- posstr
iterations[logi] <- i
logi <- logi +1
}
}
write.csv(histlog,"simext-histlog.csv")
write.csv(spacelog,"simext-spacelog.csv")
write.csv(iterations,"simext-iterations.csv")
histlog <- read.csv("simext-histlog.csv")[,2:22]
spacelog <- read.csv("simext-spacelog.csv")[,2:16]
iterations <- read.csv("simext-iterations.csv")$x
hl2 <- histlog>0
spacetab <- table(unlist(spacelog))/20
spacetab2 <- rep(0,21)
for(i in 1:length(spacetab))
{
print(i)
print(spacetab[i])
spacetab2[as.numeric(rownames(spacetab))[i]+1] <- spacetab[i]
}
pdf("extremes.pdf",width=9,height=6)
par(mfrow=c(1,1),las=1)
matplot(0:20 ,t(histlog),type="p",pch=16,xaxt="n",cex=.85,ylim=c(0,120),
xlab="Number of positive beliefs",ylab="Number of agents",
bty="L",col="grey",main="Development of extreme consensus opinions")
points(0:20,colMeans(histlog),type="l",lwd=3,col="#3182BD")
points(0:20,colMeans(hl2)*100,type="l",lwd=2,col="navy",lty=2)
points(0:20,spacetab2,type="l")
axis(1,0:20)
legend(6,80,c("Starting Distribution","Mean Outcome","Convergence Distribution"),
lty=c(1,1,2),lwd=c(1,2,3),col=c("black","#3182BD","navy"),
bty="n")
dev.off()
####################################################
## Adjacency matrix figure
####################################################
##Compute an adjacency matrix for all 15 points:
library(sna)
seed <- 3336
print(seed)
set.seed(seed)
numstates <- 15
numfeatures <- 20
kl <<- setupKnowledgeLattice(numstates,numfeatures,type="polar")
dall.1 <- as.matrix(dist(kl))^2
d2.1 <- as.matrix(dist(kl))^2<2
d3.1 <- as.matrix(dist(kl))^2<3
d4.1 <- as.matrix(dist(kl))^2<5
gall.1 <- gplot(dall.1,gmode="graph",label=rowSums(kl), edge.lwd=.2,mode="mds")
g2.1 <- gplot(d2.1,gmode="graph",label=rowSums(kl), edge.lwd=5)
#g3.1 <- gplot(d3.1,gmode="graph",edge.lwd=2,interactive=T,label=rowSums(kl))
g4.1 <- gplot(d4.1,gmode="graph",edge.lwd=.5,edge.col="grey")
lab.1 <- rowSums(kl)
set.seed(3338)
kl <<- setupKnowledgeLattice(numstates,numfeatures,type="spread")
dall.2 <- as.matrix(dist(kl))^2
d2.2 <- as.matrix(dist(kl))^2<2
d3.2 <- as.matrix(dist(kl))^2<3
d4.2 <- as.matrix(dist(kl))^2<5
gall.2 <- gplot(dall.2,gmode="graph",label=rowSums(kl), edge.lwd=.2,mode="mds")
g2.2 <- gplot(d2.2,gmode="graph",label=rowSums(kl), edge.lwd=5)
#g3.2 <- gplot(d3.2,gmode="graph",edge.lwd=2,interactive=T,label=rowSums(kl))
g4.2 <- gplot(d4.2,gmode="graph",edge.lwd=.5,edge.col="grey")
lab.2 <- rowSums(kl)
set.seed(3339)
kl <<- setupKnowledgeLattice(numstates,numfeatures,type="centered")
dall.3 <- as.matrix(dist(kl))^2
d2.3 <- as.matrix(dist(kl))^2<2
d3.3 <- as.matrix(dist(kl))^2<3
d4.3 <- as.matrix(dist(kl))^2<5
gall.3 <- gplot(dall.2,gmode="graph",label=rowSums(kl), edge.lwd=.2,mode="mds")
g2.3 <- gplot(d2.2,gmode="graph",label=rowSums(kl), edge.lwd=5)
#g3.3 <- gplot(d3.2,gmode="graph",edge.lwd=2,interactive=T,label=rowSums(kl))
g4.3 <- gplot(d4.2,gmode="graph",edge.lwd=.5,edge.col="grey")
lab.3 <- rowSums(kl)
svg("adjacency.svg",width=16,height=15)
par(mfrow=c(3,1))
gplot(d4.2,coord=g3.2,gmode="graph",edge.lwd=.5,edge.col="grey",
suppress.axes=T,ylim=c(-5,5))
gplot(d3.2,coord=g3.2,gmode="graph", edge.lwd=2,edge.col="grey40",new=F )
gplot(d2.2,coord=g3.2,gmode="graph",edge.lwd=8,vertex.cex=3.5,boxed.labels=F,displaylabels=F,
label.pos=99,vertex.sides=32,new=F,vertex.col="#3182BD")
text(g3.2[,1],g3.2[,2],lab.2)
legend(-5,-7,c("1 apart","2-3 apart", "4 apart"),col=c("black","grey40","grey"),lwd=c(5,2,1),bty="n")
gplot(d4.3,coord=g3.3,gmode="graph",edge.lwd=.5,edge.col="grey",
suppress.axes=T,ylim=c(-5,5))
gplot(d3.3,coord=g3.3,gmode="graph", edge.lwd=2,edge.col="grey40",new=F )
gplot(d2.3,coord=g3.3,gmode="graph",edge.lwd=8,vertex.cex=3.5,boxed.labels=F,displaylabels=F,
label.pos=99,vertex.sides=32,new=F,vertex.col="#3182BD")
text(g3.3[,1],g3.3[,2],lab.3)
gplot(d4.1,coord=g3.1,gmode="graph",edge.lwd=.5,edge.col="grey",
suppress.axes=T,ylim=c(-5,5))
gplot(d3.1,coord=g3.1,gmode="graph", edge.lwd=2,edge.col="grey40",new=F )
gplot(d2.1,coord=g3.1,gmode="graph",edge.lwd=8,vertex.cex=3.5,boxed.labels=F,displaylabels=F,
label.pos=99,vertex.sides=32,new=F,vertex.col="#3182BD")
text(g3.1[,1],g3.1[,2],lab.1)
dev.off()
#####################################################################
######## Simulation 3: Knowledge distribution is centered at p=.5.
########
########
########
########
numstates <<- 15
numfeatures <<- 20
numagents <<- 100
mu <<- .3 ##convergence parameter
numspaces <- 2 #50
numtries <- 2 #20
histlog <- matrix(0,numspaces*numtries,(numfeatures+1))
spacelog <- matrix(0,numspaces*numtries,(numstates))
iterations <- rep(0,numspaces*numtries)
numendstates <- rep(0,numspaces*numtries)
logi <- 1
for(space in 1:numspaces)
{
##set it up so it is polarized
kl <<- setupKnowledgeLattice(numstates,numfeatures,type="centered")
##compute precedence
precedence <- matrix(0,numstates,numstates)
for(i in 1:numstates)
for(j in 1:numstates)
{
precedence[i,j] <- precedes(kl[i,],kl[j,])
}
pedit <- editprec(precedence)
posstr <- rowSums(kl)
lasty <- rep(F,numstates)
lasty2 <- rep(F,numstates)
lastyrun <- 0
sameas <- 0
sameas1 <- 0
sameas2 <- 0
for(try in 1:numtries)
{
setup(kl,numagents)
cont <- 1
i <- 1
while(cont)
{
dat <- run(100,T)
strength <- posstr[beliefs]
par(mfrow=c(1,2))
y <- hist(beliefs,breaks=1:(numstates+1)-.5,ylim=c(0,numagents),main=paste(lastyrun,sameas1*1,sameas2*1,sameas*1))
x <- hist(strength,breaks=0:(numfeatures+1)-.5,ylim=c(0,numagents),
main=paste(space,try,cont,i,"\n",logi,"of",numspaces*numtries))
ydist <- y$counts>5
# ydist2 <- y$counts>0
sameas1 <-sum(lasty==ydist)==length(ydist)
# sameas2 <- sum(lasty2==ydist2)==length(ydist2)
sameas <- sameas1#|sameas2
if(sameas)
{
lastyrun <- lastyrun+1
}else{
lastyrun <- 0
}
num <- sum(ydist)
i <- i + 1
if(num==1 |lastyrun >25 | i>125)
{
cont<-0
}
lasty <- ydist
# lasty2 <- ydist2
}
histlog[logi,] <- x$counts
spacelog[logi,] <- posstr
iterations[logi] <- i
numendstates[logi] <- num
logi <- logi +1
}
}
if(1)
{
write.csv(histlog,"sim3b-histlog.csv")
write.csv(spacelog,"sim3b-spacelog.csv")
write.csv(iterations,"sim3b-iterations.csv")
write.csv(numendstates,"sim3b-numendstates.csv")
}
histlog <- read.csv("sim3-histlog.csv")[,2:22]
spacelog <- read.csv("sim3-spacelog.csv")[,2:16]
iterations <- read.csv("sim3-iterations.csv")$x
numendstates <- read.csv("sim3-numendstates.csv")$x
sum1 <- rowSums(histlog>5)
histlog.trim <- histlog[iterations<=125,]
postscript("sim3-dist.eps",width=6,height=6)
pdf("sim3-dist.pdf",width=6,height=6)
par(las=1)
labs <- c(sum(iterations>125),table(numendstates))
x <- barplot(labs/sum(labs)*100,col="#3182BD", xaxt="n",
main="",ylab="Percent of simulations ",
ylim=c(0,125),las=1,xlim=c(0,7), yaxt="n",
xlab="Number of groups at convergence")
text(x,(labs/nrow(histlog.trim)*100)+10,paste("n =",labs))
it.tab <- table(iterations)
it.x <- as.numeric(rownames(it.tab))
points(c(0,it.x)*6/125,100+cumsum(c(0,it.tab))/sum(it.tab)*100/4,type="l",lwd=3)
axis(3, 0:5*6/5,labels=0:5*25)
mtext("Cycles to convergence (x1000)",3,line=2)
mtext("Percentile",2,line=3,at=110,las=0)
axis(1,1:6*1.2-.5,labels=0:5)
axis(2,0:9*10)
axis(2,100+0:5*5,0:5*20)
abline(100,0)
dev.off()
postscript("sim3-dist2.eps", width=9,height=6)
par(mfrow=c(1,1))
plot(0:20,colMeans(histlog>0),bty="n",type="l",ylab="Proportion of trials",las=1,xaxt="n",
ylim=c(0,.5),lwd=2,
col="#3182BD",xlab="Issues in agreement with extreme view",pch="2")
tmp <- choose(20,0:20)
tmpsum <- sum(tmp)
points(0:20,tmp/tmpsum,lwd=3,col="grey",type="l")
axis(1,0:20)
axis(2,0:5*20,las=1)
tmp2 <-table(as.numeric(as.matrix(spacelog)))/nrow(spacelog)/15
points(as.numeric(rownames(tmp2)),tmp2,type="o",lwd=3,pch=16)
legend(2,.5,c("Initial","Binomial(20,n)","Converged"),
lwd=c(2,2,3), pch=c(NA,NA,16),
col=c("#3182BD","grey","black"),
lty=1,bty="n")
dev.off()
hist(numendstates$x[iterations$x<=125],breaks=0:5+.5)
histlog.keep <- histlog[iterations$x <= 125,]
colMeans(histlog.keep>=5)[2:22]
n2 <- histlog[numendstates$x==2,]
plot(table(as.numeric(as.matrix(spacelog[,2:16])))/nrow(spacelog),type="o",lwd=3,pch=16,xlab="Degree of agreement",ylab="percent")
#####################################################################
######## Simulation 4: Empirical distribution
########
########
########
numstates <<- 15
numfeatures <<- 20
numagents <<- 100
mu <<- .3 ##convergence parameter
numspaces <- 10000
spacelog <- matrix(0,numspaces,(numstates))
for(i in 1:numspaces)
{
kl <- setupKnowledgeLattice(numstates,numfeatures,type="spread")
spacelog[i,] <- rowSums(kl)
}
plot(table(as.numeric(as.matrix(spacelog)))/nrow(spacelog),type="o",lwd=3,pch=16,xlab="Degree of agreement",ylab="percent")
####################################################
## Other stuff
####################################################
histlog <- matrix(0,numstates,ncol(dat))
for(i in 1:ncol(dat))
{
histlog[,i] <- hist(dat[,i],breaks=0:numstates+.5,plot=F)$counts
}
cols <- rep(c("#3182BD","grey"),20)
cols<-1:25
#histlog <- rl
#cumhist <-rbind(0,apply(histlog,1,cumsum))
cumhist1 <- rlc1 ##from experiment 1
cumhist1 <- rbind(rep(0,ncol(cumhist1)),
cumhist1,
rep(100,ncol(cumhist1)))
cumhist2 <- rlc2 ##from experiment 1
cumhist2 <- rbind(rep(0,ncol(cumhist2)),
cumhist2,
rep(100,ncol(cumhist2)))
#postscript("agenthist.eps",width=10,height=5)
pdf("agenthist.pdf",width=8,height=6)
par(mfrow=c(2,1),mar=c(4,3,1,1))
cols <- brewer_pal("seq")(9)
matplot(t(cumhist1),type="l",xlab="",ylab="",bty="n",
lty=1,col="black",las=1)
title(xlab="Interactions (x1000)",line=2,cex.lab=.8)
title(ylab="Number of agents",line=2.2, cex.lab=.8)
colid <- 2
for(i in 1:(nrow(cumhist1)-1))
{
area <- sum(cumhist1[i+1,]-cumhist1[i,])
if(area>=0)
{
polygon(c(1:ncol(cumhist1),ncol(cumhist1):1),c(cumhist1[i,],rev(cumhist1[i+1,])),lwd=.5,col=cols[colid])
colid <- (colid%%9 + 1)
}
}
matplot(t(cumhist2),type="l",xlab="",ylab="",bty="n",
lty=1,col="black",las=1)
title(xlab="Interactions (x1000)",line=2,cex.lab=.8)
title(ylab="Number of agents",line=2.2, cex.lab=.8)
colid <- 2
for(i in 1:(nrow(cumhist2)-1))
{
area <- sum(cumhist2[i+1,]-cumhist2[i,])
if(area>=0)
{
polygon(c(1:ncol(cumhist2),ncol(cumhist2):1),c(cumhist2[i,],rev(cumhist2[i+1,])),lwd=.5,col=cols[colid])
colid <- (colid%%9 + 1)
}
}
# matplot(t(cumhist),type="l",xlab="Interactions (x10)",ylab="Number of agents",bty="n",
# lty=1,col="black",add=T,lwd=.5)
dev.off()
PK «©ML•ˆp
¤ codemeta.jsonPK «©ML;ˆ(î î ¤, docs/documentation.txtPK «©ML;ˆ(î î ¤N code/documentation.txtPK «©ML–HJÿ` ` ! ¤p code/MuellerKnowledgeSpaceModel.RPK Âo