Two Trees
In our previous tutorial, we explored patch-based modeling. In this next step, we will look at modeling where we are simulating individual organsims.
Overview
You are already becoming an expert in Josh! This next example is just using what you already learned but making it a little bigger. Specifically, we are going to run a simulation with lots of agents with each representing a single tree. This can sometimes be useful to include dynamics which are difficult to describe at the patch.
As with the other tutorials, this is a toy model using fake species. The use of specific location is just so that we can use external data.
Data
This time we are going to use data from Cal-Adapt which provides its data as netCDF files. Again, a later tutorial will go into the details of preprocessing for speed. For now, go ahead and grab our preprocessed temperature data and precipitation data. Then, delete the data from the prior tutorial and upload these new files.
Simulation
The simulation here will look pretty similar. Let's go ahead and define our boundaries and time range.
start simulation Main # Specify where the simluation runs grid.size = 200 m grid.low = 36.52 degrees latitude, -118.68 degrees longitude grid.high = 36.42 degrees latitude, -118.45 degrees longitude # Specify when the simulation runs (2024 to 2054) steps.low = 0 count steps.high = 30 count startYear.init = 2024 year.init = startYear year.step = prior.year + 1 # Indicate ecological and spatial limits constraints.maxOccupants = 30 count constraints.minViableTemperature = 273 K constraints.maxViableTemperature = 373 K # Indicate exports exportFiles.patch = "memory://editor/patches" end simulation
As before, we have some constants that we will use across multiple organisms.
Init
Let's start with initializing occupancy of these patches. We want to have this step use the maximum occupancy from the simulation definition.
start patch Default # Randomize starting occupancy countOccupancy.init = sample uniform from 0 count to meta.constraints.maxOccupants initialACount.init = (sample uniform from 0% to 100%) * countOccupancy initialBCount.init = countOccupancy - initialACount # Create initial occupants SpeciesA.init = create initialACount of SpeciesA SpeciesB.init = create initialBCount of SpeciesB end patch
Death
Unlike our ForeverTree, these trees can die. We will set a property on our agents called dead and remove those from the patch that have dead == true.
start patch Default
# ... prior code here ...
# Remove dead occupants at start of simulation step
SpeciesA.start = prior.SpeciesA[prior.SpeciesA.dead == false]
SpeciesB.start = prior.SpeciesB[prior.SpeciesB.dead == false]
speciesACount.start = count(SpeciesA)
speciesBCount.start = count(SpeciesB)
end patch
This looks look mostly similar to what you've seen before but you may notice that there is also a start event handler type in addition to init and step like for the "newCount.end" handler. Each step is actually broken into three phases: start, step, and end. This lets you run multiple calculations on a variable within a single step if needed where you can control the order of those calculations. Also, this is the first time we see a slice. For example, "SpeciesA[prior.SpeciesA.dead == false]" at the start of each step we remove dead trees by filtering for those that aren't dead and saving the resulting distribution "slice" back to the SpeciesA attribute on Default patch.
Extremes
Next, we will also check for extreme temperature events which our organisms will respond to.
start patch Default
# ... prior code here ...
# Determine if in an extreme temperature event
inExtremeTemperature.step = {
const currentTemperature = external temperature
const inExtremeCold = currentTemperature < meta.constraints.minViableTemperature
const inExtremeHot = currentTemperature > meta.constraints.maxViableTemperature
return inExtremeCold or inExtremeHot
}
end patch
Don't worry! You will see the organisms (agents) use this soon.
Birth
Next, we also want to take care of the birth. Let's look at making new trees as well.
start patch Default
# ... prior code here ...
# Calculate num new occupants end step (to have their first step on next step)
totalCount.start = speciesACount + speciesBCount
newCount.end = {
const probabilityNewRaw = (-1 count * (totalCount - 5 count) ^ 2 / 20 count) + 110 %
const probabilityNew = limit probabilityNewRaw to [0%, 100%]
const maxNewCount = (meta.constraints.maxOccupants - totalCount) * probabilityNew
return sample uniform from 0 count to maxNewCount
}
# If empty, split evenly
newSpeciesACount.end = {
if (totalCount > 0) {
return floor((speciesACount / totalCount) * newCount)
} else {
return newCount * 50%
}
}
newSpeciesBCount.end = {
if (totalCount > 0) {
return floor((speciesACount / totalCount) * newCount)
} else {
return newCount * 50%
}
}
# Make the new agents and concatenate to ("|" or add to end of) the current list.
SpeciesA.end = prior.SpeciesA | create newSpeciesACount of SpeciesA
SpeciesB.end = prior.SpeciesB | create newSpeciesBCount of SpeciesB
# Export counts
export.newSpeciesACount.end = newSpeciesACount
export.newSpeciesBCount.end = newSpeciesBCount
end patch
Here, we are getting a closer look at distributions in Josh. These are collections of numbers similar to lists or arrays in other languages but they may also be formal distributions like normal or uniform. We can combine two distributions together using the pipe character ("|") and we can filter the elements within a distribution using square brackets like we saw in our discussion of slices.
Before finishing up this step, note that SpeciesA is a variable or attribute of the Default patch but, by convention, it uses the same name as the organism. You could also call it "SpeciesAPresent" and things would work out.
SpeciesA
That was most of the new stuff! Let's go ahead and use some of our prior tools to define these two tree species.
Each of these will work by looking at by how much growing conditions deviate from species-specific optimal values. These are defined through a quadratic for temperature and sigmoid for precipiation. We use these to determine what percent of optimal growth each agent sees in each time step. This is done by multiplying two "impact" values together as a scale variable.
start organism SpeciesA
# Start a zero height except when if this is the simulation start
dead.init = false
height.init = {
const isFirstYear = meta.startYear == meta.year
const heightSample = sample normal with mean of 20 meters std of 5 meters
const heightSampleLimit = limit heightSample to [0 meters, 40 meters]
return heightSampleLimit if isFirstYear else 0 meters
}
# At each step, determine impact of temperature
temperatureImpact.step = map external temperature
from [270 K, 330 K]
to [0%, 100%]
quadratic
# At each step, determine impact of precipitation
precipitationImpact.step = map external precipitation
from [300 mm, 500 mm]
to [0%, 100%]
sigmoid
# At each step, determine overall growth
newGrowth.step = 1 m * temperatureImpact * precipitationImpact
height.step = prior.height + newGrowth
# At each step end, determine death if over survival probability
dead.end = {
const deadByChance = (sample uniform from 0% to 100%) > 30%
const deadByExtremeTemp = here.inExtremeTemperature
return deadByChance or deadByExtremeTemp
}
end organism
There are two things to highlight in this code. First, we can see now where the variable in the slice used above gets used. Second, like how we did a map with a sigmoid before, we use a quadratic here. By default, quadratic puts the center of the range at the maximum (upside down parabola) but this can be modified by using "quadratic(false)" instead.
SpeciesB
In practice, different species will probably vary quite a bit from each other in their defintion. However, for this toy example, SpeciesB will be about the same as SpeciesA just with some different parameters. Let's take a look!
start organism SpeciesB
# Start a zero height except when if this is the simulation start
dead.init = false
height.init = {
const isFirstYear = meta.startYear == meta.year
const heightSample = sample normal with mean of 20 meters std of 5 meters
const heightSampleLimit = limit heightSample to [0 meters, 40 meters]
return heightSampleLimit if isFirstYear else 0 meters
}
# At each step, determine impact of temperature
temperatureImpact.step = map external temperature
from [290 K, 350 K]
to [0%, 100%]
quadratic
# At each step, determine impact of precipitation
precipitationImpact.step = map external precipitation
from [200 mm, 400 mm]
to [0%, 100%]
sigmoid
# At each step, determine overall growth
newGrowth.step = 1 m * temperatureImpact * precipitationImpact
height.step = prior.height + newGrowth
# At each step end, determine death if over survival probability
dead.end = {
const deadByChance = (sample uniform from 0% to 100%) > 30%
const deadByExtremeTemp = here.inExtremeTemperature
return deadByChance or deadByExtremeTemp
}
end organism
Everything is about the same here except for the parameterization of those important temperature and precipitation impact variables.
Wrap
Let's wrap this up! We just need our units. This time, we will also demonstrate use of unit conversions. Try changing one of the temperature values above to celcius to give it a shot.
start unit K alias Kelvins alias Kelvin end unit start unit C alias Celcius K = 272.15 + current end unit start unit mm alias millimeters alias millimeter m = current / 1000 end unit
Go ahead and give it a run. This one might take a little longer! If you have Josh Cloud, this might be a good time to give it a shot. When ready, let's move forward to a tutorial on data preprocessing.