How to nest setting seeds to compare simulations when testing

[k-doering-NOAA]

We did something like this for SSMSE, but it was unwieldy and there may be a simpler way. We made a list object with seeds for each scenario and iteration. Here's the basic function where we made a list of seeds, just in case it is helpful. : https://github.com/nmfs-fish-tools/SSMSE/blob/302098fd69b9d25988aa5c8460567a467bd3d206/R/utils.R#L671

[Bai-Li-NOAA]

Thanks, Kathryn. The link is very helpful. I added some background information and questions below for further discussion.

Issue Tim mentioned that one thing might be challenging for test comparisons for simulation results is that changes to the order of different calls to simulation within TMB (or R) will change the simulated values and then tests may fail even though it's just because different random numbers are used or the order of the simulation changes through model development.

Potential solutions Tim: This is not a problem with the code. The simulated data from the model and the test just need to be carefully compared or the seed could be checked after each simulation component as a check to make sure that is the issue (and not model coding error).

Also, I think the issue is about the order of sub simulation modules from one iteration run. For milestone 1, we are testing the estimation skill of FIMS by fitting the model to simulated data from the tests, so it may not be an issue for milestone 1. However, it would be good to start thinking about the potential solutions.

Questions @k-doering-NOAA and @nathanvaughan-NOAA, did you nest set seeds based on a global seed value for each iteration? In that way each sub module gets a specific seed value and it always has the same seed value in that particular scenario and iteration. Do you have a link showing how to set seeds for individual sub modules? Thanks!

Task

• [x] Add testing design related suggestions to the FIMS Developer Handbook
• [x] Add software design related suggestions to the software design specification - Milestone 1

[nathanvaughan-NOAA]

Thanks for posting the link for the seed code Kathryn. So what we did for SSMSE was just add a fixed value to the iteration seed that was unique for each random number draw within that iteration ie seed<-iter_seed+fixed_value where the fixed value is something like 1234 or some combo of a fixed value and the year of simulation. The final submodule approach was pretty hacky and could be better particularly if FIMS is well modularized but it got us up and running.

https://github.com/nmfs-fish-tools/SSMSE/blob/302098fd69b9d25988aa5c8460567a467bd3d206/R/runSSMSE.R#:~:text=seed%20%3D%20(iter_seed%5B%5B%22iter%22%5D%5D%5B1%5D%20%2B%20123456)

https://github.com/nmfs-fish-tools/SSMSE/blob/302098fd69b9d25988aa5c8460567a467bd3d206/R/runSSMSE.R#:~:text=seed%20%3D%20(iter_seed%5B%5B%22iter%22%5D%5D%5B1%5D%20%2B%20234567%20%2B%20yr)

[nathanvaughan-NOAA]

An alternative may be to just have a TRUE/FALSE test input to functions that sets the seed to 1 or something for every random draw if it is TRUE which would only be in tests? This would make things easy for testing and have no effect on actual use runs?

[nathanvaughan-NOAA]

One more issue with the set.seed() in R is that it only accepts a 32bit integer input that equates to about 2 billion seeds. That is a lot but could be exceeded in a complex MSE scenario if every (scenario, iteration, submodule) input wants a unique seed. This is only important if the intention is to add MSE capacity into FIMS in the future. Hopefully, R adds 64bit integers to base one day.

[kellijohnson-NOAA]

Do we need to have exactly reproducible results for 500+ iterations or just be able to test things? I vote for the latter because I worry that if we get overly complex about setting seeds and reproducing results that we will lose the ability to accurately capture variance present in the model. Sorry if I am missing the point here.

[nathanvaughan-NOAA]

I think you're on point Kelli and agree. Our concern for SSMSE was being able to replicate and debug failed runs that may have been caused by extreme sample draws, but our approach was probably overkill and has not proven necessary to date. I think having a testing solution is probably sufficient, I just mentioned the seed limitation as a warning in case we went down the complex seed road.

[Bai-Li-NOAA]

Thanks, Nathan and Kelli. Both nested seeds and adding a TRUE/FALSE parameter in each module for setting up testing seeds seem work. I agree that the later one sounds easier to be implemented and maintained, also saves us from the 32bits/62bit integer issue.

