Introduce the explicit two-sex model for the FIMS catch-at-age population dynamics

[szu-yun-ko]

Notes on the Project Workflow

This is part of the "Add sex-structure to FIMS statistical catch-at-age model" Google Summer of Code project. Earlier on in PR #1581 we introduced partitioning infrastructure into the FIMS codebase that can be used for storing partitioned quantities, including sex-specific values. While that PR is still under review, I realized that I can work on the two-sex models in parallel, and after both components (the infrastructure and the two-sex model) are complete, we can then wire them together while addressing issue #1587. I think this workflow can get the project going at a faster pace and hopefully leave us time for implementing the implicit two-sex model as well.

Summary

Introduce the explicit two-sex model for the FIMS catch-at-age population dynamics. In this formulation, life-history processes can differ by sex and evolve over time, while the model keeps a single pooled abundance state $N_{a,y}$ plus a dynamic female proportion $p_{a,y}^{\mathrm{fem}}$.
This issue covers the model formulation, state-transition logic, sex-specific life-history parameters, and user-facing model selection. It does NOT include demand-driven partitioned likelihood wiring (NeedPartitionedLandings, sex-specific expected landings/index). That work stays blocked on the partition infrastructure PR and will be a follow-up addressed in issue #1587.

Current Status

In PR #1581, we've already registered fleet-level observation quantities (e.g. landings/index_numbers_at_age and landings/index_weight_at_age by partition. These are built on the infrastructure added in PR #1581, including a default sex-partition initialized in CatchAtAge::Initialize().

Implementation Details

Model Formulation

State Variables:

• $N_{a,y}$: Total numbers at age $a$, year $y$
• $p_{a,y}^{\mathrm{fem}}$: Proportion female at age $a$, year $y$

Sex-specific life-history parameters:

Parameter	Current Codebase	Explicit Two-Sex
Natural mortality	$M_{a,y}$ via log_M	$M_{a,y}^{\mathrm{fem}}$, $M_{a,y}^{\mathrm{male}}$
Fishing mortality	$F_{f,a,y} = F_{f,y}, f_y, S_f(a)$	Sex-specific $S_f^{\mathrm{sex}}(a) \rightarrow F_{a,y}^{\mathrm{sex}}$
Total mortality	$Z_{a,y} = M_{a,y} + F_{a,y}$	$Z_{a,y}^{\mathrm{sex}} = M_{a,y}^{\mathrm{sex}} + F_{a,y}^{\mathrm{sex}}$

Transition dynamics:

$$ \begin{aligned} N^{\mathrm{fem}}_{a,y} &= p^{\mathrm{fem}}_{a-1,y-1} N_{a-1,y-1} e^{-Z^{\mathrm{fem}}_{a-1,y-1}} \\ N^{\mathrm{male}}_{a,y} &= \left( 1 - p^{\mathrm{fem}}_{a-1,y-1} \right) N_{a-1,y-1} e^{-Z^{\mathrm{male}}_{a-1,y-1}} \\ N_{a,y} &= N^{\mathrm{fem}}_{a,y} + N^{\mathrm{male}}_{a,y} \\ p^{\mathrm{fem}}_{a,y} &= \frac{N^{\mathrm{fem}}_{a,y}}{N_{a,y}} \end{aligned} $$

Initial conditions: users supply $p_{a,0}^{\mathrm{fem}}$ only (age-specific at the initial year). Later $p_{a,y}^{\mathrm{fem}}$ are computed by the update above. This differs from the sex-ratio-at-age and implicit models, where $p_a^{\mathrm{fem}}$ is fixed across years.

Architectural Approach

• Store sex-specific intermediate states in population-level partitioned derived quantities indexed by IndexLayout::i_stratum_age_year() (e.g. numbers_at_age_by_partition, mortality_M_by_partition, mortality_F_by_partition, mortality_Z_by_partition).
• Keep numbers_at_age as the pooled total $N_{a,y}$.
• Extend proportion_female storage to support year-varying values for explicit mode (e.g. proportion_female_at_age_year of size n_ages × n_years).
• Branch in CalculateNumbersAA, CalculateMortality, and related functions on a new sex_structure model option.
• Fleet-level *_by_partition containers (already registered) will later read from population partitioned state when demand wiring lands.

TODOs

Model option and API surface

• Define SexStructure enum / string constant (sex_ratio_at_age, implicit_two_sex, explicit_two_sex) on Population or CatchAtAge.
• Expose via PopulationInterface + src/rcpp_population.cpp + initialize_modules.R.
• Document option in roxygen and default to sex_ratio_at_age.
• Add validation: explicit mode rejects incompatible configurations with a clear error.

Register population-level partitioned life-history containers

• Register population derived quantities (size n_strata × n_years × n_ages):
  • numbers_at_age_by_partition
  • mortality_M_by_partition
  • mortality_F_by_partition
  • mortality_Z_by_partition
