cgrandin / csas-latex Goto Github PK

@grinnellm
Can you please update as much as you can for Section Data summaries for major SARs

Simplify where you can (remove statements that involve calculations like % increase in mean weight age)

Flag sentences that you need me to add to.

Eg, for HG:
The survey index in 2017 is approximately xx tonnes, declining from xx tonnes in 2016 (Figure \ref{figSpawnIndex}). The 2017 survey index is the lowest survey index since the 1960s (check).
Mean weight at age for ages 2:6 has been stable (check) and increasing slightly since 2010 (check).
The proportions-at-age were dominated by age 4 fish (Figure \ref{figProportionAge}) and ...something about older age classes (and ref that figure)
There were xx seine-test vessel samples (each of approximately 100 fish) used to calculate mean weight at age and proportions at age data in 2017 (Appendix \ref{chap:input}).

Please make the rows for Tables 8-12:
SB_0
0.30SB_0
SB2017
SB2017/SB0
SB2018

and please add this table for AM1 - either same table or 5 more tables.

I can't generate the document on my mac at home; having the same issues that Jaclyn and I were having on her mac. Not sure why this is yet: a mac issue, or something else?

Couple edits:

-please change percentiles to 5th, 50th, and 95th. (it was a past herring subcommittee decision to report 90th credible intervals, not 95th)
-I compared the spawning biomass values to the 2016 assessment and they are very similar, with the exception of the last couple years (I will investigate this).
-Please remove the 2018 projection year from the table (the final row)
-Table header: please change "relative spawning biomass" to "Depletion (sbt/sbo)"
-The depletion values (percentiles)... these are puzzling. Please check the calculation- should be sbt/ sbo (not bo)
-the MPD and median sbt values match :)
-the MPD and median depletion values do not match...

Here's what I have so far, plotting MPDs for AM2 only, spawning biomass vs recruitment. Still needs the lines.

can we add a time series plot of sbo values?
I don't think this is currently an 'output' saved from the mcmc.

@grinnellm Can you please change \Sexpr{fish.name} to capitalize Pacific Herring?
This as per Lisa/ CSAP instructions. Thx.

Need to get the citations from the Bridging Analysis into bibtex. These will be used in the rest of the document too.

Can we create a table of maturity vectors used for the Maturity sensitivities?
Added Maturity-vectors.csv file to the data folder.

Working on this now. It'll be a table in the input data appendix. The only thing I don't like about this table is that the data comes from the biosamples database, not the stock assessment input or output files. Because some folks don't have access to the biosamples database I'll make some *.csv's with the number of biosamples and put them in the models/ directory. Then the script ResearchDocumentV2.R will be able to grab the info from those csv's, which should work for everyone. Not ideal, but it'll only have to be updated once, before we do next year's assessment.

The document can now bring in two appendices: app/A_biology.rnw and app/B_equations.rnw (these were causing errors). I think the errors were because these appendices were trying to load files for tables and equations, but those files aren't being generated yet (hopefully the required ones will be generated when we have the output for multiple stocks).
In A_biology, I commented out the R code to make 3 plots (fig:vonb, fig:lw, and fig:matage) -- Jaclyn should know which, if any, of these are required.
In B_equations, I commented out the code to make 4 tables (tab:variables, tab:steady-state-model, and tab:catch-age-model).

column 2 header: P(SB_2018 < LRP of 0.30SB0)
(the LRP is new for herring so wanting to be clear that 0.30SB0 is the LRP)

