[ericward-noaa]

This PR is related to this FIMS discussion, NOAA-FIMS/FIMS#1873, and this r4ss issue r4ss/r4ss#1104 (comment), the need to generate composition plots.

Copilot summary:
This pull request adds a comprehensive set of tools for visualizing and diagnosing composition data (age/length) in fisheries stock assessment, including new plotting functions, example data, and extensive documentation. The main changes are the introduction of new plotting functions for composition data and their documentation, the addition of a simulated example dataset, and updates to the package exports and documentation to support these features.

New composition plotting and diagnostics functionality:

• Added new plotting functions for composition data: plot_comp, plot_comp_bubble, plot_comp_mean, plot_comp_residuals, and plot_comp_sex_ratio, along with a data preparation function prepare_comp_data. These functions provide a variety of diagnostic plots for observed and expected age/length compositions, including bubble plots, mean trends, Pearson residuals, and sex ratio plots. [1] [2]

• Added man page documentation for all new functions, describing their usage, arguments, outputs, and examples: plot_comp.Rd, plot_comp_bubble.Rd, plot_comp_mean.Rd, plot_comp_residuals.Rd, plot_comp_sex_ratio.Rd, and prepare_comp_data.Rd. [1] [2] [3] [4] [5] [6]

Example data and reproducibility:

• Added a new simulated dataset example_comp_data (with generation code and documentation) to demonstrate the composition plotting functions. This data includes sex-specific age compositions for a fishery and a survey, with known patterns of misfit for testing diagnostics. [1] [2] [3]

Documentation and navigation updates:

• Updated the package documentation site navigation to include a new article on composition data diagnostics, making it easier for users to find and learn about these new features.

These changes significantly expand the package's capabilities for composition data analysis and provide users with tools and guidance for diagnosing model fit to age/length compositions.

[Schiano-NOAA]

@ericward-noaa is this ready for review?

[ericward-noaa]

Yes @Schiano-NOAA I think so

[Schiano-NOAA]

We can set up a space in the csv that holds captions and alternative text for figures that can be called and made into an RDA which is a cundle of the plot, caption, and alternative text that can be incorporated into a report or other thing.

[Schiano-NOAA]

This works great for the skate example data and provides a good foundation for the comp plots. I have a series of comments and questions and comments that need to be answered. If you don't have the bandwidth to make this work with the standard output, then I can take over from here. I think some of the issue is due to a divergence of the FIMS output and standard output which is not our intention and will be resolved by the end of the year.

I also am wondering how standard some of these approaches are. I apologize for my ignorance on this topic and would be interested to hear if this is done across the board or just the decision made for comp plots for FIMS or specific regions.

There are also some other minor things like function naming and locations that can easily be adjusted.

I thank you so much for this important contribution to {stockplotr} and cannot stress enough the help this brings to our team 😄

[ericward-noaa]

Thanks for the review @Schiano-NOAA -- I think most of these are pretty minor, and I"ll start to push changes to address them. I'll tag you when ready for re-review

[ericward-noaa]

@Schiano-NOAA I think I've dealt with most of the comments above, but there's 2 lingering ones that I think you can help with:

1. You say "We can set up a space in the csv that holds captions and alternative text for figures that can be called and made into an RDA which is a cundle of the plot, caption, and alternative text that can be incorporated into a report or other thing."

2. Related to the other PR request around data, I'm not wedded to the Big Skate data but used it because of compatibility with r4ss. It seems like if the example_data contains observed and expected counts by bin, we wouldn't need the Big Skate data any more.

Given the overlap with r4ss, it'd be great to have @iantaylor-NOAA or @Rick-Methot-NOAA take a look and see how these can be improved

bigskate_composition_plots.pdf

[iantaylor-NOAA]

@ericward-noaa, these are excellent! It's amazing to see how quickly you were able to replicate years worth of work. We live in amazing times.

Overall, I think the plots are great and are more than adequate to replace the r4ss diagnostics.

Here are a few additional changes that could be considered, none of which are vital and all could be added down the road. I'm happy to help with any of these things if folks think it's a good idea to add them.

1. Numeric values in the corner of the bar plots: You've replicated the r4ss behavior of printing the adjusted input sample size and the effective sample size associated with the McAllister-Ianelli weighting. We've moved away from using the McAllister-Ianelli weighting thanks to a Punt paper showing other methods perform better so I'm not sure there's value in printing that number. For FIMS, there's just the single input sample size. When Dirichlet-multinomial likelihoods are added to FIMS, then we could print the input sample size and the adjusted sample size (scaled by the D-M parameter). In theory the numbers could be model dependent--maybe your code already does this. In that spirit, I'm not sure the Adjusted input sample size vs Effective sample size plot is providing value for anyone any longer.
2. Potential improvements to Andre plots: There's a list of potential changes in a slide from an old talk with comments below here: refinements to CAAL diagnostic plots r4ss/r4ss#613.
3. Colors for cohorts in age comps: the hake assessment colors cohorts in age comps using code in https://github.com/pacific-hake/hake-assessment/blob/main/R/plot-age-comp-fit.R (example below pasted from https://www.pcouncil.org/documents/status-of-the-pacific-hake-whiting-stock-in-us-and-canadian-waters-in-2025/). I think @jimianelli-NOAA was the first person I saw using this approach. Note that In the case of hake it's a single-sex model, so there's no need to distinguish females from males and also the fits are based on Bayesian posterior samples so uncertainty is available.

