Vector Tutorial (NetLogo section) - NetLogo/GIS-Extension GitHub Wiki

If you had any issues getting any of these GIS steps to work, I have provided final versions of all the files you'll need to import into NetLogo so you can continue on with the tutorial if you wish.

Now that we're finally done with data preprocessing, we can get to work on actually using the GIS extension. Go ahead an open up a new NetLogo file and write extensions [gis] at the top of the code tab. This tells NetLogo that we want to be using the GIS extension and will allow us to use all the primitives defined by the GIS extension if we preface them with gis:.

To get our feet wet, let's start by importing our states shapefile to use as a background/basemap. Note that like all other file import methods in NetLogo, the filenames we use will be relative to the .nlogo file itself. I've set up my project so that all the shapefiles I am going to be using are in the same directory as the .nlogo file, but you might have to change the file paths if you've set up your project differently for some reason.

extensions [gis]           ; load the GIS extension
globals [states-dataset]   ; create a global variable to store our states dataset

to load-datasets
	clear-all
	set states-dataset gis:load-dataset "states_continental.shp"
end

This code simply loads in the states_continental.shp file and keeps it in a global variable so we can use it anywhere.

While most NetLogo models use a setup and go paradigm, when working with GIS datasets, I like to split up my reset-the-model code from my load-in-my-data code, just because when you're working with large datasets, re-parsing the datasets can be time consuming and slow down your setup procedure. As such, I use clear-all in the load-datasets procedure and selectively use clear-xxxx primitives in my setup procedure like so. The one thing you have to remember when structuring your code in this way is to always manually reset any global variables you may be using aside from your dataset variables.

to custom-clear
  reset-ticks
	clear-turtles
	clear-patches
	clear-drawing
	clear-all-plots
	;; reset any other global variables you might be using here
end

to setup
  custom-clear
  ; We'll be adding our map drawing code here shortly
  ...
end

We're almost ready to draw our states, but first we need to tell the GIS Extenison one last thing: what part of our globe we want to display. The term for this is "envelope", but if it helps you can just substitute our "envelope" for "bounding box" in your head. Either way, it is a set of two points that describe the extents of what we want to show. While you can certainly set these numbers manually, it is much simpler to just tell the extension to just "zoom to fit all of this dataset" and it'll calculate the appropriate values accordingly. So before we draw anything, we should tell NetLogo to set its world envelope to the envelope of the states dataset by adding this to the setup block:

gis:set-world-envelope gis:envelope-of states-dataset

Now we can finally draw our states. First we set a drawing color (using the NetLogo colors you are familiar with), then we use the gis:draw command, supplying the dataset we want to draw and a line width as arguments.

gis:set-drawing-color yellow
gis:draw states-dataset 1

Finally, add a load-datasets button and a setup button to the interface and set the world dimensions to 32 x 16, and while we're here, change the view updates from "continuous" to "on ticks". At this point you should be able to click load-datasets and then setup and see a map drawn like this.

Before we go on, here's all the code so far in one place.

extensions [gis]

globals [states-dataset ]

to load-datasets
  clear-all
  set states-dataset gis:load-dataset "states_continental.shp"
end

to custom-clear
  reset-ticks
	clear-turtles
	clear-patches
	clear-drawing
	clear-all-plots
	;; reset any other global variables you might be using here
end

to setup
  custom-clear
  
  gis:set-world-envelope gis:envelope-of states-dataset ;; set up the mapping between map/gis space to NetLogo world space

  gis:set-drawing-color yellow
  gis:draw states-dataset 1   ;; draw the states dataset using a line with a width of 1
end

Our next step is to set up our airport agents, which will hatch airplanes and keep track of how many infected individuals are present in their vicinity. We will soon turn each point in the point dataset of airports into its own NetLogo agent, but first let's just import the dataset and draw it using the gis:draw command to ensure it is being loaded correctly.

As before, to load a dataset you need to create a global variable to store the dataset and then set the value of that global variable to be the result of gis:load-dataset.

...
globals [states-dataset airports-dataset]