column 3 header: Med(SB_2018/ 0.30SB0)
(pls check the equation matches. This is intended to show the magnitude of difference between SB_2018 and 0.30SB. Ie, if the median ratio is 1.0 then SB_2018 ~= 0.30SB0, whereas if the median ratio is 2.5, then SB_2018 is 2.5x greater than 0.30SB0 (and so on).

column 4: P(SB_2018 < 10,700 t)

column 5: Med(SB_2018/ 10,000 t)

cols 6:9 - great as is.

for AM1: same edit to header needed for col 3

From the final Arrowtooth Res Doc, we would like to include:
Tables 1 and 2- NO b/c these are covered in Input Data Appendix (by MG)
Table 3- NO (survey inputs covered in Input Data Appendix (by MG)
Table 4- YES pls - Table label should state "initial values" (from par file)
Table 5, 6- Yes pls
Table 7- YES pls (single sex model, not female SSB), and ADD Depletion (as per Table 6 of 2014 Res Doc (./github/csas-latex/doc/CSAS 2014 Working paper.docx)
Table 8- YES pls- this needs to be age-2 recruits (not age-1)
Table 9- NO for now. I'm not sure how to add F's across gear types.
Table 10- YES but for herring U_t is calculated a bit differently- ask JC
Table 11- NO- b/c adding depletion to Table 7
Table 12, 13- leave out for now
Table 14- Decision tables- JC will add a separate "issue"

Figure 1- NO (MG already added BC stock map)
Fig 2- NO (already added as 5-panel fig)
Fig 3, 4- NO (we've included age-prop plots - by MG)
Fig 5- NO
Fig 6- NO for now
Fig 7, 8, 9, 10, 11- YES - for Fig 8 - make the prior a darker black lline (lwd=2) - it's not showing up when printed
Fig 12- NO (b/c appears in storyboard)
Fig 13-16- NO
Fig 17- YES
Fig 18- NO (b/c only including MSY ref pts in sensitivity section)- as a table
Fig 19, 20- NO
Fig 21- NO (b/c in storyboard)
Fig 22-42- No for now. Haven't decided yet on figs for sensitivity cases

I am unable to compile this morning and I am receiving the following error:
! LaTeX Error: File `pdfpages.sty' not found.
==> Fatal error occurred, no output PDF file produced!

Was the \usepackage{pdfpages} just added to the res-doc.rnw or have I lost the package somewhere?

I have added an appendix that brings in "BridgingAnalysis.pdf" (now should be in "csas-latex/doc"). The page numbers get added automatically, as are headers. We need to remove page numbers and headers from "BridgingAnalysis.pdf" -- you'l see when you compile it.

Chris: a couple of edits please to Table 5 (Posterior 5th percentile, median, 95th percentile) and MPD estimates of key parameters...): I think it makes sense for the list of estimated parameter values in Table 5 to be the same as the parameters listed in Table 3 (Ro, Steepness, M, Rbar, Rinit_bar, rho, kappa, sigma, tau, q1, q2)- also matches Fig17 (Prior probability distributions and posterior histograms) Thanks, ~J

*note we may need to calculate sigma and tau. Let me know...

a plot to show distribution of sbt in 2017 (histogram), with 0.30SBo (LRP) so that we can "see" probability of being below LRP in 2017

as discussed: age proportion bubble plots (model estimates) for each gear type (x3) x 5 stocks (AM2 only)

Table 8 SB_2018 value for HG AM2 appears to be different than the value plotted in Fig 7 (d).

My hunch is that the value plotted in Fig 7 (d) is what is habitually used to reflect "projected spawning biomass assuming no fishing, SB_2018". This is calculated from the TAC=0 rows in iscammcmc_proj_Gear1.csv, and uses "years" specified in the pfc file to calculate mean weight-age and average M.

My hunch is Table 8 SB_2018 values are median values of sbt1_2018 from iscam_sbt_mcmc.csv, which is calculated from the BH S/R assuming long-term average weight-age and M.

If my "hunches" are correct, can the SB_2018 values in Table 8 be changed to the values from iscammcmc_proj_Gear1.csv.

@cgrandin I used something like this to calculate the median values and 90%iles:

p_sbt <- cbind(iscammcmc_proj_Gear1$TAC,iscammcmc_proj_Gear1$B2018)
p_sbt_median <- tapply(p_sbt[,2],p_sbt[,1],quantile,prob=0.5,na.rm=TRUE)
...
@grinnellm probably has a fancier R routine

@grinnellm
Can you please fix the Intro to make the citation for the PA zone figure work? (search CITE PAzone figure). Thx

sensitivity case #3: changing the sd for prior on q

It seems the old (i.e., previous version with bo, not sbo) "*mcmc.csv" files are present in the folder, in the /mcmc subfolders. This is getting in the way of my script -- it thinks that these are the correct files, and doesn't move the new files to their correct place. That may not make sense. Can I delete all the /mcmc subdirectories from 2017mcmc? And the .tex and .png files too? That'll make it smaller too.

Pls use _SB__0
In the bridging doc and my writing so far, we've used _SB_0 to denote unfished spawning biomass (italicize SB and not 0). Andy will cringe but in this case I think consistency with the code makes it easier to not mistake this for B_0.

I get an error when I try to compile the PDF (RStudio):
Quitting from lines NA-488 (./maindoc/3_modelling.rnw)
I get the new q-priors runs, and ran delete.rdata.files() and build(), but no luck.

Obs_proc_error_sensitivities.pdf

Do we need to show all stocks?
There are 2 table (AM1/AM2) for each stock but they all show the same thing that regardless of initial rho and kappa values the estimated values for all leading parameters are the same if estimation is on for both rho and kappa.

We will need a table with the Fmsy values estimated for the maturity schedule base case and the sensitivity case.

Jaclyn, this is what we wanted to present for the results here correct?

In order to quickly get the source("ResearchDocumentV2.R") to work, I hacked the function GetPars(). I just did this for loading the recruitment output because the number of columns output in iscam_rt_mcmc.csv is 2 less than the length of yrRange. I haven't had time to figure it out and hoped Matt could help.

The change I made breaks Matt's age-2 recruitment plot. I also noticed that the natural mortality plot is wrong, probably due to something wrong with the M output file iscam_m_mcmc.csv. We can work together to fix these issues.. This is on commit bd32677

This is the change in GetPars():
raw <- fread( input=file.path(SAR, model, "mcmc", fn) ) %>% as_tibble( ) %>% select( 1:length(ifelse(varName=="Recruitment", yrRange - 2, yrRange)) )

can you make the prop-age barplot into a two panel (stacked) figure that includes bubbleplot for gear2

Chris- my take from our meeting last Friday is that you have the structure already in place from Arrowtooth to put multiple sensitivity cases on one plot. (pls correct me if I'm wrong)

sensitivity case #1: compare estimated time varying natural mortality (base) and constant estimated M

For each stock we want to plot:

(1) model fit to survey data (MPD estimates): 4 lines (constantM, timevaryingM, AM1, AM2)
(2) SB ~ t: 4 lines (constantM, timevaryingM, AM1, AM2)
(3) recruitment deviations ~ t: 4 lines (constantM, timevaryingM, AM1, AM2)
(4) **q ~ t (**MPD estimates): 4 lines (constantM, timevaryingM, AM1, AM2)
(5) nage2 ~ t: 4 lines (constantM, timevaryingM, AM1, AM2)
(6) mortality ~ t: 4 lines (constantM, timevaryingM, AM1, AM2)

And needing to see the MCMC diagnostics, so please add plots of:
priors vs. posteriors
mcmc trace plots
pairs plots

We won't need to include them all in the Res Doc but we can't look at them without first including them...

Sarah (thanks!!!) has rerun the mcevals for all mcmc. The new files are in Herring/2017mcmc/models directory -- sbo is now in the iscam_mcmc.csv file.

The build is failing because a variable called namesNBio in C_input.rnw hasn't been put into custom-knitr-variables.r.

@cgrandin - there are still two Tables one for each of AM1 and AM2 that use labels that look like this:

xlabel = "tab:wcvi-spawn-biomass"

Is it safe to assume the citations that are still hard-coded into the introduction need to be worked into all.bib? I've started on this already..

Chris,
When you get a chance could you please change the Fishing Mortality plots to the estimated Natural Mortality plots for the Natural Mortality sensitivities?
And can you confirm what the CIs for the inputs are based on in the index fits plots?

Thank you!

I can't remember what we decided for this plot of SSB/(SSB+C) by year:

• SSB is the MPD, or
• SSB is the median, and 90% CI,
• other?

Hey Matt, there are no mcmc outputs in the AM1Mat and AM2Mat directories so it stops on each of them. It would be good if it just continued

Table 3. Estimated and fixed parameters and prior...

row 6: change to "Variance ratio, rho ()"
row 7: change to: "Inverse total variance, kappa ()"
row 8: remove row "Survey age 50% selectivity" - we don't use this
row 9: change to: "Fishery age at 50% logistic selectivity ()
row 10: remove row "Survey SD of logistic selectivity" - we don't use this

q-priors - remove from this table and include separate (new) table to show q priors for AM1 and AM2 (from ctl files)

note : With the q-prior removed, this table is the same for all 5 stocks and AM1 and AM2. Can we consolidate to 1-table and reflect this in the header?

Chris, FYI:
Matt and I are going to meet tomorrow AM to review model diagnostics and determine if there's anything else we need. Off the top of my head, I'm recalling:

-age proportion bubble plots for each gear type (x3)
-a plot to show distribution of sbt in 2017, with 0.30SBo (LRP) so that we can "see" probability of being below LRP in 2017

Tweak Rob's R code to generate production plots using iSCAM V2 output.

Chris, my inability to compile the pdf has become inefficient. Do you have time Wed AM to try to troubleshoot this with Matt and I?

Arrg - Sorry! The TVM and constantM legend items are all opposite. I'm fixing it but wanted to get this up in case you are working on it. I hope it doesn't mean rewriting whole sections

@cgrandin
I found two options for rmdir
rm -R knitr-cache
rm -iR knitr-cache

The latter gives a prompt to ask 'are you sure'. Is the object to delete everything in the knitr-cache, including the folder? Including __packages? If yes, then seems the first option is best. Just wanting to check given it's final.

I had forgotten to add the code to bring in the iscam_m_mcmc.csv file and calculate the quantiles which are necessary for the q sensitivity 4-panel plots, so unfortunately this has to be done.

I thought you had planned to run them but I don't see them in the Herring repository.

I've just made a commit (8de595c) which adds the directory names for the natural mortality sensitivities. You need to rename and copy the directories over into the models directory and name them with the following format:

HG-natural-mortality

Then you can source("ResearchDocumentV2.r") again and that will make the mcmc directories for you.

I'm working on getting the code to load these in but it's taking quite some time as there are now 20 more models to load.

Please add:
Appendix D "Assessment of Minor Stocks"
Appendix E "Bridging Analysis" - this will be a PDF that we need to merge/ fix format for later... Chris suggests to make this the final Appendix to make page numbering easier.

Question: Is the order of the tables determined by the order in which they are cited in the text, or set in latex?

Table 1, 2- can you confirm these are 90%iles
Table 2- can you please modify into two tables:

1. col headers: Stock | Model | SB0 | 0.30*SB0 | Lower | Median | Upper
2. col headers: Stock | Model | SB2018 | Proportion age 3 in 2018 | Proportion age 4+ in 2018
   *cols 4 and 5 to include Lower | Median | Upper (90th percentiles)

Please add these tables for AM1 (and keep the AM2 tables too)