[iantaylor-NOAA]

Here's one more comment that I forgot to include.
4. Sorting panels by column instead of row: This was actually the original reason we didn't use the lattice package for r4ss back in 2008 because sorting by column wasn't available. I think the benefits of the column sorting are clear in the hake example I pasted above, where the march of the cohorts are more easily visible when sorted by column due to the bars also being vertical (e.g. pink 1980 cohort, purple 1984 and 1999 cohorts).

[ericward-noaa]

Great comments, thanks for the fast thoughts @iantaylor-NOAA ! For the plots breaking out sex, r4ss sticks with red-blue -- should we stick to the same convention here?

[iantaylor-NOAA]

For the plots breaking out sex, r4ss sticks with red-blue -- should we stick to the same convention here?

I just found this post: https://www.datawrapper.de/blog/gendercolor which makes me think that we can do better than red and blue. I don't have strong opinions on what alternative to use, however. I personally like the Telegraph 2018 purple and green, but this seems like something that the stockplotr team can decide however they see fit. Thanks to the flexibility of ggplot, I'm assuming that users could theoretically change the colors however they wish (although hopefully a good default will be good enough for everyone).

[Rick-Methot-NOAA]

This all looks great to me. Good work.
I'll reiterate suggestion for an interactive viewer which allows users to select and re-arrange plot elements to create a mind's eyeful collection that facilitates finding patterns in residuals that warrant additional model flexibility to resolve.

[Schiano-NOAA]

@ericward-noaa @iantaylor-NOAA @Rick-Methot-NOAA thanks all for your input and contributions here! I will make sure the standard output has everything needed and if not fix it, but keep the skate example in for now.

In terms of colors, we will definitely need to set a default, but I'd be interested to hear what @sbreitbart-NOAA thinks because she has looked into colors for accessibility which could be our ultimate decision maker if we are between the two options presented above.

For interactive plots, we are looking into it and have a plan to have plots be interactive in the future, particularly for integration into a dashboard where all the plots in {stockplotr} can be viewed nicely 😄

I will be re-reviewing this PR and making edits for the caption and alternative text along with the standard output later today!

[ericward-noaa]

I like @Rick-Methot-NOAA 's suggestion of the visualization tool.

One model that I think works well for this is ShinyStan https://mc-stan.org/shinystan/, which is deployed as a standalone R package for visualizing output from Bayesian models fit with Stan. It's a great tool to look through a bunch of plots quickly and help diagnose what part of a model might be misbehaving. Perhaps we could make a shinyFIMS or broader with stockplotr

[sbreitbart-NOAA]

@ericward-noaa - thank you for taking the time and effort to make such a substantial contribution to stockplotr. At a quick glance, these plots look fantastic.

Re: the color discussion: I suggest we choose two colors from the viridis palette for accessibility. I'm open to different colors as long as they're high-contrast and from a trusted source like viridis. I think this high-contrast green/purple pairing would work well:

scales::show_col(viridis::viridis(2, option = "viridis", end = 0.75))

[ericward-noaa]

@Schiano-NOAA I have brought in most of the comments above, including the default 2-colors proposed by @sbreitbart-NOAA and these ideas from @iantaylor-NOAA

1. Removed numeric values in the corner of the bar plots
2. Allowed CAAL plots to be flipped
3. Allow panels to be sorted by column or row

I didn't touch the caption thing, and left the last suggestion of Ian's (colors for cohorts in age comps) to think about

[Schiano-NOAA]

@ericward-noaa This all sounds good! I am going to do a final check through later today. Additionally, I can potentially draft up a general caption and alternative text for this that I will ask you to review so I can incorporate it into this PR. Then I think it will be ready to merge.

[Schiano-NOAA]

I am approving this for merge so it is available for FIMS to use, but I will be opening new issues in order to add in RDAs for all of the necessary function along with compatibility for all other models.

@ericward-noaa I will request review/drafting captions and alternative text for these plots. Before I do so, do all of the plots made in this PR go into an assessment report and if not, which ones do?

[iantaylor-NOAA]

do all of the plots made in this PR go into an assessment report and if not, which ones do

100% of these plots are used by assessment authors at NWFSC for model diagnostics and model selection which is an essential part of the assessment process. In recent years I would guess that maybe 25-50% of these end up in the final report depending on the authors preferences, with the rest included in an "electronic appendix". The electronic appendix is at minimum a zip file containing all the plots shared with reviewers and once finalized, posted to the council's website along with the model files and assessment report. Some authors also create a github.io page to host the diagnostic figures during the review process since it's easier to access than a zip. As we move toward more standardized reports via {asar}, I think assessment reports could include an even smaller fraction of these plots, but that doesn't remove the need for the diagnostics to be available.