• Register year-varying proportion_female_at_age_year (size n_ages × n_years) for dynamic $p_{a,y}^{\mathrm{fem}}$.
• Wire dimensions in rcpp_models.hpp using existing IndexLayout.
  Note: fleet-level *_by_partition containers already exist, this step adds the upstream population state they will consume.

Sex-specific parameter inputs

• Natural mortality: allow sex-specific log_M. If the user still supplies a pooled log_M but wishes to use the explicit two-sex model, then we will duplicate the pooled values into the sex-specific log_M containers so that both female and male log_M use the pooled values, i.e. $M_{a,y}^{\mathrm{fem}}$ = $M_{a,y}^{\mathrm{male}} = M_{a,y}^{\mathrm{(pooled)}}$.
• Selectivity / F: sex-specific selectivity per flee, update CalculateMortality to fill mortality_F_by_partition.
• Growth / weight: sex-specific weight-at-age.
• Update Prepare() transformations for any new estimated parameters.

Explicit transition dynamics

• Implement explicit branch in CalculateNumbersAA with formulas above, including plus-group handling for both sexes before summing to $N_{a,y}$.
• Implement dynamic $p_{a,y}^{\mathrm{fem}}$ update after abundance calculation.
• Age-0 / recruitment: define how $p_{0,y}^{\mathrm{fem}}$ enters (initial user input at $y=0$; align recruitment sex-split rule with mentors).
• Update CalculateSpawningBiomass to use year-varying $p_{a,y}^{\mathrm{fem}}$ and sex-appropriate weight in explicit mode.

stratum_split_factor integration

• Update split-factor call sites to read dynamically updated $p_{a,y}^{\mathrm{fem}}$ when sex_structure == explicit_two_sex.
• Ensure split factors are consistent with numbers_at_age_by_partition / numbers_at_age after each time step.

[iantaylor-NOAA]

@szu-yun-ko, thank you for all your working moving the two-sex feature forward in FIMS!

Can you (or anyone else involved in this effort) help me understand the reasons for choosing the state variables as

• $N_{a,y}$: Total numbers at age $a$, year $y$
• $p_{a,y}^{\mathrm{fem}}$: Proportion female at age $a$, year $y$

as opposed to having the state variables be the numbers in each partition?:

• $N_{a,y}^{\mathrm{fem}}$ and $N_{a,y}^{\mathrm{male}}$ or more generally $N_{a,y,p}$

I think the storage requirements would be the same, but I presume there other considerations that I'm not aware of.

[nathanvaughan-NOAA]

Hey @iantaylor-NOAA the idea is that with this approach you are always tracking the aggregate numbers but can turn on the proportional division between sexes only when needed rather than tracking all combinations in all years/ages. For example, most assessments only have sex specific data in recent years so this approach is designed to allow the model to ignore sex in the early years. You are right that it wouldn't save any storage if the proportion was tracked for all ages/years. It seemed easier to manage the dynamic scaling this way.

[iantaylor-NOAA]

@nathanvaughan-NOAA, thanks for the quick reply. I'm don't see the lack of availability of sex-specific data in early years allowing us to ignore sex for the early years because the only reason to use a two-sex model is if there are differences in growth, mortality, or selectivity between sexes, and I think those differences would need to be modeled throughout the timeseries to accurately estimate spawning biomass throughout the time series. On the other hand, if we're thinking about partitioning numbers among seasons within the year, it's easier for me to imagine that there would be benefit to ignoring seasonal dynamics in the early years, so it seems reasonable to try out this approach for sexes first.

[JonBrodziak]

Here's a concise reply that builds on Ian's point while focusing on model structure and computational tradeoffs. I tend to agree with Ian that, if a two-sex model is being used because the sexes differ in growth, mortality, selectivity, or maturity, those dynamics should generally be represented consistently across the full assessment time horizon. Otherwise, the population dynamics change structurally through time, which could introduce unnecessary complexity and inconsistencies in quantities such as spawning biomass. From a computational perspective, I'm also not convinced the aggregate-then-split approach provides a performance advantage. On a floating-point operations basis, that approach still requires the aggregate population dynamics calculations, followed by the additional calculations to partition numbers by sex. In contrast, a conventional two-sex implementation performs the population dynamics directly for each sex. While this doubles the state dynamics calculations relative to a one-sex model, it avoids the additional bookkeeping associated with reconstructing sex-specific quantities from aggregate numbers and proportions. Overall, I prefer the simpler implementation of storing the sufficient statistics for the numbers-at-age state variables in separate arrays by sex. It provides direct access to numbers at age for each sex, rather than requiring two references—one to the aggregate numbers-at-age array and another to the sex proportion array. The resulting code is simpler to read, maintain, and extend, and the data flow is more transparent.