[ChristineStawitz-NOAA]

Discussion points from modeling team meeting

• Do we need a separate R OM implemented in R in addition to simulating from the TMB functions? If so is this exported as part of MSE functionality or internal testing use only? For milestone 1 we can use model comparison om code
• There are some issues with using a single iteration of an estimation model as a simulation model because process noise is added randomly without taking into account the correlation structure e.g. temporally or between different processes (@Andrea-Havron-NOAA pers comm with Andre)

[msupernaw]

To do this from C++, we just need to implement our own Rcpp function to call R:

void set_seed(unsigned int seed) { Rcpp::Environment base_env("package:base"); Rcpp::Function set_seed_r = base_env["set.seed"]; set_seed_r(seed); }

On Mon, Jan 24, 2022 at 1:41 PM Christine Stawitz - NOAA < @.***> wrote:

Discussion points from modeling team meeting

• Do we need a separate R OM implemented in R in addition to simulating from the TMB functions? If so is this exported as part of MSE functionality or internal testing use only? For milestone 1 we can use model comparison om code
• There are some issues with using a single iteration of an estimation model as a simulation model because process noise is added randomly without taking into account the correlation structure e.g. temporally or between different processes @.*** https://github.com/Andrea-Havron-NOAA pers comm with Andre)

— Reply to this email directly, view it on GitHub https://github.com/NOAA-FIMS/FIMS-planning/issues/25#issuecomment-1020421639, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABFUSEDSDRMHK6NDXGXZK3TUXWMMVANCNFSM5MLJRYYA . Triage notifications on the go with GitHub Mobile for iOS https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675 or Android https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3Dgithub.

You are receiving this because you are subscribed to this thread.Message ID: @.***>

-- Matthew Supernaw Scientific Software Developer National Oceanic and Atmospheric Administration Office Of Science and Technology NOAA Fisheries | U.S. Department of Commerce Phone 248 - 396 - 7797

[JonBrodziak]

Hello Bai,

I think one can use the "rstream" package to do the random number generation efficiently and with a very long period.

One good source for long period generators is the Rngstreams package of L'Ecuyer et al.: http://www.iro.umontreal.ca/~lecuyer/myftp/streams00

The basic generator for Rngstreams is named "Mrg32k3a" and can be set to be the RNG for R by doing the following:

First check the default generator if curious, should be Mersenne Twister

RNGkind() [1] "Mersenne-Twister" "Inversion" "Rejection"

Load "rstream"

install.packages("rstream") library("rstream")

Now set the RNG in R to be the Mrg32k3a generator, which requires an initial seed of a vector of six integers:

rstream.RNG(new("rstream.mrg32k3a", seed = c(1,2,3,4,5,6)))

Check first value for this seed vector

runif(1) [1] 0.001009498

Check RNG setting in R if curious, should be user specified

RNGkind() [1] "user-supplied" "Inversion" "Rejection"

The period of the Mrg32k3a generator is about 2^(191) or about 10^(57.5) which should be large enough for most of our applications in FIMS.

On Thu, Jan 20, 2022 at 6:51 AM Bai Li - NOAA @.***> wrote:

Thanks, Nathan and Kelli. Both nested seeds and adding a TRUE/FALSE parameter in each module for setting up testing seeds seem work. I agree that the later one sounds easier to be implemented and maintained, also saves us from the 32bits/62bit integer issue.

— Reply to this email directly, view it on GitHub https://github.com/NOAA-FIMS/FIMS-planning/issues/25#issuecomment-1017708345, or unsubscribe https://github.com/notifications/unsubscribe-auth/ACVXIZSG6HX5J7W4N5YBBYTUXA4QNANCNFSM5MLJRYYA . Triage notifications on the go with GitHub Mobile for iOS https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675 or Android https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3Dgithub.

You are receiving this because you are subscribed to this thread.Message ID: @.***>

-- Jon Brodziak, Ph.D. NOAA Inouye Regional Center Pacific Islands Fisheries Science Center 1845 Wasp Boulevard, Building 176, NMFS/PIFSC/FRMD Mail Room 2247 Honolulu, Hawaii 96818 USA Phone: 808-725-5617 Email: @.***

“Wherever my travels may lead, paradise is where I am.” ~ Voltaire

