nmfs-ost / ss3-source-code

The source code for Stock Synthesis (SS3).
https://nmfs-ost.github.io/ss3-website/
Creative Commons Zero v1.0 Universal
37 stars 16 forks source link

[Bug]: Report 19 mature_bio and mature_num == 0 #382

Closed kellijohnson-NOAA closed 1 year ago

kellijohnson-NOAA commented 1 year ago

Describe the bug

Looking back since the start of 2022, two columns (mature_bio and mature_num) in report 19 have been populated with 0s or incorrect values, where before that they were equal, which is also incorrect. Here is a snippet of Report.sso for hake using 3fa641e

SPAWN_RECRUIT report:19  Function: 3  RecDev_method: 2   sum_recdev: 0.821395
...
Yr SpawnBio exp_recr with_regime bias_adjusted pred_recr dev biasadjuster era mature_bio mature_num raw_dev
S/Rcurve 2.04785e+06 2.17373e+06
Virg 2.04785e+06  2.17373e+06 2.17373e+06 2.17373e+06 2.17373e+06 - 0 Virg 4.09571e+06 4.09571e+06 0.0 
Init 2.04785e+06  2.17373e+06 2.17373e+06 2.17373e+06 2.17373e+06 - 0 Init 4.09571e+06 4.09571e+06 0
1946 2.04785e+06 2.17373e+06 2.17373e+06 815827 594181 -0.317019 1 Init_age 4.09571e+06 4.09571e+06 -0.317019
1947 2.04785e+06 2.17373e+06 2.17373e+06 815827 747241 -0.0878141 1 Init_age 4.09571e+06 4.09571e+06 -0.0878141

With 1486a3d both mature_bio and mature_num are 0.

SPAWN_RECRUIT report:19  Function: 3  RecDev_method: 2   sum_recdev: 0.821395
...
Yr SpawnBio exp_recr with_regime bias_adjusted pred_recr dev biasadjuster era mature_bio mature_num raw_dev
S/Rcurve 2.04785e+06 2.17373e+06
Virg 2.04785e+06  2.17373e+06 2.17373e+06 2.17373e+06 2.17373e+06 - 0 Virg 0 0 0.0 
Init 2.04785e+06  2.17373e+06 2.17373e+06 2.17373e+06 2.17373e+06 - 0 Init 0 0 0
1946 2.04785e+06 2.17373e+06 2.17373e+06 815827 594181 -0.317019 1 Init_age 0 0 -0.317019
1947 2.04785e+06 2.17373e+06 2.17373e+06 815827 747241 -0.0878141 1 Init_age 0 0 -0.0878141

Finally, here is what it looks like with e229a3f, which is v 3.30.20, where mature_bio is zero and mature_num is populated with the same value for every year.

SPAWN_RECRUIT report:19  Function: 3  RecDev_method: 2   sum_recdev: 0.821403
Yr SpawnBio exp_recr with_regime bias_adjusted pred_recr dev biasadjuster era mature_bio mature_num raw_dev
S/Rcurve 2.04785e+06 2.17373e+06
Virg 2.04785e+06  2.17373e+06 2.17373e+06 2.17373e+06 2.17373e+06 - 0 Virg 0 8.08753e+06 0.0 
Init 2.04785e+06  2.17373e+06 2.17373e+06 2.17373e+06 2.17373e+06 - 0 Init 0 8.08753e+06 0
1946 2.04785e+06 2.17373e+06 2.17373e+06 815826 594181 -0.317018 1 Init_age 0 8.08753e+06 -0.317018
1947 2.04785e+06 2.17373e+06 2.17373e+06 815826 747241 -0.0878141 1 Init_age 0 8.08753e+06 -0.0878141

To Reproduce

The input files for Pacific Hake (files used in 2021 assessment) are located here 2022.01.10_base_v2.zip

Expected behavior

I expected mature_num to not be equal to mature_bio.

Screenshots

No response

Which OS are you seeing the problem on?

Windows

Which version of SS3 are you seeing the problem on?

