So I am currently working on my master thesis on evolutionary games theory. I am running a simulation of population that plays evolutionary trust game in R. I have two populations of investors and trustees which interact between each other by the means of trust game. Nevertheless I am currently gettin some counterintuitive results from the simulation (strategy X that is being dominated by the strategy Y strategy performs better than strategy Y) and one strategy type does not evolve at all. What I do is that I run the simulation of population of investors and trustees (both of size 500) for 1000 generations and in every generation investor meets 24 random trustees. I suspect that the issue is that I generate all of the random numbers using runif (I tried sample function as well) and that after several simulations the random numbers follow same pattern. Every generation (1000 in total) i generate 24*500 random numbers. I therefore wanted to ask how to generate true random numbers time efficiently in R. I tried random package but it is too slow and after 4 generations the function refuses to generate random numbers (Error in randomNumbers(1, 1, 300, 1, 10, TRUE) : random.org suggests to wait until tomorrow). I am currently trying to implement inverse sampling method but I don't expect it to be much faster.
Thanks for help.
The code I use currently, if anyone is interested
#### VECTORIZED VERSION ####
### loading libraries ###
require(plyr)
require(grDevices)
require(erer)
require(random)
#### defining trustees ####
# [1] honor/abuse
# [2] information about previous interaction - relevant to investors who buy the information
# [3] payoff
# [4] type of trustee
trustee1 <- c(1,0,0,1)
trustee2 <- c(0,0,0,2)
#### defining investors ####
# [1] buy/not buy
# [2-4] investment decision if buying information in case 1(TH), case 2(TA), case 3(NT)
# [5] aggregated payoff
# [6] type of the investor
investor1 <- c(1,1,1,1,0,1)
investor2 <- c(1,1,0,1,0,2)
investor3 <- c(1,1,0,0,0,3)
investor4 <- c(1,1,1,0,0,4)
investor5 <- c(1,0,0,0,0,5)
investor6 <- c(1,0,1,1,0,6)
investor7 <- c(1,0,1,0,0,7)
investor8 <- c(1,0,0,1,0,8)
investor9 <- c(0,0,0,0,0,9)
investor10 <-c(0,1,1,1,0,10)
### defining trust game as a R function ###
BFX <- 10
trustgame2 <- function(investor_string,trustee_string,cost){
BF <- 10
#investor <- get(investor_string, envir = globalenv())
#trustee <- get(trustee_s tring, envir = globalenv())
investor <- investors[[investor_string]]
trustee <- trustees[[trustee_string]]
investordecision <- NULL
trusteedecision <- trustee[1]
if (investor[1]==0) {investordecision <- investor[4]}
if (investor[1]==1){
if (trustee[2]==1) {investordecision <- investor[2]
investor[5]<- investor[5] - cost}
if (trustee[2]==2) {investordecision <- investor[3]
investor[5]<- investor[5] - cost}
if (trustee[2]==3) {investordecision <- investor[4]
investor[5]<- investor[5] - cost}
if (trustee[2]==0) {investordecision <- investor[4]
investor[5]<- investor[5] - cost}
}
if ((investordecision==1) && (trustee[1]==1)) trustee[2] <- 1
if ((investordecision==1) && (trustee[1]==0)) trustee[2] <- 2
if ((investordecision==0) && (trustee[1]==1)) trustee[2] <- 3
if ((investordecision==0) && (trustee[1]==1)) trustee[2] <- 3
if ((investordecision==1) && (trusteedecision==1)){
trustee[3] <- trustee[3] +3 + BF
investor[5] <- investor[5] + 3 + BF
}
if ((investordecision==1) && (trusteedecision==0)){
trustee[3] <- trustee[3] + 5 + BF
investor[5] <- investor[5] + BF
}
if ((investordecision==0) && (trusteedecision==0)){
trustee[3] <- trustee[3] +1 + BF
investor[5] <- investor[5] + 1 + BF
}
if ((investordecision==0) && (trusteedecision==1)){
trustee[3] <- trustee[3] +1 + BF
investor[5] <- investor[5] +1 + BF
}
#assign(trustee_string, trustee, envir = globalenv())
#assign(investor_string, investor, envir = globalenv())
investors[[investor_string]] <<- investor
trustees[[trustee_string]] <<- trustee
}
### function for generating the first population of investors ###
investorspop <- function(n)
{investors <- list()
for (i in 1 : n)
{
if ((i %% 10) == 1) investors[[i]] <- investor1
if ((i %% 10) == 2) investors[[i]] <- investor2
if ((i %% 10) == 3) investors[[i]] <- investor3
if ((i %% 10) == 4) investors[[i]] <- investor4
if ((i %% 10) == 5) investors[[i]] <- investor5
if ((i %% 10) == 6) investors[[i]] <- investor6
if ((i %% 10) == 7) investors[[i]] <- investor7
if ((i %% 10) == 8) investors[[i]] <- investor8
if ((i %% 10) == 9) investors[[i]] <- investor9
if ((i %% 10) == 0) investors[[i]] <- investor10
}
shuffleI(investors)
}
### generating and setting up the names for the population of investors###
investors <- investorspop(10)
names(investors) <- paste0("investor",seq_along(investors))
### function for generating the first population of trustees###
trusteespop <- function(n)
{ trustees <- list()
for (i in 1 : n)
{if ((i %% 2) == 1) trustees[[i]] <- trustee1
if ((i %% 2) == 0) trustees[[i]] <- trustee2
}
names(trustees) <- paste0("trustee",seq_along(trustees))
shuffleT(trustees)
}
### generating and setting up the names for the population of trustees###
trustees <- trusteespop(10)
names(trustees) <- paste0("trustee",seq_along(trustees))
### prepares populations of investors and trustees with n individuals ###
prepare_population <- function(popsizeI,popsizeT)
{investors <<- investorspop(popsizeI)
trustees <<- trusteespop(popsizeT)
}
### simulating one generation of the model ###
simulategeneration <- function(ninteractions=24,investorspopulation,trusteespopulation,cost)
{ for (i in 1 : ninteractions )
{ shuffleI(investorspopulation)
vector <- as.numeric(randomNumbers(300,1,300,1,10,TRUE))
for (j in 1 : length(investorspopulation)){
trustgame2(paste0("investor",j),paste0("trustee",vector[j]),cost)}
}
}
### evolution of the new generation ###
# total number of points earned in one generation #
summedpayoff <- Reduce("+",investors)[5]
### function that gets the fractions of the investors for the next generation ###
get_fractionsI <- function(investors)
{
fractions <- rep(0,10)
fractions1 <- NULL
slots <- NULL
j <- 0
payoffs <- NULL
investors_table <- data.frame(matrix(unlist(investors),nrow=length(investors),byrow=T))
colnames(investors_table) <- c("B","TH","TA","NT","PAYOFF","TYPE")
payoffs <- aggregate(investors_table$PAYOFF, by=list(Category=investors_table$TYPE), FUN=sum)$x
slots <- count(investors_table,'TYPE')$TYPE
summedpayoff <- Reduce("+",investors)[5]
fractions1 <- payoffs / summedpayoff
fractions[slots] <- fractions1
fractions
}
get_fractionsT <- function(trustees)
{fractions <- rep(0,2)
fractions1 <- NULL
slots <- NULL
j <- 0
payoffs <- NULL
trustees_table <- data.frame(matrix(unlist(trustees),nrow=length(trustees),byrow=T))
colnames(trustees_table) <- c("H/A","INFORMATION","PAYOFF","TYPE")
payoffs <- aggregate(trustees_table$PAYOFF, by=list(Category=trustees_table $TYPE), FUN=sum)$x
slots <- count(trustees_table,'TYPE')$TYPE
summedpayoff <- Reduce("+",trustees)[3]
fractions1 <- payoffs / summedpayoff
fractions[slots] <- fractions1
fractions
}
count_investors <- function(investors)
{ number <-rep(0,10)
number1 <- NULL
slots <- NULL
j <- 0
investors_table <- data.frame(matrix(unlist(investors),nrow=length(investors),byrow=T))
colnames(investors_table) <- c("B","TH","TA","NT","PAYOFF","TYPE")
number1 <- count(investors_table,'TYPE')$freq
slots <- count(investors_table,'TYPE')$TYPE
number[slots] <- number1
number
}
count_trustees <- function(trustees)
{ number <-rep(0,2)
number1 <- NULL
slots <- NULL
j <- 0
trustees_table <- data.frame(matrix(unlist(trustees),nrow=length(trustees),byrow=T))
colnames(trustees_table) <- c("H/A","INFORMATION","PAYOFF","TYPE")
number1 <- count(trustees_table,'TYPE')$freq
slots <- count(trustees_table,'TYPE')$TYPE
number[slots] <- number1
number
}
new_generationI <- function(fractionsI,popsizeI)
{ investors <- list()
cumsumsI <- cumsum(round(fractionsI*(popsizeI)))
for (i in 1 : popsizeI)
{ if (i < cumsumsI[1] + 1) investors[[i]] <- investor1
if ((i > cumsumsI[1]) & (i < cumsumsI[2]+1)) investors[[i]] <- investor2
if ((i > cumsumsI[2]) & (i < cumsumsI[3]+1)) investors[[i]] <- investor3
if ((i > cumsumsI[3]) & (i < cumsumsI[4]+1)) investors[[i]] <- investor4
if ((i > cumsumsI[4]) & (i < cumsumsI[5]+1)) investors[[i]] <- investor5
if ((i > cumsumsI[5]) & (i < cumsumsI[6]+1)) investors[[i]] <- investor6
if ((i > cumsumsI[6]) & (i < cumsumsI[7]+1)) investors[[i]] <- investor7
if ((i > cumsumsI[7]) & (i < cumsumsI[8]+1)) investors[[i]] <- investor8
if ((i > cumsumsI[8]) & (i < cumsumsI[9]+1)) investors[[i]] <- investor9
if ((i > cumsumsI[9]) & (i < cumsumsI[10]+1)) investors[[i]] <- investor10
}
if (length(investors) < popsizeI) for (j in (length(investors)+1) : popsizeI) investors[[j]] <- get(paste0("investor",round(runif(1,0.5,10.5))))
names(investors) <- paste0("investor",seq_along(investors))
shuffleI(investors)
}
new_generationT <- function(fractionsT,popsizeT)
{trustees <- list()
cumsumsT <- cumsum(round(fractionsT*popsizeT))
for (i in 1 : popsizeT)
{
if (i < cumsumsT[1]+1) trustees[[i]] <- trustee1
if ((i > cumsumsT[1] ) && (i < cumsumsT[2]+1)) trustees[[i]] <- trustee2
}
if (length(trustees) < popsizeT) for (j in (length(trustees)+1) : popsizeT) trustees[[j]] <- get(paste0("trustee",round(runif(1,0.5,2.5))))
names(trustees) <- paste0("trustee",seq_along(trustees))
shuffleT(trustees)
}
shuffleI <- function(investors)
{investors <- sample(investors)
names(investors) <- paste0("investor",seq_along(investors))
investors
}
shuffleT <- function(trustees)
{trustees <- sample(trustees)
names(trustees) <- paste0("trustee",seq_along(trustees))
trustees
}
mutation <- function(investors,prob)
{for (q in 1 : length(investors))
{ random <- runif(1,0,1)
if (random < prob) investors[[q]] <- get(paste0("investor",round(runif(1,0.5,10.5))))
}
investors}
### simulating one generation of the population ###
simulategeneration(ninteractions = 24,investorspopulation = investors,trusteespopulation = trustees,cost=0)
### RUNNING THE WHOLE SIMULATION OF THE POPULATION ###
popsizeI <- 300
popsizeT <- 300
ninteractions <- 24
investors <- investorspop(popsizeI)
trustees <- trusteespop(popsizeT)
cost <- 0.2
numberI <- list()
numberT <- list()
fractionsil <- list()
fractionstl <- list()
for (i in 1:500)
{ numberI[[i]] <- count_investors(investors)
numberT[[i]] <- count_trustees(trustees)
simulategeneration(ninteractions,investors,trustees,cost)
fractionsI <- get_fractionsI(investors)
fractionsT <- get_fractionsT(trustees)
fractionsil[[i]] <- fractionsI
fractionstl[[i]] <- fractionsT
investors <- new_generationI(fractionsI,popsizeI)
trustees <- new_generationT(fractionsT,popsizeT)
mutation(investors,0.01)
}
### plotting the evolution of the investors ###
colourschemeI <- NULL
colourschemeT <- NULL
colourschemeI[1] <- "turquoise"
colourschemeI[2] <- "red"
colourschemeI[3] <-"black"
colourschemeI[4] <- "violet"
colourschemeI[5] <- "azure4"
colourschemeI[6] <- "darkgoldenrod3"
colourschemeI[7] <- "green"
colourschemeI[8] <- "darkblue"
colourschemeI[9] <- "mediumpurple4"
colourschemeI[10] <- "violetred3"
numbersI_table <- data.frame(matrix(unlist(numberI),nrow=length(numberI),byrow=T))
plot(numbersI_table[,1],col=colourschemeI[1],type='l',xlab="Generation",ylab="Size of the population",ylim=c(0,popsizeI),main=paste0("Simulation of population of ",popsizeI," investors and ",popsizeT," trustees. \n Price of information is ",cost,". \n The background fitness parameter is ",BFX))
for (i in 2 : 10) lines(numbersI_table[,i],col=colourschemeI[i] ,type='l',xlab="Generation",ylab="Size of the population")
colourschemeT[1] <- "darkgray"
colourschemeT[2] <- "lightsalmon2"
helpvectorI <- c("BTTT","BTNT","BTNN","BTTN","BNNN","BNTT","BNTN","BNNT","NNNN","NTTT")
numbersT_table <- data.frame(matrix(unlist(numberT),nrow=length(numberT),byrow=T))
#plot(numbersT_table[,1],col=colourscheme[1],type='l',xlab="Generation",ylab="Size of the population",ylim=c(0,500),main=paste0("Simulation of population of ",popsizeT," trustees. Price of information is ",cost,". \n The background fitness parameter is ",BFX))
for (j in 1:2) lines(numbersT_table[,j],col=colourschemeT[j] ,type='o',xlab="Generation",ylab="Size of the population")
legend(425,315,paste0("trustee",c("H","A")),lty=c(1,1),lwd=c(2.5,2.5),col=colourschemeT,cex=0.6)
legend(150,415,paste0("investor",helpvectorI),lty=c(1,1),lwd=c(2.5,2.5),col=colourschemeI,cex=0.6)
write.list(numberI, "investors.txt", t.name = NULL, row.names = FALSE)
write.list(numberT, "trustees.txt", t.name = NULL, row.names = FALSE)
investor1 <- c(1,1,1,1,0,1)
investor2 <- c(1,1,0,1,0,2)
investor3 <- c(1,1,0,0,0,3)
investor4 <- c(1,1,1,0,0,4)
investor5 <- c(1,0,0,0,0,5)
investor6 <- c(1,0,1,1,0,6)
investor7 <- c(1,0,1,0,0,7)
investor8 <- c(1,0,0,1,0,8)
investor9 <- c(0,0,0,0,0,9)
investor10 <-c(0,1,1,1,0,10)
### testing, this part is not necessary to check
prepare_population(10,10)
simulategeneration(24,investors,trustees,0)
fractionsI <- get_fractionsI(investors)
fractionsT <- get_fractionsT(trustees)
investors<- new_generationI(fractionsI,10)
trustees<- new_generationT(fractionsT,10)
## cooperativeI - 1-4 --> 100, rest 15 and 20
## uncooperativeI 5-8 --> 100, rest 15 and 20
## cooperativenoI 1-4,10 --> 60, rest 40
trustgame2("investor10","trustee10",0)
trustgame2("investor9","trustee9",0)
trustgame2("investor8","trustee8",0)
trustgame2("investor7","trustee7",0)
trustgame2("investor6","trustee6",0)
trustgame2("investor5","trustee5",0)
trustgame2("investor4","trustee4",0)
trustgame2("investor3","trustee3",0)
trustgame2("investor2","trustee2",0)
trustgame2("investor1","trustee1",0)
investors <- shuffleI(investors)
Aucun commentaire:
Enregistrer un commentaire