…

On Fri, Jul 10, 2026 at 10:15 AM Ian Taylor ***@***.***> wrote: *iantaylor-NOAA* left a comment (NOAA-FIMS/FIMS#1625) <#1625 (comment)> @nathanvaughan-NOAA <https://github.com/nathanvaughan-NOAA>, thanks for the quick reply. I'm don't see the lack of availability of sex-specific data in early years allowing us to ignore sex for the early years because the only reason to use a two-sex model is if there are differences in growth, mortality, or selectivity between sexes, and I think those differences would need to be modeled throughout the timeseries to accurately estimate spawning biomass throughout the time series. On the other hand, if we're thinking about partitioning numbers among seasons within the year, it's easier for me to imagine that there would be benefit to ignoring seasonal dynamics in the early years, so it seems reasonable to try out this approach for sexes first. — Reply to this email directly, view it on GitHub <#1625?email_source=notifications&email_token=ACVXIZQ5VZWP2ZE5D3SFC5T5EFFFZA5CNFSNUABFM5UWIORPF5TWS5BNNB2WEL2JONZXKZKDN5WW2ZLOOQXTIOJTHEYTSMBSHA2KM4TFMFZW63VKON2WE43DOJUWEZLEUVSXMZLOOSWGM33PORSXEX3DNRUWG2Y#issuecomment-4939190284>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ACVXIZSYDEUQEACGBVXQXGL5EFFFZAVCNFSNUABFKJSXA33TNF2G64TZHM2DKMRQGEZDAMBUHNEXG43VMU5TIOBVG4ZDOMJXGE4KC5QC> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

--

------------------------------ Jon Brodziak Mathematical Statistician Pacific Islands Fisheries Science Center NOAA Fisheries | U.S. Department of Commerce Office: (808) 725-5617 www.fisheries.noaa.gov | PIFSC Stock Assessment <https://pifscstockassessments.github.io/>

[nathanvaughan-NOAA]

@JonBrodziak the eventual goal is to extend what we develop here with sexes to other model dimensions like seasons/areas/species etc which provide much less contentious use cases for dynamic dimension tracking through time than if you look at this purely for tracking sex fractions. Like @iantaylor-NOAA mentioned regarding seasons.

[JonBrodziak]

I appreciate the broader motivation here. Scaling this up to seasons, areas, or species makes the intent behind dynamic dimension tracking much clearer, and it's a compelling goal for FIMS flexibility. My primary concern is how the generalized representation handles complexity once multiple dimensions interact. While the aggregate-plus-proportion approach reduces the number of primary state variables (saving memory), it effectively trades memory for code complexity. Tracking nested conditional proportions (e.g., sex within a season within an area by fleet) introduces heavy bookkeeping and lots of index reconstruction. From an architectural standpoint, explicit biological state variables, even if sparse, are often much easier to read, maintain, and extend. This the approach taken in the MAS project (so I have some cognitive bias favoring it as a preference). That said, it’s a promising concept, and I think it would be worth prototyping both approaches using a representative multi-dimensional example so we can directly compare code readability and performance.

…

On Fri, Jul 10, 2026 at 12:31 PM Nathan Vaughan ***@***.***> wrote: *nathanvaughan-NOAA* left a comment (NOAA-FIMS/FIMS#1625) <#1625 (comment)> @JonBrodziak <https://github.com/JonBrodziak> the eventual goal is to extend what we develop here with sexes to other model dimensions like seasons/areas/species etc which provide much less contentious use cases for dynamic dimension tracking through time than if you look at this purely for tracking sex fractions. Like @iantaylor-NOAA <https://github.com/iantaylor-NOAA> mentioned regarding seasons. — Reply to this email directly, view it on GitHub <#1625?email_source=notifications&email_token=ACVXIZSH4RMJ7DQFAEW5GM35EFVFLA5CNFSNUABFM5UWIORPF5TWS5BNNB2WEL2JONZXKZKDN5WW2ZLOOQXTIOJUGAYDKNBYHE42M4TFMFZW63VHNVSW45DJN5XKKZLWMVXHJLDGN5XXIZLSL5RWY2LDNM#issuecomment-4940054899>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ACVXIZUCIM25VBAWJDKFBOT5EFVFLAVCNFSNUABFKJSXA33TNF2G64TZHM2DKMRQGEZDAMBUHNEXG43VMU5TIOBVG4ZDOMJXGE4KC5QC> . You are receiving this because you were mentioned.Message ID: ***@***.***>

--

------------------------------ Jon Brodziak Mathematical Statistician Pacific Islands Fisheries Science Center NOAA Fisheries | U.S. Department of Commerce Office: (808) 725-5617 www.fisheries.noaa.gov | PIFSC Stock Assessment <https://pifscstockassessments.github.io/>