This notebook tracks progress on the development of GMSE software R package, for game-theoretic management strategy evaluation, and related issues surrounding the development and application of game theory for addressing questions of biodiversity and food security.


Contents:

Project updates



Project updates:

Update: 03 AUG 2021

A new version of GMSE 0.7.0.0 has been released on CRAN now. The changelog provides a summary of the changes and bug fixes. One new feature was written by Adrian Bach, and now GMSE includes the ability to use trajectory-based prediction to determine whether or not intervention on the part of the manager. The new feature will be used in a manuscript in preparation (note that on the GMSE website, I have now included a list of publications that use GMSE). This includes three new arguments:

Regarding the last bullet, I’ve made an update to GitHub v0.7.0.1, which fixes a minor bug that I just noticed this morning that unfortunately made its way into the CRAN version. I don’t think that it will actually affect any simulations though, so I don’t see a need to immediately update on CRAN. See Issue #76 for more details.

There are new options introduced by Jeroen for incrementing user and manager budgets using yield from the landscape.

Improvements to the gmse_gui have also been made by Jeroen, which can be accessed online here. Note that I have not yet included the trajectory strategy options that Adrian developed into the GUI.

After some conversation, I decided to change the movement operations when resources feed more than once (times_feeding > 1) or require a threshold amount of food for survival (consume_surv > 0) or reproduction (consume_repr > 0). Previously, resources would feed, then move times_feeding times. After this moving and feeding, resources would move again at the end of the time step. Somewhat confusingly and arbitrarily, resources would therefore get two moves in a row under these conditions before the timestep ends. In GMSE v0.7, resources now only move once after feeding.

The new GMSE version also fixes some minor bugs, most notably getting rid of the rare cases in which resources were able to feed once on the landscape after dying. There are also new error catches for gmse and gmse_apply functions to make sure that arguments are structured correctly (e.g., quickly flags when a list or vector is used instead of a scalar value).

Update: 21 MAY 2020

Several new features have been added, which will also make their way into the next version of GMSE v0.6.0.x. Two new arguments usr_yld_budget and man_yld_budget are now available in gmse and gmse_apply. These arguments allow an increment of user and manager budget from yield obtained from the landscape. Total yield from one time step of for a user is multiplied by usr_yld_budget, and this product is then added to whatever the user’s baseline budget was. Similarly, the mean total yield of all users is multiplied by man_yld_budget, and this product is then added to whatever the manager’s baseline budget was. There was a bit of an issue when integrating the code with the action threshold budget bonus, so the solution was to add two new columns to the agent array. The budget for any agent is then just the sum of the basline and each of these columns. Hence future development affecting the budget bonus and yield to budget effects can be more easily done independently.

A new vignette has also be written by Adrian, which will be on the updated GMSE website soon.

Update: 18 MAY 2020

A bug has been found and fixed, which was introduced into gmse_apply through the changes to the perceptions of actions for users and managers in commit 86c298eddd1c778defaac04d30ab342eef58fb50. The problem was that the internal elements of the paras vector were changed to reflect columns where perceptions of different actions were to be located in gmse. But in gmse_apply, the same change was not made, so that all of these elements were interpreted as low integers; it’s remarkable that a segfault did not occur in any of the tests. Jeroen pointed out a major difference between the dynamics in runs of gmse and gmse_apply given identical parameter values, which led the investigation into what was causing the problem. Finally, we found the aforementioned problem, which was resolved in commit 7a97eb08171147f688036851be5e3783ab631495. A subsequent fix was made to the testthat function paras vectors where needed. The dynamics in gmse and gmse_apply are once again the same. On Friday, we have tentative plans for a crash party with ConFooBio members to download the newest version of GMSE on the rev branch and try to crash it.

Update: 17 MAY 2020

Over the last two days, I have made some changes to the code and done some additional testing. Changes are summarised below.

Fixed memory leak

Since I had been doing some substantial development in C, I thought it would be wise to run valgrind to check for memory leaks. I found one, due to my forgetting to free some memory in the new landscape.c file. The memory is now correctly freed in the file, and running gmse and gmse_apply in valgrind produces no memory leaks.

==33279== 
==33279== HEAP SUMMARY:
==33279==     in use at exit: 67,085,047 bytes in 12,107 blocks
==33279==   total heap usage: 896,509 allocs, 884,402 frees, 557,718,069 bytes allocated
==33279== 
==33279== LEAK SUMMARY:
==33279==    definitely lost: 0 bytes in 0 blocks
==33279==    indirectly lost: 0 bytes in 0 blocks
==33279==      possibly lost: 0 bytes in 0 blocks
==33279==    still reachable: 67,085,047 bytes in 12,107 blocks
==33279==                       of which reachable via heuristic:
==33279==                         newarray           : 4,264 bytes in 1 blocks
==33279==         suppressed: 0 bytes in 0 blocks
==33279== Reachable blocks (those to which a pointer was found) are not shown.
==33279== To see them, rerun with: --leak-check=full --show-leak-kinds=all
==33279== 
==33279== For counts of detected and suppressed errors, rerun with: -v
==33279== ERROR SUMMARY: 2736 errors from 4 contexts (suppressed: 0 from 0)

I have also initialised a new file in the notebook for functions that we will want to periodically run through valgrind.

Perceptions of actions for users and managers

Following a conversation with Jeroen and Nils, we decided that another change would be beneficial for GMSE v0.6.0.x. In all previous versions of GMSE, how managers and users perceive the consequences of scaring actions is hidden. The consequences were calculated in a reasonable way by GMSE and used within the fitness functions of the genetic algorithm to determine which strategies of setting policy (for manager) or acting (for users) was best. The relevant values for the perceived effects of scaring, culling, castration, feeding, helping offspring, tending crops, and killing crops were set partly in the paras vector, and partly calculated internally in C. Now, instead, they are calculated for each agent (manager and users) and placed directly in the AGENTS array.

The AGENTS array has therefore increased by seven columns. Rows in each of these columns hold the perceived effects of each action for each agent. I’m not sure why I did not think of this before, but there are several benefits to doing it this way. First, anyone running gmse or gmse_apply can now directly alter how users perceive the effects of their own actions using the perceive_scare, perceive_cull, perceive_cast, perceive_feed, perceive_help, perceive_tend, perceive_kill arguments, which set how perceived utility will change given a single instance of scaring, culling, castration, feeding, helping offspring, tending crops, and killing crops, respectively. By default, all of these arguments will be calculated in a realistic way, as has been described already for the developing version, so the arguments can be ignored and GMSE will run as normal. But if, e.g., there is reason to model the users’ perceived efficacy of culling as being much higher than that of scaring (e.g., based on parameterisation from empirical social data), then more appropriate values for a specific system can be used. Additionally, with gmse_apply, the agents array can be changed directly, meaning that manager perception of action effects can be changed, and that perceptions among users can potentially be varied (e.g., if some users are more convinced that scaring will work, while others are more sceptical). I also like that the perceived effects are being returned in the agents array, so the whole parameterisation of the fitness function is much more transparent now.

The code now runs with no memory leaks, and GMSE builds on the rev branch with no notes, warnings, or errors.

Update: 15 MAY 2020

Minor bug fix

Jeroen and I have fixed a minor bug in the GMSE code that was causing a crash given the very restricted situation in which the package is built using the Rstudio drop-down menu. In such builds (but not, curiously, when building from the console), a fatal crash can be caused when looping the gmse_apply function with variable names that are identical to argument names (e.g., stakeholders = 8, followed by gmse_apply(stakeholders = stakeholders)). Just in case, to avoid this issue entirely, I have inserted a short function into gmse_apply that produces an error when attempting to set variable names equal to argument names. This can be completely circumvented in gmse_apply by setting the hidden argument my_way_or_the_highway equal to anything (i.e., not NULL). I have also included a garbage collect at the end of gmse and gmse_apply functions.

Consistency between gmse and gmse_apply

In previous versions of GMSE, the functions gmse and gmse_apply had a few slightly different defaults. This sometimes led to slightly different dynamics when running simulations and comparing between gmse and looped gmse_apply under default arguments. To avoid any cofusion, Jeroen and I have decided to make gmse and gmse_apply default arguments completely identical. The equivalence of the two functions given the same initial conditions has been confirmed through simulations.

Update: 12 MAY 2020

Updated vignettes

I have now checked vignettes 5-7 to see if updates were needed. Only an update to vignette 7 was necessary. Note that even where the code for vignettes remains unchanged, it would still be a good idea to re-knit everything for the website as a test to make sure that everything runs as expected with new defaults for v0.6.0.x.

Uneven land allocation

Jeroen has incorporated some code and a new ownership_var argument to adjust the variation in land ownership among users. I’ve re-worked the way that public land is handled in the new landscape owner building function to account for this, and everything appears to work as intended.

Update: 11 MAY 2020

Bug fix to public land allocation

A small bug was noticed by Jeroen this morning relating to public land allocation. Essentially, when the proportion of public land amounted to fewer cells on the landscape than allocated to a single user, public land would default to zero cells. We needed some way for there to be nonzero public land that was also quite small. Given such situations, a new function small_public_land is called. Public land is placed in the centre of the landscape, thereby cutting into at least two users’ owned cells equally. First a rectangle proportional to the landscape is produced, followed by more public land cells attached to each side of this rectangle as necessary to achieve approximately the correct number of cells. It is important to recognise that simulation parameter values need to be chosen with care. If the total number of landscape cells cannot be evenly distributed among users and public land, then of course deviations from gmse and gmse_apply arguments will occur. The updates cause small deviations to occur for larger sections of landscape when possible. For example, the below produces a landscape with four equal quadrants, three for users and one for public land.

sim <- gmse(time_max = 4, stakeholders = 3, public_land = 0.25,
            land_dim_1 = 100, land_dim_2 = 100, land_ownership = TRUE);

The below is a situation in which the amount of public land is not small. A total of 2000 public land cells are called for, with the remaining 8000 to be divided evenly amongst nine users. This is impossible, and the shortest-splitine algorithm cannot give that level of precision and make the landscapes clean.

