OR/16/022 Appendix 1 - NetLogo associated scripts

Beriro, D, Cave, M, Wragg, J and Hughes, A. 2016. Agent Based Modelling: Initial assessment for use on soil bioaccessibility. British Geological Survey Internal Report, OR/16/022.

NetLogo Script for the Northampton model

NetLogo Script

extensions [gis]

globals [elevation-dataset]

patches-own [elevation]

turtles-own [exposure group]

to setup

let xHm [-7 6 -5 3]

let yHm [6 5 -5 -6]

set-default-shape turtles "person"

set elevation-dataset gis:load-dataset "C:/Users/mrca/My Documents/Nthpton_As.asc"

gis:paint elevation-dataset 0

gis:apply-raster elevation-dataset elevation

let min-elevation gis:minimum-of elevation-dataset

let max-elevation gis:maximum-of elevation-dataset

ask patches

[ if (elevation <= 0) or (elevation >= 0)

[set pcolor scale-color black elevation min-elevation max-elevation]

]

ask patches

[ if pcolor = black

[set elevation 0]

]

create-turtles numturtles [

set exposure 0

set xcor -7

set ycor 6

set group 1

set color red

]

create-turtles numturtles [

set exposure 0

set xcor 6

set ycor 5

set group 2

set color green

]

create-turtles numturtles [

set exposure 0

set xcor -5

set ycor -5

set group 3

set color blue

]

create-turtles numturtles [

set exposure 0

set xcor 3

set ycor -6

set group 4

set color yellow

]

end to go

ask turtles with [group = 1][

ifelse distancexy -7 6 > dist-to-home

[facexy -7 6]

[rt random 360]  ; set random heading

forward 1  ; advance one step

]

ask turtles with [group = 2][

ifelse distancexy 6 5 > dist-to-home

[facexy 6 5]

[rt random 360]  ; set random heading

forward 1  ; advance one step

]

ask turtles with [group = 3][

ifelse distancexy -5 -5 > dist-to-home

[facexy -5 -5]

[rt random 360]  ; set random heading

forward 1  ; advance one step

]

ask turtles with [group = 4][

ifelse distancexy 3 -6 > dist-to-home

[facexy 3 -6]

[rt random 360]  ; set random heading

forward 1  ; advance one step

]

ask patches [

let my-turtles turtles-here

if any? my-turtles [

let lvls elevation / 100

print lvls

ask my-turtles

[set exposure (exposure + lvls)]

]

]

End.

R Script for the NetLogo Northampton Model

library(RNetLogo)

# this part starts up the NetLogo language

nl.path <-"C:/Program Files/NetLogo 5.3/app"

NLStart(nl.path,gui=TRUE, nl.obj=NULL, is3d=FALSE)

# now get the file path for the model we ant to run

model.path<-file.choose()

NLLoadModel(model.path)

# This part sets up the model setting the number of turtles at each site

# and the maximum distance from home the agents are allowed to move

NLCommand("clear-all")

NLCommand("set numturtles 100")

NLCommand("set dist-to-home 3")

NLCommand("setup")

# Now we need to set up some dataframes to store the arsenic exposure

# for each agent in each group

df1<-data.frame(nms=as.character(1:100))

df2<-data.frame(nms=as.character(1:100))

df3<-data.frame(nms=as.character(1:100))

df4<-data.frame(nms=as.character(1:100))

# this loop runs the NetLogo Program over 100 cycles collecting the exposure

# from each agent in each group at each cycle

for(i in 1:100){

NLDoCommand(1, "go")

Exp1<- NLGetAgentSet("exposure","turtles with [group = 1]")

df1<-cbind(df1,Exp1)

Exp2<- NLGetAgentSet("exposure","turtles with [ group = 2 ]")

df2<-cbind(df2,Exp2)

Exp3<- NLGetAgentSet("exposure","turtles with [ group = 3 ]")

df3<-cbind(df3,Exp3)

Exp4<- NLGetAgentSet("exposure","turtles with [ group = 4 ]")

df4<-cbind(df4,Exp4)

}

# now the results can be plotted out

###################################

# This plots the cumulative expusure for each agent in each group

# along with an average value for each group

par(mfrow=c(2,2))

avg1<-sapply(df1[,2:101],FUN=mean)

plot(1:100,avg1,type="l", ylim=c(0,3.5),lty=2,lwd=4,

xlab="time",ylab="Exposure", main="Group 1")

for( i in 1:100){

lines(1:100,df1[i,2:101],col=83+i)

}

lines(1:100,avg1,lty=2,lwd=3)

avg2<-sapply(df2[,2:101],FUN=mean)

plot(1:100,avg2,type="l", ylim=c(0,3.5),lty=2,lwd=2,

xlab="time",ylab="Exposure", main="Group 2")

for( i in 1:100){

lines(1:100,df2[i,2:101],col=83+i)

}

lines(1:100,avg2,lty=2,lwd=3)

avg3<-sapply(df3[,2:101],FUN=mean)

plot(1:100,avg3,type="l", ylim=c(0,3.5),lty=2,lwd=2,

xlab="time",ylab="Exposure", main="Group 3")

for( i in 1:100){

lines(1:100,df3[i,2:101],col=83+i)

}

lines(1:100,avg3,lty=2,lwd=3)

avg4<-sapply(df4[,2:101],FUN=mean)

plot(1:100,avg4,type="l", ylim=c(0,3.5),lty=2,lwd=2,

xlab="time",ylab="Exposure", main="Group 4")

for( i in 1:100){

lines(1:100,df4[i,2:101],col=83+i)

}

lines(1:100,avg4,lty=2,lwd=3)

###################################################

# This does a box and whisker plot for the final exposures over

# 100 cycles par(mfrow=c(1,1))

alldat<-cbind(df1[,101],df2[,101],df3[,101],df4[,101])

boxplot(alldat, xlab="Group",ylab="As Exposure")

###################################################

# This does a probability density function of the final exposure in each

# group

plot(density(alldat[,1]),main="Exposure PDF",

xlab="Cumulative As Exposure",col="green",

ylim=c(0,2),xlim=c(1.5,4.5),lwd=2)

lines(density(alldat[,2]),col="red",lwd=2)

lines(density(alldat[,3]),col="blue",lwd=2)

lines(density(alldat[,4]),col="purple",lwd=2)

legend("topleft",c("Group1","Group2","Group3","Group4"),

col=c("green","red","blue","purple"),bty="n",lty=1,lwd=2)

#################################################

# this quits the NetLogo program

NLQuit()

.