3.30.19 and 3.30.20 (didn't check previous versions)

Additional Context

This is likely tied to the saving of the report file via the report flag not being turned on. Potentially tied to issues #248 and #70.

Rick-Methot-NOAA commented 1 year ago

Kelli, This change occurred because the concept of mature biomass and mature numbers do not explicitly occur in a model that uses empirical wt-at-age. The input of empirical fecundity-at-age encompasses maturity, but does not explicitly address maturity separately from fecundity. I recommend that we update this report to display NA when the model is using empirical wt-at-age. Previously, SS3 was using the parameter-based maturity info to produce these quantities and there were no internal checks to assure that those maturity vectors were consistent with the empirical fecundity-at-age input in the wtatage.ss input file. I'm open to alternative suggestions. Rick

kellijohnson-NOAA commented 1 year ago

Thanks @Rick-Methot-NOAA for the clarification. I think I will default to @iantaylor-NOAA for whatever values works best with r4ss and we can always include a flag that leaves a comment in Report.sso explaining it or put something in the manual.

iantaylor-NOAA commented 1 year ago

This is related to issue #348, which was focused on the BIOLOGY table without thinking about the impact on mature_bio in the SPAWN_RECRUIT table. I note that for purposes of the BIOLOGY, @Rick-Methot-NOAA said that it wasn't possible to report NA values because they were numeric, so 0.5 was going to be used as the placeholder instead: https://github.com/nmfs-stock-synthesis/stock-synthesis/issues/348#issuecomment-1182614539. The BIOLOGY issue applies for both empirical weight at age models at those that use an age-based maturity function. I'm assuming that the mature_bio column in SPAWN_RECRUIT is accurate for age-based maturity and is only an issue for EWAA models. I think either NA or 0 would be fine as a placeholder in SPAWN_RECRUIT and would presumably just require some additional logic to replace the numeric value when WTage_rd == 1.

There is lots of output in Report.sso that doesn't apply to EWAA models, such as all the growth parameters, but it seems most efficient to just focus on areas like this one where the number is just not right and yet not obviously a placeholder.

Rick-Methot-NOAA commented 1 year ago

Let's use this current issue to fix up tables that are misleading with EWAA. Maybe add another issue to use the EWAA flag to turn off more reports now that we have specific table on/off capability

iantaylor-NOAA commented 1 year ago

New issue posted as #383 (title was initially incomplete but later edited to read "turn off some report sections when using empirical weight-at-age").

EJDick-NOAA commented 1 year ago

For non-EWAA models, is it possible to export a complete time series of mature female biomass in mt (i.e. start year to end year)? I'm looking at my most recent assessment (black rockfish in CA), and found mature_bio in the SPAWN_RECRUIT section, but it was only populated for the "Main" era. Happy to post elsewhere if this is not the best place.

Rick-Methot-NOAA commented 1 year ago

You came to the right place. We can add it to the time_series table. You or we can create an issue to that effect; there's a template for that. I checked the code and the quantity exists as: SSB_B_yr(y) += fracfemale_mult make_mature_bio(GP4(g)) natage(t, p, g); SSB_N_yr(y) += fracfemale_mult make_mature_numbers(GP4(g)) natage(t, p, g); So easy to report them anywhere. I'll defer to the r4ss gurus to pick column order in the TIME_SERIES table.

iantaylor-NOAA commented 1 year ago

@EJDick-NOAA, it looks like it's not a question of early vs late, but estimated vs fixed recruitments that lead to reporting or not of mature_bio and mature_num in the SPAWN_RECRUIT table. Either way, you don't have the values in 3.30.21.

If you could use that info prior to the next version of SS3, you could calculate them externally by taking the product of maturity at length and numbers at length by year. The R code is one way to do these calculations (added as columns to numbers-at-length tables for convenience, but the timeseries table is a better place for the SS3 solution). In the plot below the red points are those that you have from SPAWN_RECRUIT (which are only beginning of year, whereas the numbers at length table has mid-year values as well).

# read model
b1 <- r4ss::SS_output("C:/SS/black/black_2023/CA/California black rockfish pre-STAR base
     models/Northern California Black Rockfish/")

# add new columns to numbers-at-length table
b1$natlen$mature_bio <- NA
b1$natlen$mature_num <- NA

# fill in values by row (could be fancier via matrix multiplication)
for (irow in which(b1$natlen$Sex == 1)) {
  b1$natlen$mature_bio[irow] <-
    sum(b1$natlen[irow, names(b1$natlen) %in% b1$lbinspop] *
      b1$biology$Mat * b1$biology$Wt_F)
  b1$natlen$mature_num[irow] <-
    sum(b1$natlen[irow, names(b1$natlen) %in% b1$lbinspop] *
      b1$biology$Mat)
}

# confirm the match with values in SPAWN_RECRUIT
par(mfrow = c(1,2))
plot(b1$natlen$Time, b1$natlen$mature_bio)
points(b1$recruit$Yr, b1$recruit$mature_bio, col = 2, lwd = 3)

plot(b1$natlen$Time, b1$natlen$mature_num)
points(b1$recruit$Yr, b1$recruit$mature_num, col = 2, lwd = 3)

image

EJDick-NOAA commented 1 year ago

Thanks, Rick and Ian! If it's going to be added to the time series table, then it's probably worth including both mature_num and mature_bio. I will try to create an issue using the template.