[Schiano-NOAA]

@iantaylor-NOAA thanks Ian! Could you id the ones that do get into the repot? I asked the initial question to determine if pre-written captions and alternative text should be written for some or all of the plots. It would suck to write these and they never get used in the reports.

[ericward-noaa]

@Schiano-NOAA @iantaylor-NOAA I'm happy to help with captions etc. I think a lot of the text can be recycled across plots. There's a few of them that will need to be dynamic, where the axes have the option of being flipped

[Schiano-NOAA]

@ericward-noaa That would be great. You can either start drafting them in the inst/resources/captions_alt_text_template.R or just send them to me in a document. We have a way to make these dynamic, but I can work on that part if you just want to put placeholders where values should go in the captions. If you also want to draft the alternative text, here is a guide on how to write good alternative text for figures.

If you do put the draft captions in the file directly, please put them in a new PR. I am figuring out how to fix these merge conflicts and planning to merge this PR in today.

[iantaylor-NOAA]

@Schiano-NOAA that makes sense.
Since the U.S. West Coast Big Skate assessment was helpful as an example, I'm just going to tabulate what's in the actual Big Skate report from 2019: https://www.pcouncil.org/documents/status-of-big-skate-beringraja-binoculata-off-the-u-s-pacifc-coast-in-2019-october-2019/ vs. the example stockplotr figs shared in an earlier comment in this thread (direct link to PDF).

I should add two things:

1. The table below is just one example. The particular choice of which plots to include in an assessment report have depended on author preferences as well as other issues like the number of fleets (no need for aggregating across fleets if you only have that type of data from one fleet). Also, some of the plots added by @ericward-noaa are new things that we haven't had access to before (e.g. 11_residual_pattern_by_bin), and authors may want to include these if they are available.
2. any dashboard or tool to view all of these plots will also benefit from captions and presumably alternative text if we want to be able to share the dashboard publicly.

stockplotr example	assessment report	comment
NA	Figure 21: Length comp data, aggregated across time by fleet.	stockplotr fig 07 includes model fit (I'm not sure if it's necessary to have separate data-only plot, but there are fleet-specific data-only plots, just not an aggregated version)
02_bubble_obs	Figure 22: Length comp data for all years and fleets. Bubble size indicates the observed proportions, with females in red and males in blue.
20_caal_bubble_obs_WCGBTS	Figure 23 [continued in Fig 24]: Conditional age-at-length data from the current fishery.	stockplotr fig only has 3 out of 15 years
07_distribution_aggregated	Figure 36: Fits to length composition data, aggregated across time by fleet.
08_residuals	Figure 37: Pearson residuals for length composition data for all years and fleets, with females in red and males in blue. Closed bubbles are positive residuals (observed > expected) and open bubbles are negative residuals (observed < expected).
19_caal_residuals_WCGBTS	Figure 38 [continued in Fig 39]: Pearson residuals for the fit to conditional age-at-length data from the fishery. Closed bubbles are positive residuals (observed > expected) and open bubbles are negative residuals (observed < expected).	stockplotr fig only has 3 out of 15 years
NA	Figure 40: Observed sex ratios (points) from the fishery length composition data with 75% intervals (vertical lines) calculated as a Jeffreys interval based on the adjusted input sample size. The model expectation is shown in the blue line.	stockplotr example only has fleet-specific plot for one fleet, but DOES have aggregated plot 14_sexratio_aggregated which is missing from r4ss and would have been even better to include
13_sexratio_WCGBTS	Figure 41: Observed sex ratios (points) from the WCGBT Survey length composition data with 75% intervals (vertical lines) calculated as a Jeffreys interval based on the adjusted input sample size. The model expectation is shown in the blue line.

[ericward-noaa]

@Schiano-NOAA @iantaylor-NOAA I will start working on these in a new PR once everything's merged in. A few things I'm thinking about in looking at existing files:

1. There's a few arguments the composition plots don't include, but probably need to be added:
   make_rda
   figures_dir
   others?

2. There's a lot of tags / dynamic labels in the figures that could be generated within the existing plot_composition_ functions . That seems easiest to be, but could also be done in standalone functions. I think I like doing it in the function, because if those tags are generated on the fly, it means these functions could be called via Shiny (the make_rda and figures_dir could stay too). But you all have been thinking about this more, so I'm open to suggestions!

[Schiano-NOAA]

@ericward-noaa that sounds great! I already opened an issue regarding the arguments you mentioned in #451 The "rda package" involved the make_rda and figures_dir argument. I don't think there anything other than that one.

In regards to # 2, I am not sure what you mean by tags but when the quantities within the captions and alt text are extracted, they are also stored in a csv file that could be used for other means (we link it to a function to format output for SIS)