DrylandEcology / STEPWAT2

folder
4 stars 5 forks source link

Annual species biomass values intermittently 0, despite forced establishment #18

Closed kpalmqui closed 6 years ago

kpalmqui commented 7 years ago

For some years - bromus tectorum relsize and biomass values are zero, even when probability of establishment is 1.0 (which forces establishment every year). This suggests something is going wrong with either the number of viable seeds in the seedbank or the calculation of num_est (the number of seedlings to add to the plot each year). See add_annuals function in ST_resgroups.c

2hua commented 7 years ago

Redefined the seed bank. When seedbank = 0, s->seedprod[i] = RandUniRange(1, s->max_seed_estab); Else seed production = s->max_seed_estab * lastyear_relsize.

kpalmqui commented 7 years ago

We do not need to redefine the seed bank. X is the number of viable seeds in the seedbank in the current year (see function _get_annual_maxestab). What we need to redefine is how (and how many) seedlings are added to the plot (i.e. what fraction of x). As such, we need to retain x here. In those below commits it is commented out.

These commit needs to be reverted: https://github.com/Burke-Lauenroth-Lab/STEPWAT2/commit/242f3cf389b00d12a00946caf36e96c8ee60f8b0

https://github.com/Burke-Lauenroth-Lab/STEPWAT2/commit/c40cd56f9e0f3f743a3cb17786145a57bb8de6e1

2hua commented 7 years ago

Yes, I ignored the X here.... Seed bank parameter need not be redefined. I will retain this X value and use x * seedling establishment rate. What the most different from the original workflow for annual establishment is the final return value. My opinion is new_size is the final return value rather than sum_size.

2hua commented 6 years ago

All the changes will be done based on this version of annual codes: https://github.com/Burke-Lauenroth-Lab/STEPWAT2/blob/f7206482cacbe24144cb53f2d2a0141b3ed741a6/ST_resgroups.c

  1. add_annual function a.Modify newsize = x lastyear_relsize to newsize = x seedling establishment rate no matter last year resize is zero or not. b.Delete the newsize = RandUniRange(1, x) when last year resize is 0. c.Transfer newsize to an integer. As group_establish function: num_est = _add_annuals(rg, sp, Species[sp]->lastyear_relsize). 2._add_annual_seedprod function (Need to be decided) Replace Globals.currYear == 1 to Lastyear relsize ==0 then there will be enough seed bank in the plot for establishment. Which is not quite sure as if we use Globals.currYear as a decision statement, the seedbank will always be super low as the decay of each annuals, which is not accordance with the real nature about relatively higher amount of seedbank in each plot.

@kpalmqui @dschlaep

kpalmqui commented 6 years ago

1a. agree. Let's also change newsize to num_est for clarity. 1b. agree. 1c. agree.

  1. I think you are proposing to change If(Globals.currYear==1) with If(lastyear_relsize==0)? The question is whether we want to add seeds to the seedbank if there was no standing biomass in the previous year. Conceptually, you would not expect seeds to be added to the seedbank in the current year for a given species if that species was not present in the previous year. So conceptually, I think this code makes sense. In order to accumulate seeds in the seedbank for annuals, we may just need to adjust the max_seed_estab parameter - It is currently set to 20 for cheatgrass. @2hua is this too low of a value based on values you have found in the literature?
2hua commented 6 years ago

In some references the max_seed_estab could reach 1800 per year per square meter in wet year or 200 in dry year. s->max_seed_estab this parameter especially for cheatgrass need to be more precise. Obviously, 20 is too small for cheatgrass. Following around 30 refs, 1800 is the max number I have seen for cheatgrass. @kpalmqui

2hua commented 6 years ago

After these three changes, the s->max_seed_estab parameter of annual which should be a fixed value will become 0 after year 2. @kpalmqui
Commit as: https://github.com/Burke-Lauenroth-Lab/STEPWAT2/blob/annual_testing/ST_resgroups.c

kpalmqui commented 6 years ago

@2hua perhaps you can put together a small table that summarizes the max_seed_estab parameters (calculated for 1m2) you have found from the literature to guide our selection of a max_seed_estab parameter?

I am not sure I am following your second comment - the max_seed_estab is a fixed parameter and read from input, so how can it change after year 2?

2hua commented 6 years ago

Cheatgrass EIND index list was updated on Dropbox Link: https://www.dropbox.com/s/x7utfv1p7k72ynv/Cheatgrass%20EIND%20index%20list.xlsx?dl=0 The largest number I could find is 2200. @kpalmqui

2hua commented 6 years ago

@kpalmqui fc5d0ee2a0e7fa39a3d736be1838a53387c6f35c This commit modified the loop of seed production . The definition of seed production in this species structure is an array. if we use (s->max_seed_estab) * (s->lastyear_relsize) to show seed production in the current year. Should we change the array to a current year seed production value?

2hua commented 6 years ago

1a modification has been pushed for you to review. 4701c29849409ee17a603f024a7a70fce3bf9eb6 @kpalmqui @dschlaep

kpalmqui commented 6 years ago

Just tested and yes now there are some non-zero values (especially when you increase pestab values for annual species). Ready for the second commit.

