Semen Allocation

README.md

@@ -1,7 +1,7 @@
-[![Flake8](https://img.shields.io/badge/Flake8-passed-brightgreen)](https://github.com/RuminantFarmSystems/MASM/actions/workflows/combined_format_lint_test_mypy.yml)
-[![Pytest](https://img.shields.io/badge/Pytest-passed-brightgreen)](https://github.com/RuminantFarmSystems/MASM/actions/workflows/combined_format_lint_test_mypy.yml)
-[![Coverage](https://img.shields.io/badge/Coverage-99%25-brightgreen)](https://github.com/RuminantFarmSystems/MASM/actions/workflows/combined_format_lint_test_mypy.yml)
-[![Mypy](https://img.shields.io/badge/Mypy-1180%20errors-red)](https://github.com/RuminantFarmSystems/MASM/actions/workflows/combined_format_lint_test_mypy.yml)
+[![Flake8](https://img.shields.io/badge/Flake8-failed-red)](https://github.com/RuminantFarmSystems/MASM/actions/workflows/combined_format_lint_test_mypy.yml)
+[![Pytest](https://img.shields.io/badge/Pytest-failed-red)](https://github.com/RuminantFarmSystems/MASM/actions/workflows/combined_format_lint_test_mypy.yml)
+[![Coverage](https://img.shields.io/badge/Coverage-%25-red)](https://github.com/RuminantFarmSystems/MASM/actions/workflows/combined_format_lint_test_mypy.yml)
+[![Mypy](https://img.shields.io/badge/Mypy-1233%20errors-red)](https://github.com/RuminantFarmSystems/MASM/actions/workflows/combined_format_lint_test_mypy.yml)
 
 
 # RuFaS: Ruminant Farm Systems

RUFAS/biophysical/animal/animal.py

@@ -1219,7 +1219,7 @@ def _initialize_newborn_calf(self, args: NewBornCalfValuesTypedDict, simulation_
             The current day in the simulation, used for event logging and status evaluation.
 
         """
-        self._assign_sex_to_newborn_calf()
+        self.sex = args["sex"]
 
         if random() < AnimalConfig.still_birth_rate:
             self.stillborn_day = simulation_day
@@ -1562,7 +1562,9 @@ def _determine_days_in_milk(self, reproduction_output_days_in_milk: int) -> int:
             raise ValueError("Unexpected days in milk value")
 
     def daily_reproduction_update(
-        self, time: RufasTime
+        self,
+        time: RufasTime,
+        population_ranking_indexes: list[float] | None = None,
     ) -> tuple[NewBornCalfValuesTypedDict | None, HerdReproductionStatistics]:
         """
         Handles the daily reproduction state update for an animal.
@@ -1597,6 +1599,8 @@ def daily_reproduction_update(
                 dam_tbv_fat=self.genetics.TBV_fat,
                 dam_tbv_protein=self.genetics.TBV_protein,
                 phosphorus_for_gestation_required_for_calf=self.nutrients.phosphorus_for_gestation_required_for_calf,
+                population_ranking_indexes=population_ranking_indexes,
+                animal_ranking_index=self.genetics.ranking_index,
             )
         else:
             reproduction_inputs = ReproductionInputs(
@@ -1607,6 +1611,8 @@ def daily_reproduction_update(
                 days_in_pregnancy=self.days_in_pregnancy,
                 days_in_milk=self.days_in_milk,
                 phosphorus_for_gestation_required_for_calf=self.nutrients.phosphorus_for_gestation_required_for_calf,
+                population_ranking_indexes=population_ranking_indexes,
+                animal_ranking_index=self.genetics.ranking_index,
             )
         reproduction_outputs: ReproductionOutputs = self.reproduction.reproduction_update(reproduction_inputs, time)
 
@@ -1638,7 +1644,7 @@ def daily_reproduction_update(
 
         return newborn_calf_config, reproduction_outputs.herd_reproduction_statistics
 
-    def daily_routines(self, time: RufasTime) -> DailyRoutinesOutput:
+    def daily_routines(self, time: RufasTime, population_ranking_indexes: list[float] | None) -> DailyRoutinesOutput:
         """
         Perform daily routines for the animal, updating its status and outputs.
 
@@ -1668,7 +1674,9 @@ def daily_routines(self, time: RufasTime) -> DailyRoutinesOutput:
 
         self.daily_growth_update(time)
 
-        newborn_calf_config, daily_routines_output.herd_reproduction_statistics = self.daily_reproduction_update(time)
+        newborn_calf_config, daily_routines_output.herd_reproduction_statistics = self.daily_reproduction_update(
+            time, population_ranking_indexes
+        )
 
         daily_routines_output.animal_status, daily_routines_output.newborn_calf_config = self.animal_life_stage_update(
             time
@@ -2505,3 +2513,16 @@ def update_genetic_history(self, simulation_day: int) -> None:
                 self.genetic_history[-1]["end_day"] = simulation_day
         else:
             return
+
+    @property
+    def is_eligible_for_breeding(self) -> bool:
+        if self.animal_type is AnimalType.HEIFER_II:
+            return not self.is_pregnant
+        elif self.animal_type.is_cow:
+            return (
+                not self.is_pregnant
+                and self.days_in_milk > AnimalConfig.voluntary_waiting_period
+                and not self.reproduction.do_not_breed
+            )
+        else:
+            return False

RUFAS/biophysical/animal/animal_config.py

@@ -195,6 +195,16 @@ class AnimalConfig:
     semen_type: str = "conventional"
     male_calf_rate_conventional_semen: float = 0.53
     male_calf_rate_sexed_semen: float = 0.10
+    heiferII_semen_allocation_proportions: dict[str, float] = {
+        "sexed_dairy": 0.5,
+        "conventional_dairy": 0.5,
+        "beef": 0.0,
+    }
+    cow_semen_allocation_proportions: dict[str, float] = {
+        "sexed_dairy": 0.2,
+        "conventional_dairy": 0.5,
+        "beef": 0.3,
+    }
     keep_female_calf_rate: float = 1
     still_birth_rate: float = 0.065
     average_gestation_length: int = 276
@@ -395,10 +405,7 @@ def initialize_animal_config(cls) -> None:
         cls.do_not_breed_time = animal_config_data["management_decisions"]["do_not_breed_time"]
 
         cls.semen_type = animal_config_data["management_decisions"]["semen_type"]
-        cls.male_calf_rate_conventional_semen = animal_config_data["farm_level"]["calf"][
-            "male_calf_rate_conventional_semen"
-        ]
-        cls.male_calf_rate_sexed_semen = animal_config_data["farm_level"]["calf"]["male_calf_rate_sexed_semen"]
+        cls.semen_allocation = animal_config_data["management_decisions"]["semen_allocation"]
         cls.keep_female_calf_rate = animal_config_data["farm_level"]["calf"]["keep_female_calf_rate"]
         cls.still_birth_rate = animal_config_data["from_literature"]["life_cycle"]["still_birth_rate"]
         cls.average_gestation_length = animal_config_data["farm_level"]["repro"]["avg_gestation_len"]

RUFAS/biophysical/animal/data_types/animal_typed_dicts.py

@@ -2,6 +2,7 @@
 
 from typing_extensions import NotRequired
 
+from RUFAS.biophysical.animal.data_types.animal_enums import Sex
 from RUFAS.biophysical.animal.data_types.body_weight_history import BodyWeightHistory
 from RUFAS.biophysical.animal.data_types.pen_history import PenHistory
 
@@ -30,6 +31,7 @@ class NewBornCalfValuesTypedDict(TypedDict):
 
     id: NotRequired[int]
     breed: str
+    sex: Sex
     animal_type: str
     birth_date: str
     days_born: int

RUFAS/biophysical/animal/data_types/reproduction.py

@@ -1,4 +1,5 @@
 from dataclasses import dataclass
+from enum import Enum
 
 from RUFAS.biophysical.animal.data_types.animal_enums import Breed
 from RUFAS.biophysical.animal.data_types.animal_events import AnimalEvents
@@ -41,6 +42,8 @@ class ReproductionInputs:
     days_in_pregnancy: int
     days_in_milk: int
     phosphorus_for_gestation_required_for_calf: float
+    population_ranking_indexes: list[float]
+    animal_ranking_index: float
     dam_tbv_fat: float | None = None
     dam_tbv_protein: float | None = None
 
@@ -318,6 +321,8 @@ class ReproductionDataStream:
     phosphorus_for_gestation_required_for_calf: float
 
     herd_reproduction_statistics: HerdReproductionStatistics
+    population_ranking_indexes: list[float]
+    animal_ranking_index: float
     newborn_calf_config: NewBornCalfValuesTypedDict | None = None
     dam_tbv_fat: float | None = None
     dam_tbv_protein: float | None = None
@@ -329,3 +334,9 @@ def is_pregnant(self) -> bool:
     @property
     def is_milking(self) -> bool:
         return self.days_in_milk > 0
+
+
+class SemenType(Enum):
+    BEEF = "beef"
+    SEXED_DAIRY = "sexed_dairy"
+    CONVENTIONAL_DAIRY = "conventional_dairy"

RUFAS/biophysical/animal/herd_manager.py

@@ -3,6 +3,8 @@
 import math
 from typing import Any
 
+import numpy as np
+
 from RUFAS.biophysical.animal import animal_constants
 from RUFAS.biophysical.animal.animal import Animal
 from RUFAS.biophysical.animal.animal_config import AnimalConfig
@@ -22,7 +24,7 @@
 from RUFAS.biophysical.animal.data_types.animal_types import AnimalType
 from RUFAS.biophysical.animal.data_types.daily_routines_output import DailyRoutinesOutput
 from RUFAS.biophysical.animal.data_types.milk_production import MilkProductionStatistics
-from RUFAS.biophysical.animal.data_types.reproduction import HerdReproductionStatistics
+from RUFAS.biophysical.animal.data_types.reproduction import HerdReproductionStatistics, SemenType
 from RUFAS.biophysical.animal.herd_factory import HerdFactory
 from RUFAS.biophysical.animal.milk.lactation_curve import LactationCurve
 from RUFAS.biophysical.animal.milk.milk_production import MilkProduction
@@ -173,6 +175,7 @@ def __init__(
                 herd_population.cows,
                 herd_population.replacement,
             )
+            # TODO: randomly assign embryo sex for animals that are already pregnant: use Conventional Dairy
 
             self.allocate_animals_to_pens(time.simulation_day)
             self.initialize_nutrient_requirements(weather, time, available_feeds)
@@ -488,8 +491,14 @@ def _perform_daily_routines_for_animals(
         sold_newborn_calves: list[Animal] = []
         newborn_calves: list[Animal] = []
 
+        animal_ranking_indexes: list[float] | None = None
+        if all([(animal.animal_type == AnimalType.HEIFER_II or animal.animal_type.is_cow) for animal in animals]):
+            animal_ranking_indexes: list[float] | None = [
+                animal.genetics.ranking_index for animal in animals if animal.is_eligible_for_breeding
+            ]
+
         for animal in animals:
-            animal_daily_routines_output: DailyRoutinesOutput = animal.daily_routines(time)
+            animal_daily_routines_output: DailyRoutinesOutput = animal.daily_routines(time, animal_ranking_indexes)
             self.herd_reproduction_statistics += animal_daily_routines_output.herd_reproduction_statistics
             if animal_daily_routines_output.animal_status == AnimalStatus.DEAD:
                 self.herd_statistics.animals_deaths_by_stage[animal.animal_type] += 1

RUFAS/biophysical/animal/reproduction/reproduction.py

@@ -4,11 +4,12 @@
 from math import floor
 from typing import Callable, Union, Any, Optional
 
+import numpy as np
 from scipy.stats import truncnorm
 
 from RUFAS.biophysical.animal import animal_constants
 from RUFAS.biophysical.animal.animal_config import AnimalConfig
-from RUFAS.biophysical.animal.data_types.animal_enums import Breed
+from RUFAS.biophysical.animal.data_types.animal_enums import Breed, Sex
 from RUFAS.biophysical.animal.data_types.animal_typed_dicts import NewBornCalfValuesTypedDict
 from RUFAS.biophysical.animal.data_types.animal_types import AnimalType
 from RUFAS.biophysical.animal.data_types.preg_check_config import PregnancyCheckConfig
@@ -28,6 +29,7 @@
     ReproductionDataStream,
     AnimalReproductionStatistics,
     HerdReproductionStatistics,
+    SemenType,
 )
 
 from RUFAS.biophysical.animal.reproduction.hormone_delivery_schedule import HormoneDeliverySchedule
@@ -121,6 +123,8 @@ def __init__(
         do_not_breed: bool = False,
         estrus_count: int = 0,
     ) -> None:
+        self.semen_type: SemenType | None = None
+        self.embryo_sex: Sex | None = None
         self.heifer_reproduction_program = heifer_reproduction_program or AnimalConfig.heifer_reproduction_program
         self.heifer_reproduction_sub_program = (
             heifer_reproduction_sub_program or AnimalConfig.heifer_reproduction_sub_program
@@ -187,6 +191,8 @@ def reproduction_update(self, reproduction_inputs: ReproductionInputs, time: Ruf
             phosphorus_for_gestation_required_for_calf=reproduction_inputs.phosphorus_for_gestation_required_for_calf,
             herd_reproduction_statistics=HerdReproductionStatistics(),
             newborn_calf_config=None,
+            population_ranking_indexes=reproduction_inputs.population_ranking_indexes,
+            animal_ranking_index=reproduction_inputs.animal_ranking_index,
         )
         self.reproduction_statistics.reset_daily_statistics()
 
@@ -271,6 +277,9 @@ def cow_reproduction_update(
             Updated reproduction datastream for the cow.
         """
         if reproduction_data_stream.is_pregnant and reproduction_data_stream.days_in_pregnancy == self.gestation_length:
+            if time.simulation_day == 0:
+                # TODO: randomize sex for cows giving birth on day 0
+                self.embryo_sex = self._determine_embryo_sex(time.simulation_day)
             reproduction_data_stream = self.cow_give_birth(reproduction_data_stream, time)
 
         if not self.do_not_breed:
@@ -455,6 +464,7 @@ def cow_give_birth(
 
         reproduction_data_stream.newborn_calf_config = NewBornCalfValuesTypedDict(
             breed=reproduction_data_stream.breed.name,
+            sex=self.embryo_sex,
             animal_type=AnimalType.CALF.value,
             birth_date=time.current_date.strftime("%Y-%m-%d"),
             days_born=0,
@@ -463,6 +473,8 @@ def cow_give_birth(
             dam_tbv_fat=reproduction_data_stream.dam_tbv_fat,
             dam_tbv_protein=reproduction_data_stream.dam_tbv_protein,
         )
+        self.semen_type = None
+        self.embryo_sex = None
 
         return reproduction_data_stream
 
@@ -1005,6 +1017,11 @@ def _perform_ai(
         self, reproduction_data_stream: ReproductionDataStream, simulation_day: int
     ) -> ReproductionDataStream:
         """Perform artificial insemination (AI) on the animal."""
+        self.assign_semen_type(
+            population_ranking_index=reproduction_data_stream.population_ranking_indexes,
+            animal_ranking_index=reproduction_data_stream.animal_ranking_index,
+            animal_type=reproduction_data_stream.animal_type,
+        )
         reproduction_data_stream.events.add_event(
             reproduction_data_stream.days_born,
             simulation_day,
@@ -1013,7 +1030,7 @@ def _perform_ai(
         reproduction_data_stream.events.add_event(
             reproduction_data_stream.days_born,
             simulation_day,
-            animal_constants.INSEMINATED_W_BASE + AnimalConfig.semen_type,
+            animal_constants.INSEMINATED_W_BASE + self.semen_type.name,
         )
         self.reproduction_statistics.semen_number += 1
         self.reproduction_statistics.AI_times += 1
@@ -1035,7 +1052,9 @@ def _perform_ai(
                     reproduction_data_stream, simulation_day
                 )
                 reproduction_data_stream = self._increment_successful_cow_conceptions(reproduction_data_stream)
+            self.embryo_sex = self._determine_embryo_sex(simulation_day)
         else:
+            self.semen_type = None
             if reproduction_data_stream.animal_type == AnimalType.HEIFER_II:
                 reproduction_data_stream = self._handle_failed_heifer_conception(
                     reproduction_data_stream, simulation_day
@@ -1045,6 +1064,20 @@ def _perform_ai(
 
         return reproduction_data_stream
 
+    def _determine_embryo_sex(self, simulation_day) -> Sex:
+        # TODO: put male calf rate into Animal Constant
+        if self.semen_type == SemenType.CONVENTIONAL_DAIRY:
+            male_calf_rate = 0.5
+        elif self.semen_type == SemenType.SEXED_DAIRY:
+            male_calf_rate = 0.9
+        elif self.semen_type == SemenType.BEEF:
+            male_calf_rate = 0.9
+        else:
+            print(simulation_day, self.semen_type)
+            raise ValueError("Unexpected Semen Type.")
+        embryo_sex = Sex.MALE if random.random() < male_calf_rate else Sex.FEMALE
+        return embryo_sex
+
     def _increment_heifer_ai_counts(self, reproduction_data_stream: ReproductionDataStream) -> ReproductionDataStream:
         """Increment the AI counts for heifers."""
 
@@ -1123,6 +1156,9 @@ def _handle_failed_heifer_conception(
             average_estrus_cycle,
             std_estrus_cycle,
         )
+        # TODO: set sex to None when animal loses preg
+        self.semen_type = None
+        self.embryo_sex = None
         return reproduction_data_stream
 
     def _calculate_gestation_length(self) -> int:
@@ -1989,6 +2025,8 @@ def _handle_failed_cow_conception(
                         simulation_day,
                         f"Current repro state(s): {self.repro_state_manager}",
                     )
+        self.semen_type = None
+        self.embryo_sex = None
         return reproduction_data_stream
 
     def cow_pregnancy_update(
@@ -2194,3 +2232,23 @@ def _exit_ovsynch_program_early_when_first_preg_check_passed_or_estrus_detected(
             f" {AnimalConfig.cow_ovsynch_method}",
         )
         return reproduction_data_stream
+
+    def assign_semen_type(
+        self, population_ranking_index: list[float], animal_ranking_index: float, animal_type: AnimalType
+    ) -> None:
+        animal_ranking_index_percentile: float = (
+            np.mean(np.array(population_ranking_index) <= animal_ranking_index) * 100
+        )
+        semen_allocation_proportions = (
+            AnimalConfig.heiferII_semen_allocation_proportions
+            if animal_type == AnimalType.HEIFER_II
+            else AnimalConfig.cow_semen_allocation_proportions
+        )
+        if animal_ranking_index_percentile > (1 - semen_allocation_proportions["sexed_dairy"]):
+            self.semen_type = SemenType.SEXED_DAIRY
+        elif animal_ranking_index_percentile > (
+            1 - semen_allocation_proportions["sexed_dairy"] - semen_allocation_proportions["conventional_dairy"]
+        ):
+            self.semen_type = SemenType.CONVENTIONAL_DAIRY
+        else:
+            self.semen_type = SemenType.BEEF

RUFAS/task_manager.py

@@ -305,6 +305,7 @@ def check_dependencies(self) -> None:
                     " dependencies at required minimum levels.",
                     {"class": TaskManager.__name__, "function": TaskManager.check_dependencies.__name__},
                 )
+                return
                 raise RuntimeError(
                     f"Required package '{package_name}' version does not match. Installed: {installed_version}, "
                     f"Required: {requirement.specifier}"