The views expressed in this message are my own and do not necessarily reflect any position of NOAA.

[Bai-Li-NOAA]

Thanks everyone for your suggestions. Below is the summary of the solutions. I will move bullet points 2 and 3 to FIMS developer handbook - chapter 7. I suggest that we include 2.i and 5 in the Appendix of the software design specification doc.

1. For FIMS milestone 1, we will test the estimation skill of FIMS by fitting the model to simulated data from the tests, so no worries there.

2. Once we start developing simulation modules, there are two ways that help compare simulated data from FIMS and a test.

   1. Add a TRUE/FALSE parameter in each FIMS simulation module for setting up testing seed. When testing the module, use the parameter=TRUE to fix the seed number in R and conduct tests.
   2. If approach 2.i does not work, then carefully check simulated data from each component and make sure it is not a model coding error.

3. FIMS will use set.seed() from R to set seed. “rstream” package will be investigated if one of the requirements of FIMS simulation module is to generate multiple streams of random numbers to associate distinct streams of random numbers with different sources of randomness. rstream was specifically designed to address the issue of needing very long streams of pseudo-random numbers for parallel computations. Please see rstream paper and RngStreams for more details.

4. It is not necessary to develop a separate R operating model for internal testing. If FIMS has a modularized simulation engine, developers can use unit tests and functional tests to test FIMS simulation sub modules.

5. The process noise needs to be added to the FIMS simulation model with consideration of correlation structure because process noise is added randomly without taking into account the correlation structure if a single iteration of an estimation model is used as a simulation model.

[JonBrodziak]

On 3), please note that rstream was specifically designed to address the issue of needing very long streams of pseudo-random numbers for parallel computations and is an upgrade over L'Ecuyer's earlier package.

On Wed, Jan 26, 2022 at 6:26 AM Bai Li - NOAA @.***> wrote:

Thanks everyone for your suggestions. Below is the summary of the solutions. I will move bullet points 2 and 3 to FIMS developer handbook - chapter 7. I suggest that we include 2.i and 5 in the Appendix of the software design specification doc.

1.

For FIMS milestone 1, we will test the estimation skill of FIMS by fitting the model to simulated data from the tests, so no worries there. 2.

Once we start developing simulation modules, there are two ways that help compare simulated data from FIMS and a test.

1. Add a TRUE/FALSE parameter in each FIMS simulation module for setting up testing seed. When testing the module, use the parameter=TRUE to fix the seed number in R and conduct tests.
   2. If approach 2.i does not work, then carefully check simulated data from each component and make sure it is not a model coding error.

2. FIMS will use set.seed() from R to set seed. “rstream” package will be investigated if one of the requirements of FIMS simulation module is to generate multiple streams of random numbers to associate distinct streams of random numbers with different sources of randomness.

3. It is not necessary to develop a separate R operating model for internal testing. If FIMS has a modularized simulation engine, developers can use unit tests and functional tests to test FIMS simulation sub modules.

4. The process noise needs to be added to the FIMS simulation model with consideration of correlation structure because process noise is added randomly without taking into account the correlation structure if a single iteration of an estimation model is used as a simulation model.

— Reply to this email directly, view it on GitHub https://github.com/NOAA-FIMS/FIMS-planning/issues/25#issuecomment-1022367150, or unsubscribe https://github.com/notifications/unsubscribe-auth/ACVXIZSVYHWRYH5PAQVRB5LUYAOD3ANCNFSM5MLJRYYA . Triage notifications on the go with GitHub Mobile for iOS https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675 or Android https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3Dgithub.

You are receiving this because you commented.Message ID: @.***>

-- Jon Brodziak, Ph.D. NOAA Inouye Regional Center Pacific Islands Fisheries Science Center 1845 Wasp Boulevard, Building 176, NMFS/PIFSC/FRMD Mail Room 2247 Honolulu, Hawaii 96818 USA Phone: 808-725-5617 Email: @.***

“Wherever my travels may lead, paradise is where I am.” ~ Voltaire

The views expressed in this message are my own and do not necessarily reflect any position of NOAA.

[Bai-Li-NOAA]

Thanks @JonBrodziak. The bullet point 3) has been updated.