dschlaep commented 6 years ago

I find a commit with for example a title "Isssue #18 1a modification" and no further text useless. For a few days, we may know by heart what '1a modification' is, but we will not remember this in a few weeks and every new person starting working on this code will have no clue. Digging through the thread of issue #18 to find the explanation of '1a modification' is a wast of time even if successful.

I strongly urge that we all attempt to write better commit (and issue) messages following the guidelines for our repositories and incorporating additional advice:

Why good messages are important: Erlang: Writing good commit messages: "Good commit messages serve at least three important purposes:

How to write good messages: Who-T: On commit messages: "A good commit message should answer three questions about a patch:

dschlaep commented 6 years ago

I updated our guidelines: https://github.com/Burke-Lauenroth-Lab/workflow_guidelines#communication-commit-messages-comments-on-issues-etc

2hua commented 6 years ago

@dschlaep Thank you Daniel for the instructions about how to write commit. Sorry for the unclear issue 18 1 a modification last time. This new commit is based on issue 18 1b. I made a precise commit: "Delete the newsize = RandUniRange(1, available seeds) when last year relative size is 0" to show a more reasonable annual establishment mechanism. 0b4a287cec6c4d952263f51b2f3d9ebe55b2f190 This commit maybe need two of you to test @kpalmqui @dschlaep .

kpalmqui commented 6 years ago

@dschlaep Yes, thank you Daniel - I appreciate you updating our guidelines and I need to do a better of this myself.

@2hua I have tested the latest commit and I am getting results similar to what is documented in our Github issues posted in June and July #26 #18 - annual biomass is occasionally non-zero, but not very frequently and often too high. Non-zero values of course increase as both pestab (probability of establishment) and eind (max number of indivs that can establish in a year) increase. However, a.cool.grass biomass consistently goes to 0 after about 20 years and remains 0 throughout the rest of the simulation (this is despite increasing pestab to 0.5 and eind to 100). This was the problem that was documented in June and July and it continues to be a problem in the annual2 branch despite a few code changes. This could either be caused by:

1) there is no viable seed (you should track the viable seedbank (x) each year to see if this is the problem). 2) how we add seedlings to the plot each year needs to be reconsidered (currently line 260 of ST_resgroups.c): newsize = x * s->seedling_estab_prob; where newsize is the number of seedlings added, x is the size of the viable seedbank, and seedling_estab_prob is the pestab for each species.

Let me know what you find with the viable seedbank. Let's discuss the cause of this and potential solutions on Friday (or earlier via Github). Let's not push any commits to deal with this problem until we establish the cause of the problem and identify the solution.

2hua commented 6 years ago

The viable seed was tracked. Same results were found as you described. As the seed production = max seed establish * last year relative size. I changed the eind value to 2200 in cheatgrass but the viable seed will be gone after 20 years as well. Daniel and I talked yesterday afternoon about this issue a little bit . One of the possible problems may be the random number for annual establishments rate. For cheatgrass, the rate is 1.5%, only 1.5% possibility could call seedbank parameter. Meanwhile 98.5% possibility of last year relative size is 0. There will be no enough seed adding into seed bank. @kpalmqui @dschlaep

kpalmqui commented 6 years ago

@2hua That is exactly why I testing increasing pestab to 0.50 (as described above). Probably best to use the input or code parameter name, (pestab), rather than "annual establishments rates" in this conversation, so we are clear.

In terms of the seedbank - what do you mean by "same results were found as you described". I did not describe tracking of the seedbank in my above comments. This would include not just seed production in the current year, but the total viable seedbank (x).

2hua commented 6 years ago

x (viable seeds) parameter for establishment year was tracked of the current annual2 branch and all the input parameters are following the current GitHub values: Current year 4 Sp name cryp seedbank = 1.25000 Current year 15 Sp name cryp seedbank = 2162.37500 Current year 17 Sp name chen seedbank = 74.00000 Current year 27 Sp name cryp seedbank = 1807.57788 Current year 35 Sp name cryp seedbank = 17980.61719 Current year 45 Sp name brte seedbank = 18460.00000 Current year 50 Sp name cryp seedbank = 17316.57227 Current year 70 Sp name brte seedbank = 0.00000 Current year 86 Sp name chen seedbank = 1480.00000 Current year 104 Sp name cryp seedbank = 39298.87891 Current year 105 Sp name cryp seedbank = 22327.46680 Current year 120 Sp name brte seedbank = 0.00000 Current year 129 Sp name chen seedbank = 709.32617 Current year 132 Sp name cryp seedbank = 8742.25586 Current year 134 Sp name cryp seedbank = 3506.61108 Current year 167 Sp name cryp seedbank = 5756.00000 Current year 169 Sp name brte seedbank = 0.00000 Current year 192 Sp name chen seedbank = 626.52502 Current year 200 Sp name chen seedbank = 311.79761 Current year 220 Sp name chen seedbank = 205.00000 Current year 232 Sp name brte seedbank = 0.00000 Current year 248 Sp name brte seedbank = 0.00000 Current year 289 Sp name chen seedbank = 338.75714 Current year 293 Sp name cryp seedbank = 5151.89990 @kpalmqui From the results, I think the x value is reasonable for the annuals establishment.