sim <- gmse(time_max = 4, stakeholders = 9, public_land = 0.20, 
            land_dim_1 = 100, land_dim_2 = 100, land_ownership = TRUE);

Instead, 1840 cells are allocated to public land, with each user being allocated rectangular blocks of 900 or 920 landscape cells. But when public land to be allocated is very low, fewer cells than given to a single user, deviations are then placed on users to make the public land proportion as accurate as possible in the centre of the landscape.

sim <- gmse(time_max = 4, stakeholders = 7, public_land = 0.03, 
            land_dim_1 = 100, land_dim_2 = 100, land_ownership = TRUE);

In the case of the above, e.g., exactly 300 public land cells are put down in the centre of the landscape, while users are given somewhere between 1303 and 1462 cells.

Update gmse GUI

The gmse_gui function has now been updated so that all options are available (and default landscape allocation is enacted) for running simulations from within the shiny GUI. I have not put this up on the shiny server yet, but on the rev branch, the new GUI can be used within R. This now concludes all of the planned updates for code changes to v0.6.0.x, and now the last two things to do will be to add vignettes and update the website. Please contact me to request any additional changes.

Update: 10 MAY 2020

Yesterday and today I implemented the fourth bullet point listed on 6 MAY 2020, which changes the default pattern of land ownership in GMSE. In over a decade of coding individual-based models, I think that this might be the first time that I have needed to apply recursion to a problem. I actually only just learned that recusion is possible in R, but decided to write the function in C just to get a bit more speed.

The goal of the new algorithm is to divide a rectangular landscape of length land_dim_1 and width land_dim_2 into stakeholder smaller rectangles that are as equal as possible in area (i.e., cell number). To do this, I used a shortest-splitline algorithm in a new file landscape.c. The relevant function is the break_land function. For a given rectangular landscape, the break_land function draws a line splitting the landscape along its longest dimension (i.e., if the rectangle has a length of 20 and a width of 10, it is split lengthwise). If the number of owners to be placed within the landscape is even, then the line splits exactly halfway down the middle of the landscaope. If the number of owners to be placed within the landscape is odd, then the line is placed in proportion to the ratio of ownership – in other words, if 3 owners need to be placed on the landscape, then the split line would be placed 2/3 of the way along the longest dimension (3/5 of the way for 5 owners, 4/7 of the way for 7 owners, and so forth). The rectangular landscape is then split into two smaller sub-rectangles. The same splitting procedure is then run on these sub-rectangles with the new number of owners (i.e., half of the original in the case of an even split, and appropriate portions such as 1 and 2, or 2 and 3, for the odd number of original owners, with the biggest rectangle getting the higher number of owners). The splitting of the landscape into smaller and smaller sub-rectangles continues until the number of owners that need to be placed on a rectangle equals one, in which case an integer representing an agent’s ID is placed on all landscape cells in the rectangle.

In the case of public land, the correct proportion public_land is put on the landscape by running the above algorithm in the break_land function with an appropriately added number of owners. If, for example, half of the land is meant to be public, then the total number of unique sub-rectangles to be made is adjusted to users / (1 - proportion_public_land). Where if proportion_public_land = 0.5, e.g., the total number of sub-rectangles to be made will be twice the number of users. Once all of the cells are assigned an integer value by the algorithm, all cells with a number higher than any user ID are then re-assigned as public land.

These changes have now been pushed to the rev branch, so as of writing this, they are available to use. Relevant testthat functions have been adjusted where necessary to account for the code changes, and the new code passes all checks and all CRAN tests with no notes, warnings, or errors.

I like the above algorithm so much that I do not see any real need to even make the old straight-line simple sort approach an option in gmse. The old approach was a simplification that always seemed a bit unrealistic, and I always wanted to come up with a better way of allocating owned landscape cells. Please contact me or raise an issue on GitHub if you really want the old way back in gmse (note that gmse_apply of course always allows custom landscape definition).

With this new feature, I have now tackled bullet points 1, 3, and 4 from 6 MAY 2020. The next thing that I will do will be to update the GUI, which I do not expect to take long. Once this is finished, I will begin writing my vignette on resource density, and updating the website accordingly. I think I might also update the vignette on default GMSE structures to briefly discuss the above algorithm and how landscape ownership is assigned. I will put out another announcement before sending to CRAN, but please let me know if there are changes that you would like to see in GMSE v0.6.0.x.

Update: 8 MAY 2020

I have made some progress on the plans from 6 MAY 2020. None of my updates are incorporated into the main R branch, so they are not permanent and I am still happy to take suggestions for improvements or other ideas in addition or in contrast to what folows.

I have completed the first bullet point of the plans from 6 MAY 2020. The order of operations is therefore changed in the resource model on the work-in-progress rev branch. Hence, calls to the default resource model proceed by first having resources consume on the landscape, then birth, death, and movement. Within resource consumption rules, feeding happens first, followed by movement for up to times_feed iterations. The landscape consumption also now calls either one function or another: one function where consume_surv > 0 or consume_repr > 0 in which consumption happens in a random order among resources for a variabl enumber of times_feed, and a second (faster) function where consumption happens once for all resources (in order of the array, as the order of consumption does not affect the resources themselves).

Despite the information in the vignette on GMSE structures, column 13 of the agents array was actually being used within the more or less unused anectdotal function, in which each agent ‘looks around’ from the cell on which they are standing to see how many resources are in the general area. This function has no purpose at the moment (I was thinking about using it to model disagreement with the manager’s estimate and the users’ inferences about population size), but writes to the column anywa in gmse (but not gmse_apply). I’m going to refer to it as written to, but unused in the vignette, and use column 14 of the agent array to record the number of cells owned.

Writing to the number of cells owned is done in a new count_agent_cells in R or the count_owned_cells function in C. I wrote the R function first, but then had second thoughts about relying on it in gmse and particularly gmse_apply. What if we later want land ownership to change over the course of a simulation? I’m not sure in this case, actually, how much extra speed we get from C (it relies on two nested for loops across the landscape, rather than a which or == in R, which I suspect is actually just as fast since the underlying C is probably better than mine). In any case, it seems to make sense to have the owned cells calculated at the beginning of the user function, so this is what I ultimately decided on. But the R function can also be used.

I have changed the fitness function for users as such:

These changes are now pushed to the rev branch, which is still under ongoing development. This effectively takes care of the third bullet point from 6 MAY 2020. Hence, bullets one and three are now done, so here is what remains.

I will try to make the above changes as soon as possible, but the changes that will have the biggest effect on simulations are now finished, so people can play around with the new options. The update in progress builds with no warnings, errors, or notes, and passes all tests in testthat.

Update: 6 MAY 2020

New version of GMSE

A new version of GMSE v0.5.0.0 is now available on GitHub. Here are the major updates:

The new version can be downloaded from GitHub, but will not be submitted to CRAN. It is likely to be the last version of GMSE that is completely backwards compatible with previous versions. In other words, if you set parameter values for previous versions of GMSE as arguments in gmse or gmse_apply and ran simulations, then you are effectively running the same code (and should get the same results) in v0.5.0.0. In the next version update, this will not necessarily be true because minor changes to the default order of operations are planned (see below). In practice, this is highly unlikely to change any results from simulations in previous versions, but it is good to be aware of what has changed if you are expecting old code to run identically on new versions of GMSE.

Plans for GMSE v0.6.0.x: what will change

While backwards compatibility has been a priority up until now with GMSE development, the most recent developments in GMSE v0.5.0.0 make it clear that the benefits of make some small changes to the default simulations outweigh the costs. This mostly has to do with the order of operations in the default resource submodel (currently: movement, feeding, birth, death). Given the new option of resource consumption affecting birth and death rates, and the consequences of farmer actions (i.e., when land_ownership == TRUE), a change in the order makes sense. When feeding affects birth and death in v0.5.0.0, resources first move once, then move again, then eat, move, eat, move, etc. for as many times_feeding as defined. Hence, resources are moving twice at the start of a time step, and potentially landing on good landscape cells, but then not eating at all. With respect to farmer actions, if resources move far in a single time step, then their efforts might have very little actual effect on their individual yield. This is because, since movement happens first, resources might be highly likely to move to another farmer’s land anyway, making the actions futile. This might not be a problem for modelling some systems (where such actions are in fact futile), but the fitness function of the genetic algorithm is designed in such a way as to have farmers assume that one scare or cull will decrease the effects of a resource accordingly. In practice, this doesn’t matter too much because the absolute values in the fitness functions are arbitrary for now (i.e., we could multiply them by the probability that a resource did not leave the farmer’s landscape and expect the same farmer actions – there’s nothing else they can do). But it could become an issue later, and we might at some point be really interested in looking at patterns among different users in a more interesting way. Therefore, in GMSE v0.6.0.x, I plan to change the order of operations in the default resource submodel (some more of my reasoning is discussed in Issue 56). Overall, I plan for the following changes:

The end result of all of this will be, I think, an improved GMSE, but also one that breaks from previous versions in the default options and potentially stakeholder behaviour. For example, land owning users will inevitably see culling as a slightly better option for increasing yield, all else being equal (i.e., given identical costs, as set by the manager). I’m keen to start this new version as soon as possible. The new version will be pushed to CRAN after these changes have been implemented. I will attempt this by the end of the month, or early June at the latest. Please contact me either via email or GitHub if you have any thoughts or concerns about the next version.

Update: 23 APR 2020

With commit 6d1e5e6ca41bc6bb45c62cb174f34480c3c5add1, I believe that I now have resolved Issue #56 and created a version of GMSE that now allows resource death and reproduction to be tied to resource consumption on the landscape. Initial testing shows that this works in both gmse and gmse_apply. I have not had a chance to code in the suggestion from Jeroen to have it be the probability of survival or reproduction that changes with increasing consumption. Rather, given that a resource has consumed a sufficient quantity of resources consume_surv, the resource will survive. A separate new gmse argument consume_repr ties reproduction to landscape consumption by having resources produce floor(consumed / consume_repr) offspring given some total consumption consumed.