to load-datasets
  clear-all
  set states-dataset gis:load-dataset "states_continental.shp"
  set airports-dataset gis:load-dataset "final_airports_with_populations.shp"
end
...

To verify that the airport points were loaded like we expected them to, we can use gis:draw like we did for the polygon states dataset. This time instead of representing the width of a line, the second, numerical, input represents the radius of the circle that represents each point. When you run your updated code, you should see a number of red circles above each airport in the dataset.

to setup
  ...
  gis:set-drawing-color yellow
  gis:draw states-dataset 1
  gis:set-drawing-color red
  gis:draw airports-dataset 2.5 ;; draw each airport circle with a radius of 2.5
end

Now it is time to create turtles for each of the points in the dataset, and here is where things get a bit tricky if NetLogo is your primary experience programming. In order to create a turtle for each point in the dataset, we need to use a not very commonly used feature of NetLogo, the foreach loop.

The foreach loop runs a section of code over and over again for each item in a list. (If you have programmed before, this is similar to the for item in items: pattern in Python, the for (int num : numbers) pattern in Java or C++, or the foreach (int num in numbers) pattern in C#.) For example, if we had a list of numbers [1 2 3 42 3.14 69105] and wanted to print out each of them times two, we would write:

let numbers [1 2 3 42 3.14 69105]
foreach numbers [number -> show (2 * number)] ;; for each *number* in the list *numbers*, print out (show) *that number* times 2.

In general, the syntax works like this: foreach list-of-items [item -> (do something with that item)]. Remember that a list is different from an agentset in NetLogo. An agentset can only contain agents and it implies that there is no order among the agents within it. This is why we need a foreach loop to create turtles based on each point instead of using ask. ask only works on agentsets and the point features in a shapefile are not agents. A list, on the other hand, can contain any kind of value and implies a specific order, which, not coincidentally, is the order in which the items are "visited" in a foreach loop.

(If you're curious about the mechanics of the -> syntax, then you can read more about anonymous procedures in NetLogo in the NetLogo documentation.)

Now that we know how to use foreach loops, we can go through each point feature in our airports-dataset and create a turtle for each of them.

First, create a breed for airports like so:

breed [airports airport]
airports-own [ abbrev name population infected-population]

Second, go into your setup command, remove the gis:draw command that draws the airports and insert this bare-bones version of our foreach loop:

foreach gis:feature-list-of airports-dataset [ this-airport-vector-feature -> ; foreach VectorFeature in the airports dataset
    let location gis:location-of gis:centroid-of this-airport-vector-feature  ; get the location of its center in NetLogo world-space
    create-airports 1 [                                          ; create an airport
      set xcor item 0 location																	 ; set its xcor to the x coordinate of the airport
      set ycor item 1 location																	 ; set its ycor to the y coordinate of the airport
    ]
  ]

If you run this new code, you should see default NetLogo turtles at the location of each airport. This is probably the most complicated section of code we'll write today, so let's go through it line by line.

foreach gis:feature-list-of airports-dataset [ this-airport-vector-feature ->

If you remember the foreach syntax from before, then this code should seem mostly familiar. The only new addition is the command gis:feature-list-of. This command looks into the input dataset (final_airports_with_populations.shp in this case) and reports back a list of each individual feature within it. In this case, since the airports dataset is a dataset of point features, gis:feature-list-of will report back a list of every point in the dataset, represented by the GIS extension as a "VectorFeature" object. (the VectorFeature object is used regardless of the actual shape type. For example, if we were to run gis:feature-list-of with our states-dataset ,it would report back a list of shape features, but they would also be represented by VectorFeatures.)

let location gis:location-of gis:centroid-of this-airport-vector-feature

In this line, we create a local variable location and store the location of the center of this-airport-vector-feature within it. It may seem counterintuitive that we need two gis: commands to get the location of a point feature, but remember, all our code knows is that this-airport-vector-feature is a VectorFeature, it doesn't know if it is a point, a line, or a polygon, and it doesn't much make sense to try to get the location-of a polygon or a line. Before we can get a location, we need to take the VectorFeature and turn it into a single point, called a "Vertex" in the context of the GIS extension. To do this we use gis:centroid-of which reports the "center of mass" of a given feature, which in the case of a point, is just the center of the point. Once we have that Vertex, we can pass it into gis:location-of and finally get a two-item list representing an xcor and ycor.

create-airports 1 [
	set xcor item 0 location
	set ycor item 1 location
]

These lines simply create an airport and set its xcor and ycor to the xcor and ycor that we just stored in location (the 0th and 1st item of the list, respectively.) At This point, if you run your setup procedure, you should get something like this:

Now that we have the spatial component of our data loaded into NetLogo, it is time to plug in the second piece of our GIS dataset, the feature properties themselves. If you remember, each feature in a GIS dataset is associated with a row of data in the attribute table, and each row has a value for each column, or field of data. For our airports dataset in particular, each airport feature has a "name", an abbreviation ("abbrev"), and a population size that it services ("pop_max"). To get these values out of the dataset and into properties of our airport turtles, we use the gis:property-value command, which takes as input an individual VectorFeature and a string representing the name of the field. So to grab the population, name, and abbreviation from the dataset and add them as properties of our turtle, add these lines of code within the create-airports block:

...
set population gis:property-value this-airport-vector-feature "pop_max"
set abbrev gis:property-value this-airport-vector-feature "abbrev"
set name gis:property-value this-airport-vector-feature "name"
...

To verify that the fields were loaded correctly, you can add a set label abbrev below and, after running this new code, you should see each airport has a label that coo responds to their abbreviation code.

Before we move on to creating the NetLogo model itself, you should take a second to think back on what you have learned and pat yourself on the back. We will touch on exporting data out of NetLogo back into a GIS format at the end, but with just the skills you've learned so far, you can probably already accomplish whatever task led you to learn how to use the GIS NetLogo extension in the first place. You can download data in the proper formats, cleanup and manipulate that data in a GIS program, and import it into a NetLogo world. That's no small feat.

In version 1.3.0 of the GIS extension, a pair of new primitives were added to make the process of creating turtles from point datasets easier. Instead of a foreach loop, you can use create-turtles-from-points or create-turtles-from-points-manual to automatically create one turtle for each point in your points dataset and automatically populate their turtles-own variables with values from their GIS counterpart. In general, the process is to add variables to the breed you want to convert points into with the same names as the property/field names in your dataset and then use create-turtles-from-points with the dataset and the name of the breed to create the turtles. Like the standard create-turtles primitive, you can add in a command block to give the newly created turtles commands upon creation.

To substitute in the create-turtles-from-points primitive into the setup procedure above, you would replace the foreach loop with the following:

;; the variable names must match the field names (ignoring case)
;; to manually map GIS field names to turtle variable names,
;; use gis:create-turtles-from-points-manual
airports-own [abbrev name pop_max] 

gis:create-turtles-from-points airports-dataset "airports" [
	;; the create-turtles-from-points primitive automatically did the 
	;; `set abbrev gis:property-value this-airport-vector-feature "abbrev"` 
	;; step for us, so the value of this turtle's `abbrev` is already 
	;; set correctly.
	set label abbrev
	set shape "circle"
]

Before continuing work on the model, there are a few quick things we need to do as setup. First, add a quick line to set the pcolor of all the patches to a navy within the setup command: ask patches [ set pcolor blue - 3 ] and add a line to the create-airports block that sets each airport's shape to a circle: set shape "circle". Switch back to the interface tab and click on the Settings button and disable vertical and horizontal wrapping (otherwise our planes would be able to magically transport between the Atlantic and Pacific). While you're in the interface tab, add a new Input panel of type "Number" to the interface and name it growth. This value represents the growth rate of the disease within each city. I set mine at 0.1 to start. Finally, add a go forever button to handle running the model.

Here's our code as it stands, just in case you want to catch up:

extensions [gis]

globals [states-dataset airports-dataset]

breed [airports airport]
airports-own [ abbrev name population infected-population]

to load-datasets
  clear-all
  set states-dataset gis:load-dataset "states_continental.shp"
  set airports-dataset gis:load-dataset "airports_with_population_of_nearest_city.shp"
end

to custom-clear
  reset-ticks
	clear-turtles
	clear-patches
	clear-drawing
	clear-all-plots
	;; reset any other global variables you might be using here
end

to setup
  custom-clear
  ask patches [ set pcolor blue - 3 ]

  gis:set-world-envelope gis:envelope-of states-dataset

  gis:set-drawing-color yellow
  gis:draw states-dataset 1

  foreach gis:feature-list-of airports-dataset [ this-airport-vector-feature ->
    let location gis:location-of gis:centroid-of this-airport-vector-feature
    create-airports 1 [
      set xcor item 0 location
      set ycor item 1 location
      set population gis:property-value this-airport-vector-feature "pop_max"
      set abbrev gis:property-value this-airport-vector-feature "abbrev"
      set name gis:property-value this-airport-vector-feature "name"
      set label abbrev
      set shape "circle"
    ]
  ]
end

The next step is to set up the infection modeling within each airport (keeping in mind that we are representing the entire region that an airport is servicing within each airport turtle). Go back to the code and after your set label command within the create-airports block, go ahead and add these lines of code that set an initial infected population and setup some visual characteristics for each airport.

set infected-population 0           ; each city initially has no infected inhabitents
set color black                     ; black represents no infection
set size sqrt (population / 500000) ; this seemed to be a good way to visually scale each city by population without any city being too much larger than any other.

After that, we want to setup the air traffic routes that our airplanes will travel along, and NetLogo links are the perfect choice here. After the close of the foreach block, add in this quick command that will create links between each airport and up to 3 of its surrounding airports.

ask airports [ create-links-with up-to-n-of 3 other airports in-radius 30 ]

Next, we need to start the infection somewhere, so go ahead and randomly pick an unlucky city to have a patient zero.

ask one-of airports [ set infected-population 1 ]

That concludes our setup command, leaving us now ready to create some airplanes.

To start, create a new airplanes breed at the top of your program:

breed [planes plane]
planes-own [ destination infected ]

For our highly simplified model, we pretend like each airplane only carries one passenger and that passenger can either be infected or not infected.

To make our code a bit cleaner, we're going to encapsulate the behavior of hatching planes into a command that our airports will run called hatch-departing-planes. It will create a list of destinations reachable from the current airport (as represented by link-neighbors) and will have a 50% chance of hatching an airplane each tick. If an airplane is hatched, it will have its destination set to one of the possible destinations. It will then face that destination and finally, decide if it has any infected passengers on board. For each plane, there is a x% chance that it will be infected where x = (infected-population / total_population). Throw in a few visual changes, and we're done.

to hatch-departing-planes
  let possible-destinations link-neighbors
  hatch-planes random 2 [ ; random 2 is 0 50% of the time and 1 50 % of the time
    set label ""     ; because `hatch` transfers all properties from the hatching turtle to the hatched turtle, we need to clear out the label
    set size 1
    set color white
    set shape "airplane"
    set infected false ; not-infected is the default
    
    set destination one-of possible-destinations
    face destination      ; face means turn towards a specified turtle
    if random-float 1 < ([infected-population] of myself / [population] of myself) [ ; 'myself' is the departing airport here
      set color red
      set infected true
    ]
  ]
end

Next we want to give our planes a new command, fly so that they can, well, fly. Because we had them face their destination already, they can just fly straight, but we do need them to stop and offload their passenger when they reach their destination. To accomplish this, we have them check if they are within one unit of their destination and, if so, die and have their passenger deplane, possibly bringing their infection with them.

to fly
  if distance destination < 1 [
    if infected [                                                    ; if the plane is infected
      ask destination [                                              ; ask the destination to
        set infected-population infected-population + 1              ; add to the infected population
      ]
    ]
    die                                                              ; in either case, die
  ]
  fd 1                                                               ; if we aren't dead, fly forward one unit
end

Now that we have a behavior for the airports (hatching planes) and a behavior for the planes, we can put them together in a go procedure and check that everything is working as expected.

to go
  ask airports [ hatch-departing-planes ]
  ask planes [ fly ]
  tick
end

If you flip back to the interface tab and run your go procedure, you should see a bunch of white planes flying between different airports. Since there is still only one infected case in the whole model, it is highly highly unlikely that any of the airplanes will carry that infection to another city. However, unfortunately, viral disease is rarely so contained. The next step is to have each city's infected-population rise at a rate proportional to how large the infected population currently is. We're going to encapsulate this behavior inside a command called calculate-infections:

to calculate-infections
  set infected-population infected-population + growth * infected-population * (1 - (infected-population / population))
end

This equation is modeled off a typical logistic model of population growth and its principal characteristic is that it starts growing exponentially when there are many susceptible (not-infected) individuals it can infect but slows down once most of the population already has been infected and stops growth completely once there is no more susceptible population. Most models that model disease also have a concept of individuals that have recovered from the model and are no longer susceptible, but our simulated people are not so lucky and will remain permanently infectious.

In order to visualize this infection growth, let's use scale-color to recolor our airports according to how much of their population is infected. Lets call this command recolor-airport

to recolor-airport
  set color scale-color red (infected-population / population) 0 2
  ; use 0 as the minimum and 2 as the maximum instead of 2 because that way, 100% infected population will be 100% red, whereas if we used 2 as the maximum, 100% infected would be white. 
end

Throw these in to the ask airports block inside the go procedure and we've got our model!

to go
  ask airports [ 
    calculate-infections
    recolor-airport
    hatch-departing-planes
  ]
  ask planes [ fly ]
  tick
end

Here's our final model & code in case you need a reference, as well as a download link to a working version with all associated files:

extensions [gis]

globals [states-dataset airports-dataset]

breed [airports airport]
airports-own [ abbrev name population infected-population]
breed [planes plane]
planes-own [ destination infected ]

to load-datasets
  clear-all
  set states-dataset gis:load-dataset "states_continental.shp"
  set airports-dataset gis:load-dataset "airports_with_population_of_nearest_city.shp"
end

to custom-clear
  reset-ticks
  clear-turtles
  clear-patches
  clear-drawing
  clear-all-plots
  ;; reset any other global variables you might be using here
end

to setup
  custom-clear
  ask patches [ set pcolor blue - 3 ]
  gis:set-world-envelope gis:envelope-of states-dataset

  gis:set-drawing-color yellow
  gis:draw states-dataset 1

  foreach gis:feature-list-of airports-dataset [ this-airport-vector-feature ->
    let location gis:location-of gis:centroid-of this-airport-vector-feature
    create-airports 1 [
      set xcor item 0 location
      set ycor item 1 location
      set population gis:property-value this-airport-vector-feature "pop_max"
      set abbrev gis:property-value this-airport-vector-feature "abbrev"
      set name gis:property-value this-airport-vector-feature "name"
      set label abbrev
      set shape "circle"
      set infected-population 0
      set color black
      set size sqrt (population / 500000)
    ]
  ]

  ask airports [
    create-links-with up-to-n-of 3 other airports in-radius 30
  ]

  ask one-of airports [
    set infected-population 1
  ]
end

to go
  ask airports [ 
    calculate-infections
    recolor-airport
    hatch-departing-planes
  ]
  ask planes [ fly ]
  tick
end

to calculate-infections
  set infected-population infected-population + growth * infected-population * (1 - (infected-population / population))
end

to recolor-airport
  set color scale-color red (infected-population / population) 0 2
end

to hatch-departing-planes
  let possible-destinations link-neighbors
  hatch-planes random 2 [
    set label ""
    set size 1
    set color white
    set shape "airplane"
    set infected false
    
    set destination one-of possible-destinations
    face destination
    if random-float 1 < ([infected-population] of myself / [population] of myself) [
      set color red
      set infected true
    ]
  ]
end

to fly
    if distance destination < 1 [
      if infected [                                                    
        ask destination [                                              
         set infected-population infected-population + 1
       ]
      ]
      die                                                              
    ]
    fd 1                                                               
end