kpalmqui commented 6 years ago

This problem is caused by: 1) depletion of the seedbank because lastyear_relsize is often 0 and hence no seeds are added to the seedbank. Lastyear_relsize is often 0 due to the random number draw comparison to pestab, before we add seedlings to the plot. If we fail to add seedlings in multiple,back-to-back years (because the random number draw is not less than pestab), lastyear_relsize will be 0 and no seeds will be added to the seedbank in the next year. For a species like Bromus tectorum, which has viable seed for only 3 years, this issue will result in the viable seedbank going to 0 for the remainder of the simulation.

Solution: Remove this random number draw and instead draw a value from a distribution (range 0 to 1) with mean = pestab set in inputs and variance = x. More details on the solution to be added by @2hua.

2hua commented 6 years ago

Beta distribution could be applied for the distribution (range 0 to 1) with mean = pestab set in inputs. Meanwhile, two variance values (beta and alpha) will be set in species database structure.

dschlaep commented 6 years ago

The beta distribution is normally parametrized by the two shape parameters alpha and beta. They are not "variances" (http://mathworld.wolfram.com/BetaDistribution.html). However, the mean and variance are functions of alpha and beta:

mean(X) = alpha / (alpha + beta) var(X) = (alpha beta) / [(alpha + beta)^2 (alpha + beta + 1)]

As discussed in our meeting yesterday, we want to parametrize the random draw by the mean and variance of the beta distribution, i.e., pestab and a new input representing the variance. alpha and beta need to be calculated from those two inputs accordingly.

2hua commented 6 years ago

General beta distribution random function was created by the last commit. My thought is calling the beta random number from the Global random definition file rather than just defining the beta function in this annual establishment module. At the same time the input alpha and beta should be a function of pestab and a new input variance value. I will rethink about how to express alpha and beta by pestab and the new input var.

2hua commented 6 years ago

mean(X) = alpha / (alpha + beta) var(X) = (alpha beta) / [(alpha + beta)^2 (alpha + beta + 1)] Implement the derivations of alpha and beta by mean and var in the annual establishment. Read the var and pestab parameters from species.in. Now adding codes of the derivation of these binary quadratic equations is the next step. @kpalmqui

2hua commented 6 years ago

The roots of these binary quadratic equations are: alpha = (m^2-m^3-vm) / v beta = (m-m^3-v+vm) / v When we get the m = pestab and v (new parameter read from species.in), the alpha and beta could be calculated and the beta distribution random number draw could be applied for the annual establishment.

2hua commented 6 years ago

Output of biomass of new beta code implement was uploaded on dropbox: https://www.dropbox.com/s/lbfh3cu61ag32v9/var%20index%20for%20annual.xlsx?dl=0 The variance parameter needs to be more accurate. Now I still have no idea about a specific instruction to users about the suggest value for variance. @kpalmqui If you think the logic is right, maybe I could make a new commit based on this implement,

kpalmqui commented 6 years ago

@2hua I am not sure what your question is or what you would like me to do. Are you asking what variances should be specified for each species?

2hua commented 6 years ago

The results from dropbox shows the sensitivity of this variance parameter( e.g.: even between 0.0001 and 0.001, the results were totally different), therefore this parameter needs to be explained to the users about how they could define this parameter. And probably we need to give a default value for each species. When we plan to define this variance value range, I suggest to take Brummer et al 2016 as a reference: using boosted regression tree approach to determine the importance of climate and local biotic and biotic factors in cheatgrass establishment. Now the simplified way is just take the conclusion of this paper as consideration: July precipitation, mean temperature during the driest quarter and native grass canopy cover should be taken as trade-offs for the variance of cheatgrass. BTW, the boosted regression tree methods are here: https://github.com/dmlc/xgboost More quantitive informations will be gathered from Brummer's paper and supplemental materials. @kpalmqui

2hua commented 6 years ago

Three more quantitive numbers should be taken into consideration: 1) 10 mm July PPT 2) 15C temperature of driest annual quarter 3) 25% coverage of native grass The first 2 have positive relation and the 3rd has negative relation. The 4th factor in this paper is another invasive species crested wheat grass, which was not shown in this current model and should be taken into consideration in the further implementation. Now my opinion is using these first 3 factors to model this regression tree of variance parameter.

2hua commented 6 years ago

Forb biomass test data @kpalmqui : forb biomass

2hua commented 6 years ago

forb test2 Pestab=0.95;Eind=10;Var=0.0001;A.warm.MaxBio=4.700;A.cool.forb.MaxBio=4.429

2hua commented 6 years ago

Pestab=0.85;Eind=10;Var=0.0002;A.warm.MaxBio=4.700;A.cool.forb.MaxBio=4.429 forb test3

2hua commented 6 years ago

Pestab=0.35;Eind=10;Var=0.0001;A.warm.MaxBio=4.700;A.cool.forb.MaxBio=4.429 forb test4

2hua commented 6 years ago

Raw data: forb.test.xlsx