There is still considerable stochasticity introduced in terms of whether or not a resource will survive or reproduce. First, each resource can potentially feed more than once by using the new times_feeding argument. When times_feeding > 1 (default is 1), resources will consume the landscape, then move, then consume, and so forth until all of their times_feeding have run out. In practice, all resources get the same times_feeding value, but the code for executing it depends on a new column in the resources array. This means that there can be variation in the number of times a resource feeds if using gmse_apply (or gmse, if we later decide to sample feeding rates from a distribution). Feeding events also occur in a random order, which avoids the first resource on the list filling themselves up before moving on to the next resource (as would be a natural way to program it with a for loop). I did it this (slightly less efficient) way to more realistically model how resources would really move around simultaneously on a landscape eating as they move from place to place.

Hence, in the end, there is a lot of uncertainty surrounding whether an individual resource will have enough to survive and reproduce in a time step. But the probabilistic idea of Jeroen could also be imposed with a function of some sort that would link consumption to a probability of death (will need to think of the best way to actually do this in the code though). Testing shows that resource increases to a carrying capacity imposed by the constraints of the landscape, with a lot more stochasticity than occurs just by settting a value for global res_death_K. To get started, for example, the below demonstrates the new features.

sim <- gmse(consume_surv = 2, consume_repr = 3, times_feeding = 6, 
            res_birth_type = 0, res_death_type = 0, land_ownership = TRUE, 
            stakeholders = 12, scaring = TRUE, time_max = 40, user_budget = 1);

I set the user_budget to unity above just to show how the population rises to a natural carrying capacity when the manager and users are unable to cause any population change. We can recover a (usually) well-regulate population by settting the user_budget back to the default.

sim <- gmse(consume_surv = 2, consume_repr = 3, times_feeding = 6, 
            res_birth_type = 0, res_death_type = 0, land_ownership = TRUE, 
            stakeholders = 12, scaring = TRUE, time_max = 40, 
            user_budget = 1000);

Just to clarify, each landscape cell still puts out a value of 1. The argument consume_surv = 2 means that resources need to consume two cells worth to survive a time step (by default, resources consume half of what’s on a landscape cell with each visit). The argument consume_repr = 3 means that resources need to consume three cells worth to produce one offspring. The argument times_feeding = 6 means that resources get to feed on six cells in each time step.

Setting res_birth_type = 0 and res_death_type = 0 effectively removes all other factors in GMSE that affect resource birth and death. Hence, what is observed from running the above code is purely a consequence of the new rules limiting survival and reproduction due to landscape consumption. Note that these values need to be specified, else GMSE will assume that other defaults apply in addition to the effects of landscape consumption. In other words, the global external carrying capacity from res_death_K will still limit population growth because by default res_death_type = 2, and the mean birth number lambda will still apply because res_birth_type = 2. I don’t want to override this just in case anyone wants to include multiple types of processes affecting birth and death (I also want the new version to be backwards compatible with existing code), but I think it will be important for me to write some documentation that makes using these new features really clear. I’ll start a new vignette.

The new version is GMSE v0.5.0.0 with Adrian now included in the list of contributors. The action_thres and budget_bonus can now be run in gmse and gmse_apply, and usr_budget_rng and age_repr are also available options. The new version passes all CRAN checks with no notes, warnings, or errors. The new version should be completely backwards compatible. If you have existing code to run, it should do the same thing in the new version that it did in previous versions. Everything is now just a new feature that can be used, but is not by default.

The last thing that I think I want to do before merging with master, updating the website, and sending to CRAN, is to make a new default landscape. It would be good to make each landscape into more farm-like blocks rather than the narrow stripes that exist now. I don’t have a good algorithm for doing this for an arbitrary number of stakeholders (i.e., for N stakeholders, separate an X by Y landscape into roughly equal blocks of sizes ca (X*Y / N)). If anyone wants to have a go at this and create a function that does it, please let me know! Else I’ll try to get around to it before early next month where the plan is to send to CRAN. If anyone has time (Nils Jeroen Adrian), feel free to download from the dev branch and go.

install.packages("devtools");
library(devtools);
install_github("ConFooBio/GMSE", ref = "dev")

It should be up now.

Update: 13 APR 2020

I have now completed two major tasks in a GMSE update, first mentioned earlier in the month. The new code is up on the rev branch.

I have resolved Issue #53, so partial matching is not allowed in gmse_apply (or in the tests). I have also resolved Issue #59, meaning that it is now possible to use different action thresholds and budget bonuses in both gmse and gmse_apply (arguments action_thres and budget_bonus), and variable user budgets can now be specified in gmse, and in gmse_apply (argument usr_budget_rng – of course, manual edits to user budgets were always allowed in gmse_apply by directly editing the agents array).

In making these changes, I have also allowed gmse_apply to keep track of the current time step when old calls to the function are passed back through old_list (i.e., when looping old list as in here, the time step will also be tracked automatically through the loops, not just outside of gmse_apply through the loop index – this makes the action threshold possible, as it only applies after the first time step).

I now need to move on to Issue #56, which will more or less complete the new version, though I also might update the default landscape ownership placement.

New tests have been written in testthat, and the package builds with now errors, warnings, or notes.

Update: 1 APR 2020

I am now beginning a moderate update to GMSE, which I hope to have complete and up on CRAN by the end of the month. There are three specific goals of this update:

  1. Resolve Issue #53: Just to head off any potential issues in gmse_apply, several places in the code use the format sim_list$element_name rather than `sim_list[[“element_name”]]. The latter is more robust because it prohibits partial matching. Hence in gmse_apply.R, this should be changed.

  2. Resolve Issue #59: The code from AdrianBach needs to be incorporated into the master branch of GMSE. Doing this will mainly consist of resolving a merge conflict and a bit of refactoring, and I plan to do this early on in developing what I hope to be a new version of GMSE by the end of APR 2020.

  3. Resolve Issue #56: This is the major one. I want to add options in GMSE to allow resource dynamics to be dependent up on landscape properties. Specifically, for resources to consume the landscape cell contents, and for options to exist that allow resource birth, death, and (potentially) movement to be dependent on consumption. This will mainly involve updates to the resource model, with updates to the gmse and gmse_apply function as needed.

  4. In addressing 1-3, I will need to check to make sure that gmse_apply is handling all of the new options correctly, and new tests will need to be written as necessary.

I will continue to update as need be and make sure that any updates from other team members (e.g., Issue #58) are merged appropriately before the new GMSE version.

Update: 4 FEB 2020

Yesterday, I received an email noting a problem with the GMSE documentation. Specifically, a warning arose because ... was in the documentation but is no longer allowed as an argument in the default R functions. Notice the issue on this line of the manager.R file, which is contradicted later by the function itself. This has now been resolved.

Update: 15 MAR 2019

GMSE v0.4.0.11 is now available on GitHub. Some minor changes and bug fixes are included in the update version.

I will update once more when GMSE v0.4.0.11 is successfully submitted to CRAN.

Update: 21 AUG 2018

GMSE v0.4.0.7 is now available on GitHub, and the official GMSE repository is now transferred to the ConFooBio organisation. I am currently in the process of fixing some website issues, but most of the important stuff has transferred to the new website location.

Update: 16 MAY 2018

GMSE v0.4.0.3 is now available on CRAN. A new website for GMSE has also been launched. This website was built with the R package pkgdown, recently released on CRAN. The site contains all of the vignettes and documentation for GMSE, and also includes a link to this lab notebook. A submission of the accompanying manuscript will soon be uploaded on bioRxiv.

Update: 14 MAY 2018

A new GMSE v0.4.0.3 has now been pushed to the master branch on GitHub and has been submitted to CRAN. The biggest update in this new version is a series of vignettes, plus a minor improvement to the genetic algorithm. More updates will follow soon, including some re-organisation of the GMSE project and a new manuscript submission.

Update: 13 APR 2018

I have re-worked the way that a manager restimates how the change in their policy affects users’ actions. The new new_act function in the genetic algorithm (games.c) performs well for getting more precise cost settings. The former way of doing it was much more of a blunt instrument, and it had a ceiling issue – that is, the manager would believe that higher costs caused fewer actions even when the resulting cost was over the users’ budgets.

/* =============================================================================
 * This function updates an action based on the change in costs & paras
 *     old_cost: The old cost of an action, as calculated in policy_to_counts
 *     new_cost: The new cost of an action, as calculated in policy_to_counts
 *     paras: Vector of global parameters
 * ========================================================================== */
int new_act(double old_cost, double new_cost, double old_act, double *paras){
    
    int total_acts;
    double users, max_to_spend, acts_per_user, cost_per_user, total_cost;
    double pr_on_act, budget_for_act, mgr_budget, min_cost;
    
    users        = paras[54] - 1; /* Minus one for the manager */
    min_cost     = paras[96];     /* Minimum cost of an action */
    max_to_spend = paras[97];     /* Maximum per user budget   */
    mgr_budget   = paras[105];    /* Manager's total budget    */
    
    total_cost    = 0.0;
    if(old_cost < mgr_budget){
        total_cost    = old_act * old_cost; /* Total cost devoted to action */
    }

    cost_per_user = (total_cost / users); /* Cost devoted per user */
    pr_on_act     = cost_per_user / max_to_spend;    /* Pr. devoted to action */

    /* Assume that the proportion of the budget a user spends will not change */
    budget_for_act = max_to_spend * pr_on_act;

    /* Calculate how many actions to expect given acts per user and users */
    acts_per_user  = budget_for_act / (new_cost + min_cost);  
    total_acts     = (double) users * acts_per_user;

    return(total_acts);
}

This new way of assessing how users will act is now the function to be run in the background of all manager genetic agorithms. Very nicely, this also resolves an annoyance with the maximum allowed budgets. Previously, it was unclear why maximum budgets greater than 10000 were causing problems (managers were making bad predictions). I have now set the maximum budget to an order of magnitude higher, and there are no longer any apparent issues. A new version of GMSE will soon have this update.

Update: 21 DEC 2017

New Issue #40: Age distribution bump

Running simulations using gmse_apply, jeremycusack noticed a small but noticeable sharp decline in the population size at a generation equal to the maximum age of resources in the population (used a maximum age of 20). This decline is caused by the initial seed of resources having a uniform age distribution. In the first generation, these resources reproduce offspring that all have an age of zero, leading to an age structure in the population with many zero age individuals and a uniform distribution of ages greater than zero. The initial seed of individuals with random ages died gradually, but there were enough individuals in the initial offspring cohort that made it to the maximum age for it to have a noticeable effect in decreasing population size (i.e., all of these resources died on the maximum_age + 1 time step).

This effect can be avoided entirely given sufficient burn in generations of a model, and is less of a problem when the maximum age is low because this allows the age distribution to stabilise sooner. Further, using gmse_apply can avoid the issue by directly manipulating resources ages after the initial generation. Nevertheless, it would be useful to have a different default of age distributions instead of a uniform distribution.

One way to do this would be to find the age (\(A\)) at which a resource is expected to be alive with a probability of \(0.5\), after accounting for mortality (\(\mu\)). This is simply calculated below:

\((1 - \mu)^A = 0.5\)

The above can be re-arranged to find A,

\(A = \frac{log(0.5)}{log(1 - \mu)}\).

Note that we could use a switch function (or something like it in R) to make \(A = 0\) when \(\mu = 1\), and revert to a uniform distribution of \(\mu = 0\) (though this should rarely happen).

The value of \(\mu\) would depend on res_death_type, and be processed in make_resource, which is used in both gmse and gmse_apply. If res_death_type = 1 (density independent, rarely used), then \mu is simply equal to remov_pr. If res_death_type = 2 (density dependent), then \mu could be found perhaps using something like the following:

mu = (RESOURCE_ini * lambda) / (RESOURCE_ini + RESOURCE_ini * lambda) gi This would get a value that is at least proportional to expected mortality rate of a resource (if res_death_type = 3, then we could use the some of types 1 and 2). Overall, the documentation should perhaps recommend finding a stable set of age distributions for a particular set of parameter combinations when using gmse_appy (i.e., through simulation), then using that distribution as an initial condition. But something like the above could probably get close to whatever the stable age distribution would be, at least close enough to make the decline in population size trivial.

I will start to consider some of the above as a potential default for the next version of GMSE. The best way to do this is probably to look at how code from the res_remove function in the resource.c file can best be integrated into a function called by the R function make_resource (i.e., either use the equations, or estimates of them, or somehow call res_remove directly).

Update: 13 DEC 2017

Improved convergence criteria

I have introduced and then immediately resolved Issue #39.

The convergence criteria has now been fixed with commit f598d8e52b47ef2017cac13d09aac1fb7aa6b506. To do this, I re-configured some of the genetic algorithm code into easier to read functions for checking the fitness increase. Now two separate ways of checking the increase in fitness from one genetic algorithm generation to the next exist; one for managers and one for users. This is needed because user fitness values are greater than zero and increase as their utility is maximised, but manager fitness values are less than zero and increase toward zero as their utility is maximised. The genetic algorithm now checks for a percentage improvement in fitness.

Now the default value of converge_crit equals 1, which means it does actually play a role sometimes (or is expected to). The genetic algorithm will continue until the percent increase in fitness from the previous generation is less than one percent. In practice, this doesn’t noticeably affect much, but it does allow better strategies to be found more quickly, and without having to play with ga_mingen to find them under extreme parameter settings (e.g., huge budgets and rapid shifts in abundance).

The new fix has now been checked and built with Winbuilder into v0.3.2.0, but I am leaving this on the development branch for now in anticipation of other potential improvements to be made soon.

Update: 1 NOV 2017

CRAN ready GMSE v0.3.1.7 – more flexibility, better error messages

I have now completed some substantial coding of error messages, which will be called in both gmse and gmse_apply. Essentially, these provide some help to software users who parameterise their models in a way that does not work with GMSE. For example, if the parameter stakeholders is set to equal a negative number, an error message will be returned that informs the user that at least one stakeholder is required in the model. These error messages become a bit more important in gmse_apply, where it is possible for users to include arguments that don’t make sense (e.g., arrays of incorrect dimensions, or arguments that contradict one another).

The function gmse_apply has also been improved to make looping it easier. What had been happening during testing was that we were finding it all too easy to crash R by reading in parameters that contradicted one another (e.g., changing setting the landscape dimensions through land_dim_1 and land_dim_2 caused a crash when also trying to add in a LAND of different dimension – now this returns an error that LAND and land_dim_1 disagree about landscape size). This has been resolved in two ways. First, I have included many error messages meant to catch bad and contradictory arguments in gmse_apply (and, to a lesser extent gmse); it is still possible to crash R by setting things incorrectly, but you have to work very hard to do it – i.e., it almost has to be deliberate, as far as I can tell. Second, I have added the argument old_list to gmse_apply, which is FALSE by default, but can instead take the output of a previous full list return of gmse_apply (where get_res = Full). An element of the full list includes the basic output from which key parameters can be pulled. As a reminder, the basic gmse_apply output looks like the below.

$resource_results
[1] 1062

$observation_results
[1] 680.2721

$manager_results
         resource_type scaring culling castration feeding help_offspring
policy_1             1      NA     110         NA      NA             NA

$user_results
        resource_type scaring culling castration feeding help_offspring tend_crops kill_crops
Manager             1      NA       0         NA      NA             NA         NA         NA
user_1              1      NA       9         NA      NA             NA         NA         NA
user_2              1      NA       9         NA      NA             NA         NA         NA
user_3              1      NA       9         NA      NA             NA         NA         NA
user_4              1      NA       9         NA      NA             NA         NA         NA

An example gmse_apply used in a loop is below.

to_scare <- FALSE;
sim_old  <- gmse_apply(scaring = to_scare, get_res = "Full", stakeholders = 6);
sim_sum  <- matrix(data = NA, nrow = 20, ncol = 7);
for(time_step in 1:20){
    sim_new               <- gmse_apply(scaring = to_scare, get_res = "Full", 
                                        old_list = sim_old);
    sim_sum[time_step, 1] <- time_step;
    sim_sum[time_step, 2] <- sim_new$basic_output$resource_results[1];
    sim_sum[time_step, 3] <- sim_new$basic_output$observation_results[1];
    sim_sum[time_step, 4] <- sim_new$basic_output$manager_results[2];
    sim_sum[time_step, 5] <- sim_new$basic_output$manager_results[3];
    sim_sum[time_step, 6] <- sum(sim_new$basic_output$user_results[,2]); 
    sim_sum[time_step, 7] <- sum(sim_new$basic_output$user_results[,3]); 
    sim_old               <- sim_new;
    print(time_step);
}
colnames(sim_sum) <- c("Time", "Pop_size", "Pop_est", "Scare_cost", 
                       "Cull_cost", "Scare_count", "Cull_count");

The ouput sim_sum is shown below.

      Time Pop_size   Pop_est Scare_cost Cull_cost Scare_count Cull_count
 [1,]    1      733  839.0023         NA       110          NA         54
 [2,]    2      768  702.9478         NA       110          NA         54
 [3,]    3      824  725.6236         NA       110          NA         54
 [4,]    4      933  907.0295         NA       110          NA         54
 [5,]    5     1180  816.3265         NA       110          NA         54
 [6,]    6     1345 1224.4898         NA        10          NA        426
 [7,]    7     1114 1269.8413         NA        10          NA        425
 [8,]    8      820  884.3537         NA       110          NA         54
 [9,]    9      952  793.6508         NA       110          NA         54
[10,]   10     1101  884.3537         NA       110          NA         54
[11,]   11     1299 1111.1111         NA        12          NA        402
[12,]   12     1079  907.0295         NA       110          NA         54
[13,]   13     1227 1564.6259         NA        10          NA        431
[14,]   14      934  839.0023         NA       110          NA         54
[15,]   15     1065 1133.7868         NA        10          NA        423
[16,]   16      768  725.6236         NA       110          NA         54
[17,]   17      869  929.7052         NA       110          NA         54
[18,]   18      949  907.0295         NA       110          NA         54
[19,]   19     1049  884.3537         NA       110          NA         54
[20,]   20     1200 1020.4082         NA        64          NA         90

We can take advantage of gmse_apply to dynamically change parameter values mid-loop. For example, below shows the same code, but with a policy of scaring introduced on time step 10.

to_scare <- FALSE;
sim_old  <- gmse_apply(scaring = to_scare, get_res = "Full", stakeholders = 6);
sim_sum  <- matrix(data = NA, nrow = 20, ncol = 7);
for(time_step in 1:20){
    sim_new               <- gmse_apply(scaring = to_scare, get_res = "Full", 
                                        old_list = sim_old);
    sim_sum[time_step, 1] <- time_step;
    sim_sum[time_step, 2] <- sim_new$basic_output$resource_results[1];
    sim_sum[time_step, 3] <- sim_new$basic_output$observation_results[1];
    sim_sum[time_step, 4] <- sim_new$basic_output$manager_results[2];
    sim_sum[time_step, 5] <- sim_new$basic_output$manager_results[3];
    sim_sum[time_step, 6] <- sum(sim_new$basic_output$user_results[,2]); 
    sim_sum[time_step, 7] <- sum(sim_new$basic_output$user_results[,3]); 
    sim_old               <- sim_new;
    if(time_step == 10){
        to_scare <- TRUE;
    }
    print(time_step);
}
colnames(sim_sum) <- c("Time", "Pop_size", "Pop_est", "Scare_cost", 
                       "Cull_cost", "Scare_count", "Cull_count");

The above simulation results in the following output for sim_sum.

      Time Pop_size   Pop_est Scare_cost Cull_cost Scare_count Cull_count
 [1,]    1      745  657.5964         NA       110          NA         54
 [2,]    2      805 1111.1111         NA        12          NA        400
 [3,]    3      473  634.9206         NA       110          NA         54
 [4,]    4      504  566.8934         NA       110          NA         54
 [5,]    5      577  498.8662         NA       110          NA         54
 [6,]    6      600  430.8390         NA       110          NA         54
 [7,]    7      648  612.2449         NA       110          NA         54
 [8,]    8      714  702.9478         NA       110          NA         54
 [9,]    9      813  612.2449         NA       110          NA         54
[10,]   10      914 1020.4082         NA        64          NA         90
[11,]   11     1011 1179.1383         57        10          49        301
[12,]   12      858  725.6236         10       110         193         37
[13,]   13     1011 1043.0839         37        30           0        198
[14,]   14      989 1043.0839         57        30           0        198
[15,]   15      983 1065.7596         48        20          10        270
[16,]   16      851  839.0023         10       110         193         37
[17,]   17      962 1111.1111         38        12          58        306
[18,]   18      783  612.2449         10       110         193         37
[19,]   19      862  816.3265         10       110         193         37
[20,]   20      963  702.9478         10       110         182         38

Hence, in addition to all of the other benefits of gmse_apply, one new feature is that we can use it to study change in policy availability – in this case, what happens when scaring is introduced as a possible policy option. Similar things can be done, for example, to see how manager or user power changes over time. In the example below, users’ budgets increase by 100 every time step, with the manager’s budget remaining the same. The consequence appears to be decreased population stability and a higher likelihood of extinction.

ub          <- 500;
sim_old     <- gmse_apply(get_res = "Full", stakeholders = 6, user_budget = ub);
sim_sum     <- matrix(data = NA, nrow = 20, ncol = 6);
for(time_step in 1:20){
    sim_new               <- gmse_apply(get_res = "Full", old_list = sim_old,
                                        user_budget = ub);
    sim_sum[time_step, 1] <- time_step;
    sim_sum[time_step, 2] <- sim_new$basic_output$resource_results[1];
    sim_sum[time_step, 3] <- sim_new$basic_output$observation_results[1];
    sim_sum[time_step, 4] <- sim_new$basic_output$manager_results[3];
    sim_sum[time_step, 5] <- sum(sim_new$basic_output$user_results[,3]);
    sim_sum[time_step, 6] <- ub;
    sim_old               <- sim_new;
    ub                    <- ub + 100;
    print(time_step);
}
colnames(sim_sum) <- c("Time", "Pop_size", "Pop_est", "Cull_cost", "Cull_count",
                       "User_budget");

The output of sim_sum is below.

      Time Pop_size   Pop_est Cull_cost Cull_count User_budget
 [1,]    1     1215 1405.8957        10        292         500
 [2,]    2     1065 1224.4898        10        336         600
 [3,]    3      833  680.2721       110         36         700
 [4,]    4      936  907.0295       110         42         800
 [5,]    5     1174 1224.4898        10        401         900
 [6,]    6      887  521.5420       110         54        1000
 [7,]    7      988  680.2721       110         60        1100
 [8,]    8     1084  975.0567       110         60        1200
 [9,]    9     1208  861.6780       110         66        1300
[10,]   10     1360 1133.7868        10        520        1400
[11,]   11      975  861.6780       110         78        1500
[12,]   12     1079 1156.4626        10        560        1600
[13,]   13      597  770.9751       110         90        1700
[14,]   14      595  476.1905       110         96        1800
[15,]   15      586  612.2449       110        102        1900
[16,]   16      584  770.9751       110        108        2000
[17,]   17      557  589.5692       110        114        2100
[18,]   18      519  521.5420       110        120        2200
[19,]   19      469  521.5420       110        120        2300
[20,]   20      430  453.5147       110        126        2400

There is an important note to make about changing arguments to gmse_apply when old_list is being used: The function gmse_apply is trying to avoid a crash, so the function will accomodate parameter changes by rebuilding data structures if necessary. For example, if you change the number of stakeholders (and by including an argument stakeholders to gmse_apply, it is assumed that stakeholders are changing even they are not), then a new array of agents will need to be built. If you change landscape dimensions (or just include the argument land_dim_1 or land_dim_2), then a new landscape willl be built. This is mentioned in the documentation.

GMSE v0.3.3.7 passes all CRAN checks in Rstudio. I will make sure that the code works with win-builder, then prepare the new submission. Alternatively, as always, the newest GMSE version can be downloaded through GitHub if you have devtools installed in R.

devtools::install_github("bradduthie/GMSE")

I will soon update the manuscript for GMSE and upload it to biorXiv.

Update: 23 OCT 2017

Bug fix concerning density-based estimation

An error with density-based resource estimation (observe_type = 0) at very high values of agent_view was identified by Jeremy. When managers had a view of the landscape that encompassed a number of cells that was calculated to be larger than the actual number of landscape cells (as defined by land_dim_1 * land_dim_2), the manager would understimate actual population size. This occurred only in the manager.c file and not in the equivalent R function shown during plotting. The bug was fixed in commit a916b8f8a40041b5f08984cf73348108482dde59 with a simple if statement. This has therefore been resolved in a patched GMSE v0.3.1.3, which is now availabe on GitHub.

Update: 19 OCT 2017

Bug fix concerning resource movement

An error with the res_move_obs parameter was identified by Jeremy. This parameter was supposed to only affect resource movement during observation, but an if statement corrected in commit 5eeb88d285af57984171e7d72410659b3b441af3 was causing res_move_obs = FALSE to stop moving entirely in the resource model. This has now been resolved in a patched GMSE v0.3.1.1, which is now available on GitHub.

New option for removal of resources

A new option has been included for the argument res_death_type. By setting res_death_type = 3 in gmse or gmse_apply, resources can experience both density dependent (caused by res_death_K) and density independent (caused by remove_pr) removal simultaneously. Effects of each are independent of one another (i.e., both processes occur simultaneously, so the calculation of population size affecting removal due to carrying capacity includes resources that might experience density independent mortality).

Update: 16 OCT 2017

New group_think parameter in GMSE v0.3.1.0

A new group_think parameter has been developed by Jeremy and me, and included into an updated v0.3.1.0. This parameter is defined as FALSE by default, but when set to be TRUE will cause all users to act as a single block instead of independently. In the code, what happens is that a single user (user ID number 2) runs through the genetic algorithm, but then instead of having the resulting actions apply to only this user, they apply to all users so that the genetic algorithm only needs to be run once in the user model. This decreases simulation time, particularly when there are a lot of users to model, but at a cost of removing all variation in actions among users. The group_think parameter can be defined in both gmse() and gmse_apply(), but I have not added it as an option in gmse_gui().

Update: 13 OCT 2017

GMSE v0.3.0.0 now available with gmse_apply

The gmse_apply function is now available on a new GMSE version 0.3.0.0. (minor tweaks to other functions have also been made, but nothing that changes the user experience of gmse – mostly typos corrected in the documentation). The new function allows software users to integrate their own submodels (resource, observation, manager, and user) into GMSE, or to use their own submodels entirely within a single function.


GMSE apply function

The gmse_apply function is a flexible function that allows for user-defined sub-functions calling resource, observation, manager, and user models. Where such models are not specified, GMSE submodels ‘resource,’ ‘observation,’ ‘manager,’ and ‘user’ are run by default. Any type of sub-model (e.g., numerical, individual-based) is permitted as long as the input and output are appropriately specified. Only one time step is simulated per call to gmse_apply, so the function must be looped for simulation over time. Where model parameters are needed but not specified, defaults from gmse are used.

gmse_apply arguments


Example uses of gmse_apply

A simple run of gmse_apply() will return one generation of gmse using default submodels and parameter values.

sim <- gmse_apply();

For sim, the default ‘basic’ results are returned as below.

$resource_results
[1] 1102

$observation_results
[1] 1179.138

$manager_results
       scaring culling castration feeding help_offspring
policy      NA      10         NA      NA             NA

$user_results
        resource_type scaring culling castration feeding help_offspring tend_crops kill_crops
Manager             1      NA       0         NA      NA             NA         NA         NA
user_2              1      NA      70         NA      NA             NA         NA         NA
user_3              1      NA      75         NA      NA             NA         NA         NA
user_4              1      NA      69         NA      NA             NA         NA         NA
user_5              1      NA      74         NA      NA             NA         NA         NA

Note in the case above we have the total abundance of resources returned, the estimate of resource abundance from the observation function, the costs the manager sets for the only available action of culling, and the number of culls attempted by each user.

The above was produced by all of the individual-based functions that are default in GMSE; custom generated subfunctions can instead be included provided that they fit the specifications described above. For example, we can define a very simple logistic growth function to send to res_mod instead.

alt_res <- function(X, K = 2000, rate = 1){
    X_1 <- X + rate*X*(1 - X/K);
    return(X_1);
}

The above function takes in a population size of X and returns a value X_1 based on the population intrinsic growth rate rate and carrying capacity K. Iterating the logistic growth model by itself under default parameter values with a starting population of 100 will cause the population to increase to carrying capacity in roughly 7 generations. The function can be substituted into gmse_apply to use it instead of the default GMSE resource model.

sim <- gmse_apply(res_mod = alt_res, X = 100, rate = 0.3);

The gmse_apply function will find the parameters it needs to run the alt_res function in place of the default resource function, either by running the default function values (e.g., K = 2000) or values specified directly into gmse_apply (e.g., X = 100 and rate = 0.3). If an argument to a custom function is required but not provided either as a default or specified in gmse_apply, then an error will be returned.

To integrate across different types of submodels, gmse_apply translates between vectors and arrays between each submodel. For example, because the default GMSE observation model requires a resource array with particular requirements for column identites, when a resource model subfunction returns a vector, or a list with a named element ‘resource_vector,’ this vector is translated into an array that can be used by the observation model. Specifically, each element of the vector identifies the abundance of a resource type (and hence will usually be just a single value denoting abundance of the only focal population). If this is all the information provided, then a resource_array will be made with default GMSE parameter values with an identical number of rows to the abundance value (floored if the value is a non-integer; non-default values can also be put into this transformation from vector to array if they are specified in gmse_apply, e.g., through an argument such as lambda = 0.8). Similarly, a resource_array is also translated into a vector after the default individual-based resource model is run, should the observation model require simple abundances instead of an array. The same is true of observation_vector and observation_array objects returned by observation models, of manager_vector and manager_array (i.e., COST) objects returned by manager models, and of user_vector and user_array (i.e., ACTION) objects returned by user models. At each step, a translation between the two is made, with necessary adjustments that can be tweaked through arguments to gmse_apply when needed. Alternative observation, manager, and user, submodels, for example, are defined below; note that each requires a vector from the preceding model.

# Alternative observation submodel
alt_obs <- function(resource_vector){ 
    X_obs <- resource_vector - 0.1 * resource_vector;
    return(X_obs);
}

# Alternative manager submodel
alt_man <- function(observation_vector){
    policy <- observation_vector - 1000;
    if(policy < 0){
        policy <- 0;
    }
    return(policy);
}

# Alternative user submodel
alt_usr <- function(manager_vector){
    harvest <- manager_vector + manager_vector * 0.1;
    return(harvest);
}

All of these submodels are completely deterministic, so when run with the same parameter combinations, they produce replicable outputs.

gmse_apply(res_mod = alt_res, obs_mod = alt_obs, 
           man_mod = alt_man, use_mod = alt_usr, X = 1000);

The above, for example, produces the following output (Note that the X argument needs to be specified, but the rest of the subfunctions take vectors that gmse_apply recognises will become available after a previous submodel is run).

$resource_results
[1] 1500

$observation_results
[1] 1350

$manager_results
[1] 350

$user_results
[1] 385

Note that the manager_results and user_results are ambiguous here, and can be interpreted as desired – e.g., as total allowable catch and catches made, or as something like costs of catching set by the manager and effort to catching made by the user. Hence while manger output is set in terms of costs of performing each action, and user output is set in terms of action attempts, this need not be the case when using gmse_apply (though it should be recognised when using default GMSE manager and user functions).

GMSE default submodels can be added in at any point.

gmse_apply(res_mod = alt_res, obs_mod = observation, 
           man_mod = alt_man, use_mod = alt_usr, X = 1000)

The above produces the results below.

$resource_results
[1] 1500

$observation_results
[1] 1655.329

$manager_results
[1] 655.3288

$user_results
[1] 720.8617

If we wanted to, for example, specify a simple resource and observation model, but then take advantage of the genetic algorithm to predict policy decisions and user actions, we could use the default GMSE manager and user functions (written below explicitly, though this is not necessary).

gmse_apply(res_mod = alt_res, obs_mod = alt_obs, 
           man_mod = manager, use_mod = user, X = 1000)

The above produces the output below returning culling costs and culling actions attempted by four users (note that the default manager target abundance is 1000).

$resource_results
[1] 1500

$observation_results
[1] 1350

$manager_results
       scaring culling castration feeding help_offspring
policy      NA      10         NA      NA             NA

$user_results
        resource_type scaring culling castration feeding help_offspring tend_crops kill_crops
Manager             1      NA       0         NA      NA             NA         NA         NA
user_2              1      NA      70         NA      NA             NA         NA         NA
user_3              1      NA      70         NA      NA             NA         NA         NA
user_4              1      NA      71         NA      NA             NA         NA         NA
user_5              1      NA      73         NA      NA             NA         NA         NA

Instead of using the gmse function, we might simulate multiple generations by calling gmse_apply through a loop, reassigning outputs where necessary for the next generation (where outputs are not reassigned, new defaults will be inserted in their place, so, e.g., if we were to just loop without reassigning any variables, nothing would update and we would be running the same model, effectively, multiple times). Below shows how this might be done.

sim1      <- gmse_apply(get_res = "full", lambda = 0.3);
RESOURCES <- sim1$resource_array;
LAND      <- sim1$LAND;
PARAS     <- sim1$PARAS;
results   <- matrix(dat = NA, nrow = 40, ncol = 4);

for(time_step in 1:40){
    sim_new <- gmse_apply(RESOURCES = RESOURCES, LAND = LAND, PARAS = PARAS,
                          COST = COST, ACTION = ACTION, stakeholders = 10,
                          get_res   = "full", agent_view = 20);
    
    results[time_step, 1] <- sim_new$resource_vector;
    results[time_step, 2] <- sim_new$observation_vector;
    results[time_step, 3] <- sim_new$manager_vector;
    results[time_step, 4] <- sim_new$user_vector;
    
    RESOURCES <- sim_new$resource_array;
    LAND      <- sim_new$LAND;
    PARAS     <- sim_new$PARAS;
    COST      <- sim_new$COST;
    ACTION    <- sim_new$ACTION;
}

colnames(results) <- c("Abundance", "Estimate", "Cull_cost", "Cull_attempts");

The above results in the following output for results.

      Abundance  Estimate Cull_cost Cull_attempts
 [1,]      1195 1165.9726        10   716
 [2,]      1045  939.9167       110   461
 [3,]      1160 1160.0238        10   715
 [4,]      1056 1183.8192        10   715
 [5,]      1014  850.6841       110   468
 [6,]      1171 1237.3587        10   717
 [7,]      1026  993.4563       110   464
 [8,]      1202  957.7632       110   464
 [9,]      1394 1469.3635        10   702
[10,]      1333 1457.4658        10   702
[11,]      1277 1397.9774        10   702
[12,]      1175 1415.8239        10   702
[13,]      1088  701.9631       110   468
[14,]      1275 1207.6145        10   718
[15,]      1200 1332.5402        10   718
[16,]      1116 1029.1493        45   512
[17,]      1249 1814.3962        10   699
[18,]      1141 1273.0518        10   722
[19,]      1019  963.7121       110   455
[20,]      1216 1629.9822        10   708
[21,]      1088 1130.2796        10   708
[22,]       988 1035.0982        38   537
[23,]      1056 1029.1493        45   505
[24,]      1154  749.5538       110   463
[25,]      1344 1499.1077        10   722
[26,]      1268 1386.0797        10   712
[27,]      1165 1493.1588        10   707
[28,]      1061 1070.7912        19   633
[29,]      1019 1076.7400        17   663
[30,]       961  600.8328       110   457
[31,]      1135  874.4795       110   450
[32,]      1338 1189.7680        10   701
[33,]      1275 1600.2380        10   710
[34,]      1174 1362.2844        10   709
[35,]      1104 1112.4331        12   685
[36,]      1003 1302.7960        10   715
[37,]       828 1183.8192        10   712
[38,]       649  785.2469       110   462
[39,]       739 1023.2005        56   488
[40,]       813  910.1725       110   455

Note that managers increase the cost of culling based on the time step’s estimated abundance, and user culling attempts decrease when culling costs increase.

In addition to the flexibility of allowing user-defined submodels, gmse_apply is also useful for modellers who might be interested in simulating processes not currently available in gmse by itself. For example, if we wanted to model a sudden environmental perturbation decreasing population size, or a sudden influx of new users, after 30 generations, we could do so in the loop.

In the near future, the gmse_apply function will be included in the GMSE vignette and submitted to CRAN with the rest of v0.3.0.0 – in the mean time, I believe that all major bugs have been ironed out, but please let me know or report an issue if you are able to crash the function (i.e., if you run it and it causes R to crash – you should always get an error message before this happens).

To download the latest GMSE v0.3.0.0, simply run the below in R (make sure that devtools is installed).

devtools::install_github("bradduthie/GMSE")

I welcome any feedback, and I expect to submit an update to CRAN around late October.

Update: 12 OCT 2017

New function gmse_apply complete and tested

I have now completed the gmse_apply function, which exploits the full modularity of GMSE by allowing software users to develop their own sub-functions and string them together with any combination of GMSE default sub-functions. As a brief summary, gmse_apply includes the following features:

Any arguments for custom user functions can simply be passed along by specifying them in gmse_apply. For example, if we have a custom resource function alt_res below:

alt_res <- function(X = 1000, K = 2000, r = 1){
    X_1 <- X + r*X*(1 - X/K);
    return(X_1);
}

We can simply include the above in gmse_apply as follows to use the very simple logistic growth sub-model with the individual-based submodels that are defaults of GMSE.

sim_app <- gmse_apply(res_mod = alt_res);

The gmse_apply function simply adds in GMSE defaults for unspecified models, but we can specify them too.

sim_app <- gmse_apply(res_mod = alt_res, obs_mod = observation);

To adjust parameters in the alternative resource model, simply add in the arguments as below.

sim_app <- gmse_apply(res_mod = alt_res, X = 2000, K = 5000, r = 1.2);

The gmse_apply function will know where to place them, and update them should they be needed for other models.

I will give a more lengthy description of how to use gmse_apply tomorrow, when I push GMSE v0.3.0.0 to the master branch of GitHub and advertise the update.

Update: 6 OCT 2017

Compensation suggestion

A suggestion from Jeremy to include a compensation option for users. Users could devote some of their budget to compensation, then managers could compensate a proportion of their damaged yield. Implementing this will require consideration from the manager’s perspective with respect to the genetic algorithm – the users’ perspective will be easier because a user can remember their previous losses and assess compensation versus culling. Managers might have to think about how compensation could incentivise non-culling, but this might actually already work given the way the manager anticipates actions; more investigation into this will be useful following the finalisation of gmse_apply(), which is in progress.

Update: 28 SEP 2017

Progress has been made on the gmse_apply() function. My goal is to make this as modular as possible – to allow any four functions to be included in the GMSE framework, including arbitrary arguments to each function. The gmse_apply() function will recognise which arguments go along with which functions, and naturally string together results from one sub-function to the input of the next sub-function (though this will demand that the output from functions is labelled in a way that matches the arguments of the next function; e.g., if you have a ‘N_total’ as input for the observation model, then ‘N_total’ will either need to be labelled output of the resource model or specificied directly in gmse_apply()). Default submodels will be the IBMs used in gmse(), and where arguments are not specified by the software user in gmse_apply() (e.g., LAND) they will be built from default gmse() parameters.

Update: 27 SEP 2017

The GMSE GUI has been updated with all of the new features in version 0.2.3.0. The gmse_gui() function is likewise updated in a new patch version 0.2.3.1. I did this quickly because the GUI was actually easy to update; plans for the gmse_apply function are now also clear, and I hope to have a working function and version 0.3.0.0 by the end of the week, or by early next week.

Update: 26 SEP 2017


GMSE Version 0.2.3.0 on GitHub

I have pushed a new version 0.2.3.0 of GMSE onto the master branch of GitHub, which means that the most up-to-date version can be installed using the code below (make sure the devtools library is installed).

devtools::install_github("bradduthie/GMSE")

The new version includes multiple new features:

To run a simple default simulation, the gmse function remains unchanged.

sim <- gmse();

To plot the effort of managers and users, use the below.

plot_gmse_effort(agents = sim$agents, paras = sim$paras, 
                 ACTION = sim$action,  COST = sim$cost);

Below summarises the results more cleanly, extracting key information from sim.

gmse_summary(sim);

And as before, the GUI can be called directly from the R console.

gmse_gui();

The GUI does not yet allow you to get a vew of the plot_gmse_effort output, or a gmse_summary, but this will be a goal for future versions of GMSE.

If able, I recommend updating to version 0.2.3.0 as soon as possible. In the coming few days, I will also add the gmse_apply function, primarily for developers who will benefit from a more modular way of using GMSE, allowing for different types of submodules to be used within the broader GMSE framework. When the new apply function has been added (and possibly the GUI improved), I will submit a new version 0.3.x.x to CRAN.


Bug Fix and tweaks to agent prediction

I have now fixed a bug in the code that was causing confusion between culling and castration. After recompiling and running simulations, manager and user actions improve. I have also made some minor changes to default gmse() options. Regarding the predicted consequences of manager and user actions (i.e., the predictions from the agents’ perspective that guid their decision making), I have adjusted some things to make them more in line with what is expected in the simulation as follows (recall that managers are interested in global abundance and users are interested specifically in how abundance affects themselves):

  1. Scaring: Managers predict no change in resource abundance, while users predict a decrease of 1
  2. Culling: Managers and users predict a decrease of \(1 + \lambda\) (Note that this brings in knowledge of birth rate a priori – might want to allow for a change in this in the simulation, but it also seems realistic for agents to recognise that adults can reproduce, and a value is needed to reflect this)
  3. Castration: Managers and users predict a decrease of \(\lambda\)
  4. Feeding: Managers and users predict an increase of \(\lambda\)
  5. Help offspring: Managers and users predict an increase of 1

These values are a bit more in line with what will actually happen, so we assume that managers and users are a bit more informed now. It also allows for a bit more differentiation among actions. Overall, the model appears to perform better now – meaning that managers and users appear to be better predictors of the conseuqneces of their actions.

Before finishing the gmse_apply() function, I will push an updated version of GMSE to GitHub with these changes, plus new plotting options.

Update: 25 SEP 2017

I have written a gmse_summary function (see below), which returns a simplified list that includes four elements, each of which is a table of data: 1. resources, a table showing time step in the first column, followed by resource abundance in the second column. 2. observations, a table showing time step in the first column, followed by the estimate of population size (produced by the manager) in the second column. 3. costs, a table showing time step in the first column, manager number in the second column (should always be zero), followed by the costs of each action set by the manager (policy); the far-right column indicates budget that is unused and therefore not allocated to any policy. 4. actions, a table showing time step in the first column, user number in the second column, followed by the actions of each user in the time step; additional columns indicate unused actions, crop yield on the user’s land (if applicable), and the number of resources that a user successfully harvests (i.e., ‘culls’).

At the moment, I have not added in the actual number of resources that a user culls. This will be added shortly, after which I will post a new function. Doing so is a bit more complicated because it requires me to go into the C code and make a recording every time it happens (see how I plan to do this below the function).

gmse_summary <- function(gmse_results){
    time_steps <- dim(gmse_results$paras)[1];
    parameters <- gmse_results$paras[1,];
    #--- First get the resource abundances
    res_types    <- unique(gmse_results$resource[[1]][,2]);
    resources    <- matrix(dat  = 0, nrow = time_steps, 
                           ncol = length(res_types) + 1);
    res_colna    <- rep(x = NA, times = dim(resources)[2]);
    res_colna[1] <- "time_step";
    for(i in 1:length(res_types)){
        res_colna[i+1] <- paste("type_", res_types[i], sep = "");
    }
    colnames(resources) <- res_colna;
    #--- Next get estimates abd the costs set by the manager
    observations    <- matrix(dat  = 0, nrow = time_steps, 
                              ncol = length(res_types) + 1);
    costs   <- matrix(dat = NA, nrow = time_steps*length(res_types), ncol = 10);
    agents  <- gmse_results$agents[[1]];
    users   <- agents[agents[,2] > 0, 1];
    actions <- matrix(dat  = NA, ncol = 13,
                      nrow = time_steps * length(res_types) * length(users));
    c_row  <- 1;
    a_row  <- 1;
    for(i in 1:time_steps){
        the_res            <- gmse_results$resource[[i]][,2];
        manager_acts       <- gmse_results$action[[i]][,,1];
        resources[i, 1]    <- i;
        observations[i, 1] <- i;
        land_prod          <- gmse_results$land[[i]][,,2];
        land_own           <- gmse_results$land[[i]][,,3];
        for(j in 1:length(res_types)){
            #---- Resource abundance below
            resources[i,j+1] <- sum(the_res == res_types[j]);
            #---- Manager estimates below
            target_row <- which(manager_acts[,1] == -2 & 
                                    manager_acts[,2] == res_types[j]);
            estim_row  <- which(manager_acts[,1] ==  1 & 
                                    manager_acts[,2] == res_types[j]);
            target <- manager_acts[target_row, 5];
            adjusr <- manager_acts[estim_row,  5];
            observations[i,j+1] <- target - adjusr;
            #---- Cost setting below
            costs[c_row, 1]  <- i;
            costs[c_row, 2]  <- res_types[j];
            estim_row    <- which(manager_acts[,1] ==  1 & 
                                  manager_acts[,2] == res_types[j]);
            if(parameters[89] == TRUE){
                costs[c_row, 3] <- manager_acts[estim_row,  8];
            }
            if(parameters[90] == TRUE){
                costs[c_row, 4] <- manager_acts[estim_row,  9];
            }
            if(parameters[91] == TRUE){
                costs[c_row, 5] <- manager_acts[estim_row,  10];
            }
            if(parameters[92] == TRUE){
                costs[c_row, 6] <- manager_acts[estim_row,  11];
            }
            if(parameters[93] == TRUE){
                costs[c_row, 7] <- manager_acts[estim_row,  12];
            }
            if(parameters[94] == TRUE){
                costs[c_row, 8] <- parameters[97];
            }
            if(parameters[95] == TRUE){
                costs[c_row, 9] <- parameters[97];
            }
            costs[c_row, 10] <- manager_acts[estim_row, 13] - parameters[97];
            c_row <- c_row + 1;
            #--- Action setting below
            for(k in 1:length(users)){
                usr_acts <- gmse_results$action[[i]][,,users[k]];
                actions[a_row, 1] <- i;
                actions[a_row, 2] <- users[k];
                actions[a_row, 3] <- res_types[j];
                res_row <- which(usr_acts[,1] == -2 & 
                                     usr_acts[,2] == res_types[j]);
                if(parameters[89] == TRUE){
                    actions[a_row, 4] <- usr_acts[res_row,  8];
                }
                if(parameters[90] == TRUE){
                    actions[a_row, 5] <- usr_acts[res_row,  9];
                }
                if(parameters[91] == TRUE){
                    actions[a_row, 6] <- usr_acts[res_row,  10];
                }
                if(parameters[92] == TRUE){
                    actions[a_row, 7] <- usr_acts[res_row,  11];
                }
                if(parameters[93] == TRUE){
                    actions[a_row, 8] <- usr_acts[res_row,  12];
                }
                if(j == length(res_types)){
                    if(parameters[104] > 0){
                        land_row <- which(usr_acts[,1] == -1);
                        if(parameters[95] > 0){
                            actions[a_row, 9]  <- usr_acts[land_row, 10];
                        }
                        if(parameters[94] > 0){
                            actions[a_row, 10] <- usr_acts[land_row, 11];
                        }
                    }
                    actions[a_row, 11] <- sum(usr_acts[, 13]);
                }
                if(parameters[104] > 0){
                    max_yield <- sum(land_own == users[k]);
                    usr_yield <- sum(land_prod[land_own == users[k]]);
                    actions[a_row, 12] <- 100 * (usr_yield / max_yield);
                }
                a_row <- a_row + 1;
            }
        }
    }
    cost_col <- c("time_step", "resource_type", "scaring", "culling",
                  "castration", "feeding", "helping", "tend_crop", 
                  "kill_crop", "unused");
    colnames(costs)        <- cost_col;
    colnames(resources)    <- res_colna;
    colnames(observations) <- res_colna;
    action_col <- c("time_step", "user_ID", "resource_type", "scaring", 
                    "culling", "castration", "feeding", "helping", "tend_crop",
                    "kill_crop", "unused", "crop_yield", "harvested");
    colnames(actions) <- action_col;
    the_summary <- list(resources    = resources, 
                        observations = observations, 
                        costs        = costs, 
                        actions      = actions);
    return(the_summary);
}

To record kills, I think that the best way is to use the resource mortality adjustment column (at the moment, column 17 in C and 18 in R of the resource array). Mortality as of now is just adjusted to 1 in the event of a kill, and mortality occurs whenever a random probability is greater than or equal to 1. Hence, I can replace the 1 value with the user’s ID (for non-managers, this must be at least 1), and then the resource array will record the ID of the user that killed it at the particular time step. Note that this cannot be done for other adjustments such as growth rate or offspring production because the values are not interpreted as probabilities.

I will do the above tomorrow, which should not take too long. I will then continue work on the gmse_apply function.

Update: 22 SEP 2017

Currently, the gmse() function returns a list that includes all of the data produced by the model, some details of which are required for plotting.

sim_results <- list(resource    = RESOURCE_REC,
                    observation = OBSERVATION_REC,
                    paras       = PARAS_REC,
                    land        = LANDSCAPE_REC,
                    time_taken  = total_time,
                    agents      = AGENT_REC,
                    cost        = COST_REC,
                    action      = ACTION_REC
                );

I think that this list is fine, perhaps necessary, to keep, but the ConFooBio group has also concluded that there should be some easier to understand summary of the data. I propose that some function written, a gmse_summary(), that summarises the results in an easier to understand way would be useful. The function could just be run as below.

sim         <- gmse();
sim_summary <- gmse_summary(sim);

The output of gmse_summary() should be a list of all of the relevant information that a user might want to plot or analyse. It should include the following list elements.

  1. sim_summary$resources
  2. sim_summary$observations
  3. sim_summary$costs
  4. sim_summary$actions

More might be needed, but the above should be a good starting point that will provide four clear data tables for the user. The tables will look like the below.

1. Resource abundances over time

time_step abundance
1 100
2 104
99 116
100 108

In the above, only the resource abundance is reported to the software user, though it might also be useful to have additional columns as well eventually.

2. Observation estimates of abundance over time

time_step estimated_abundance
1 102
2 101
99 121
100 112

In the above, only the estimate from the observaiton submodel is reported to the software user. Additional columns might also be useful for things like confidence intervals, though for now I’m not sure if this is needed.

3. Costs set in each time step

time_step manager scaring castration culling feeding helping unused
1 0 40 NA 60 NA NA 0
2 0 36 NA 62 NA NA 2
99 0 0 NA 100 NA NA 0
100 0 3 NA 97 NA NA 0

In the above, the manager number is always 0 because this is the number of the agent that has that role in GMSE. All impossible actions (specificed by the simulation) are labelled NA, while the possible scaring and culling actions are given values that correspond to the cost of each action for users in each time step. Hence the table summarises policy for each time step in a way that software users can interpret more cleanly.

4. Actions in each time step

time_step user scaring castration culling feeding helping tend_crop kill_crop unused crop_yield harvested
1 1 50 NA 50 NA NA NA NA 0 90 12
1 2 59 NA 40 NA NA NA NA 1 92 9
1 3 100 NA 0 NA NA NA NA 0 89 0
2 1 44 NA 66 NA NA NA NA 0 88 16
2 2 52 NA 48 NA NA NA NA 0 94 12
2 3 98 NA 0 NA NA NA NA 2 90 0
99 1 36 NA 63 NA NA NA NA 1 79 20
99 2 40 NA 60 NA NA NA NA 0 83 18
99 3 28 NA 72 NA NA NA NA 0 88 12
100 1 35 NA 62 NA NA NA NA 3 82 18
100 2 37 NA 63 NA NA NA NA 0 84 22
100 3 23 NA 77 NA NA NA NA 0 84 13

The above action table has more rows in it than the cost table because a row is needed for each user in each time step. This gives the software user full access to each individual user’s actions, and their results. Note that as above, castration, feeding, and helping, are not options. Additionally, in this hypothetical simulation, tending or killing crops are not options, so no actions are performed. Users divide their budget between scaring and culling in each time step. The last two columns also give useful information to the software user. The first is crop yield on the user’s owned land (should probably be NA if land_ownership = FALSE), which will reflect the percentage of the total possible yield (or maybe raw yield?) for each user – hencing allowing the table to direclty correlate actions with yield. The last column is the number of resources ‘harvested,’ which I think should count successful ‘kills’ (rather than just actions devoted to culling). The realised culling might be lower than the actions devoted to culling, for example, if not enough resources are actually on the user’s land to cull. Additional statistics for each user could be added in as columns, but this seems a good place to start. This gmse_summary producing a list of the above four tables will be included in the next version of GMSE, along with the new plotting function highlighting the conflict itself, and the gmse_apply function discussed on 6 SEPT.

Update: 15 SEP 2017

Continued progress has been made on slides for an upcoming talk.

Update: 13 SEP 2017

I will be giving a talk on 19 September 2017 for the Mathematics and Statistics Group at the University of Stirling on GMSE as a general tool for management strategy evaluation. Slides for this talk will be available on GitHub.

Update: 8 SEP 2017

The alternative approach from Wednesday is being implemented smoothly. Passing user-defined functions in a modular way is possible, but inputs and outputs need to be carefully considered within gmse_apply(). The objective is to make things as easy and flexible as possible for the user, while also making sure that the function runs efficiently.

Update: 6 SEP 2017

A modular function for modellers

I am beginning work on a gmse_apply() function, which will improve the modularity of GMSE for developers. The goal behind this function will be to provide a simple tool for allowing developers to use their own resource and observation models and, with the correct inputs, take advantage of the manager and user functions. Hence, simple resource and observation models will be possible, but the flexibility of GMSE should be retained as much as possible. A few starting points include the following:

Inputs and outputs of different functions will then include the following:

As an alternative, at least to the implementation, I think that the call could be made at the level of the individual resource() and observation() R functions. This was kind of always the plan, but there’s a semi-dirty way to mix numerical resource and observation models with the full individual-based manager and user models. This can be done by adding a model option to be user defined through an if( is.function(model) == TRUE ) in the resource.R function. If the condition is satisfied, then resource() will shift to the user generated model. This can actually be done for all of the submodels very easily.

This alternative might be a better way to go. The aforementioned ‘dirty’ part of the technique might be to check to see if the output is in the correct form, then, if only a vector is returned – turn it into the correct form by making a data frame that has the same number of rows by calling make_resource. The type 1 values could correspond to vector elements. Admittedly, this could get slow for huge population sizes, but population sizes would have to be massive for R to slow down from simplying making a matrix with a lot of rows. In any case, it would at least standardise the input and output for the user of gmse_apply in a way that plays nice with everything else in GMSE.

Similarly, the observation function could also call make_resource if a vector is returned (since individual variation wouldn’t be relevant in the numerical model).

With this alternative approach, no changes to the C code need to be made – the inputs and outputs just need to be tweaked into a standardised way when a vector or scalar is returned from any user-defined model (small detail – population size needs to be an integer). This can be an option later for the user and manager models – though I’m not sure how this would work, exactly. A benefit here is that some parts of the model could concievably individual-based, with others being numerical – the trade-off being the requirement for discrete resource numbers and a very small amount of slowdown (which will almost certainly not be noticable for any resaonable model).

The gmse_apply function would then initialise a very small number of agents, and a small landscape (unless otherwise specified) in every run. The possibility of passing more options could be applied with a simple .... This would also require a sub-function build_para_vec, which would be used for the sole purpose of taking the list of options included (same as in gmse()) and passing it to the sub-function, with any functions not passed being assumed defaults (and most would be irrelevant). So the default function should then look like

gmse_apply(resource_function = "IBM", observation_function = "IBM", manager_function = "IBM", user_function = "IBM", res_number = 100, ...);

I think at least an initial population needs to be specified, but everything else can be left up to the user, with the elipses passing to the function building the parameter vector (which can also be called by gmse(), replacing some clutter). Overall, the function will run without any input if none is specified, defaulting to an IBM with a population size of 100 for one generation. All other options, including non-standard functions, are left to the user.

Additional thoughts

Working this through, I’m slightly apprehensive about the motivation for including gmse_apply function. Once you strip the mechanistic approach from the resource and observation models, all you really have are two values: (1) the population abundance or density and (2) the estimate of the population size or density. Once you include the manage_target into gmse_apply (necessary, I believe), then the genetic algorithm is really just a fancy way of getting the difference between the population estimate and the target size, and then setting a number of culling actions acceptable for users. Users then cull as much as possible because they’re assumed to want to use the resource as much as possible. Of course, we can consider other parameters that affect user actions (e.g., maximum catch, effort), but if we’re interested in learning about how these concepts affect harvesting in theory, then they can and should probably be studied using a simpler model. The real point of the genetic algorithm is that it allows for complex, multi-variable goal driven behaviour, as might occur given indirect effects (e.g., organisms on crop yield) or multiple options (e.g., culling versus scaring or growing) and spatial complexities. There seems little to be gained by calling the genetic algorithm to tell users to cull as much as possible, which can be done with a (very) simple function.

Update: 28 AUG 2017

I have finally fixed the annoyance in the shiny app of GMSE that caused the bottom of the browser to black, hence making it difficult to set parameter values in some tabs.

Additionally, by hovering over the different options in the application, the software user can now see a brief description of what each option does in the simulation.

Update: 23 AUG 2017

I am experimenting with ways of demonstrating the conflict between what a manager incentivises, and what the users actually do, in GMSE. Below are some plots that show this for a few sample simulations. The five panels in each plot correspond to the five possible actions where policy is set. Policy set by the manager is shown with the black solid line, with the thin coloured lines reflecting individual user effort expended into each action.

The right axis is fairly easy to interpret – it’s just the percentage of the user’s total budget devoted to a particular action (note, this is not necessarily the number of actions a user performs because different actions can cost different amounts – hence the term ‘effort’).

The left axis is a bit trickier – it’s how permissive of an action the manager is in practice. High values correspond to an action being highly permitted by the manager (i.e., the manager invests no effort in making these actions costly), whereas low values correspond to an action being less permitted (i.e., the manger invests highly in making these actions costly for users).

The end result is that the lines indicating manger permissiveness are typically correlated with user effort towards any particular action. In the first example below, this is true for scaring and culling (as the manager becomes more permissive of these actions, users tend to take advantage and spend more effort doing them). Note that users do not feed because they have nothing to gain by feeding the resources, even though the manager is usually permits feeding (around generation 75, the population started going way over the manager’s target).